Why Earthquake Prediction Is Impossible but Aftershock Forecasting Works: Measured on 130 Noto Aftershocks

Yesterday, on July 28, 2026 at 16

JST, an M7.1 earthquake hit Kumamoto, Japan (JMA magnitude, 16 km deep). An M6.1 followed 41 minutes later, and aftershocks were still arriving this morning. From my terminal it looks like this:

$ quake-lens recent --limit 4
time                      lat       lon   depth   mag  src   place
2026-07-29T04:17:00Z   32.700   130.800    10.0   2.1  p2p   熊本県熊本地方
2026-07-29T04:13:00Z   32.500   130.600    10.0   3.5  p2p   熊本県熊本地方
2026-07-29T04:09:00Z   32.500   130.600    10.0   3.1  p2p   熊本県熊本地方
2026-07-29T03:52:00Z   32.500   130.600    10.0   2.9  p2p   熊本県熊本地方

Nobody can tell you, looking at this list, whether a bigger one is coming. That is the official position of seismology. The USGS states in its FAQ that no scientist has ever predicted a major earthquake and that they do not expect to know how any time in the foreseeable future. Japan reached the same conclusion institutionally: in 2017, disaster planning for the Nankai Trough was rebuilt on the premise that reliable short-term prediction is not feasible.

And yet the same field, right after a large earthquake, states without hesitation: “expect earthquakes of similar size for about the next week.” The discipline that declared prediction impossible somehow gets to talk about the coming week.

There is no contradiction. “Prediction” and “forecasting” are two different words. In this post I’ll sort out the difference, then actually compute the forecasting side on 130 aftershocks of the 2024 Noto Peninsula earthquake, forecast the aftershock count for day 8, and check the answer. Everything runs on the Python standard library alone.

Prediction and forecasting are different words

Prediction means calling an individual earthquake in advance: when, where, and how big. This is the deterministic version, and it has never been scientifically established. Radon levels, animal behavior, electromagnetic noise: candidate precursors have been studied for decades, and not one reproducible prediction method survived.

Forecasting means making statistical statements about populations of earthquakes. You cannot answer “will an M6 hit tomorrow,” but you can answer “how much will this aftershock rate decay over a week” and “roughly what fraction of M4+ events will be M5+.” It is the world of probabilities and rates.

For engineers there is a one-line translation. Nobody can predict the exact second the next HTTP request arrives at your server, yet you can forecast tomorrow’s peak request rate accurately enough for capacity planning. Same thing here. Aftershock arrivals are in fact modeled as a non-homogeneous Poisson process, the same mathematics you would use on traffic analytics.

Individual events are unreadable. Population rates are readable. Seismology only speaks with confidence on the side where the law of large numbers is playing for its team.

The tools: two empirical laws, both about a century old

The machinery behind aftershock forecasting is surprisingly small: two short formulas.

The first is the Gutenberg-Richter law. The number N of earthquakes at or above magnitude M follows

log10 N = a − bM

The slope b is close to 1.0 worldwide. Each step down in magnitude means roughly 10 times more earthquakes: a region with ten M5 events gets about one M6. It is a conversion table across magnitudes.

The second is the Omori-Utsu law. The aftershock rate n(t) at t days after the mainshock decays as a power law:

n(t) = K / (t + c)^p

Fusakichi Omori published the original form in 1894, and one of the aftershock records he worked from was the 1889 Kumamoto earthquake. The tool we are pointing at yesterday’s Kumamoto M7.1 was born from Kumamoto data 130 years ago. With Tokunji Utsu’s 1961 modification, it still underpins the Japan Meteorological Agency’s official aftershock probability estimates today.

Notice that neither formula says anything about when a specific earthquake will occur. They only describe frequencies and rates, properties of the population.

Measuring 130 aftershocks of the Noto earthquake

A seismic observation station on the coast overlooking the Sea of Japan, with an aftershock waveform crossing the water (illustration)

Staring at formulas is no fun, so let’s estimate from real data. Pulling the USGS public catalog for 10 days after the Noto Peninsula earthquake (2024-01-01 16

JST, M7.5), magnitude 4 and above, gives 131 events; excluding the mainshock leaves 130 aftershocks.

First the b-value. Aki’s (1965) maximum likelihood estimator has a closed form that needs only the mean magnitude.

$ quake-lens catalog --start 2024-01-01 --end 2024-01-11 --min-mag 4.0 \
    --format json > noto.json
$ quake-lens bvalue noto.json --mc 4.5
b       = 1.0449
se      = 0.1306
n_used  = 64
mc      = 4.50

b = 1.04 ± 0.13, right in line with the global average of 1.0.

What is this good for? Use it as the conversion table. The estimate used 64 events at M4.5+. Plugging b = 1.045 into Gutenberg-Richter predicts 64 × 10^(−1.045×0.5) ≈ 19.2 events at M5.0+ and 5.8 at M5.5+. The observed counts were 16 and 7. Within a few events, two parameters capture the entire frequency structure across magnitudes. “Given this many M4-class events, about this many M5.5-class events should be mixed in”: one formula buys you that estimate.

The b-value itself is monitored in seismology. A low b means proportionally more large earthquakes relative to small ones, which is read as a signal of concentrated stress. But this too is population statistics; it does not translate into “when is this fault due.”

There is one pothole everyone hits when they try this for real. The completeness magnitude Mc is the threshold above which the catalog is assumed to record everything. Lower it from 4.5 to 3.0 on the same data and you get b = 0.27, nearly a factor of four off the global average. The culprit is not the data but the setting: the USGS global catalog does not fully capture Japan’s M3-class events. Small earthquakes that happened but were never recorded inflate the mean magnitude, and the estimator breaks. Adding “more” data (lowering Mc) makes the result worse. It is practically a textbook exercise in questioning your assumptions before your statistics.

Next, the decay. Fitting Omori-Utsu with Ogata’s (1983) maximum likelihood method reduces to a 2-D optimization over (c, p).

$ quake-lens omori noto.json --mainshock 2024-01-01T07:10:00Z
K       = 23.3856
c       = 0.0224
p       = 0.8861
n_used  = 130
window  = [0.0000, 8.8356] days

p = 0.89, the exponent controlling how fast the sequence dies down; typical values sit around 1.0. Running rates through the fitted model: one hour after the mainshock, M4+ aftershocks were arriving at a pace of 267 per day; one day later, 23 per day; by day 8, 3.7 per day. Half of everything happens in the first day, and the observations agree: 67 of the 130 aftershocks landed in day one.

Where the model holds, and where it breaks

Here are observed daily counts against model expectations:

Day	Observed	Model expected
0–1	67	72.6
1–2	16	16.6
2–3	15	10.4
3–4	13	7.7
4–5	8	6.2
5–6	4	5.2
6–7	3	4.4

The first two days match almost exactly, and the weekly shape tracks well after that. The excess on days 2 through 4 comes from larger aftershocks recruiting aftershocks of their own (secondary aftershocks), which a plain Omori fit flattens into a single decay curve.

But does it predict the future? I refit on the first 7 days only (K=24.0, c=0.019, p=0.85), forecast the count for day 8, and then checked.

Forecast: 4.3 events. Observed: 0.

A miss. Worse, the following 0.84 days produced 4 events in a burst. This is the honest face of statistical forecasting. Drawing 0 from a Poisson distribution with expectation 4.3 has probability of about 1.4%: inside the model’s support, but out on the edge. Aftershocks do not taper politely; they arrive in clumps. This is exactly why operational practice extends the model to ETAS, which treats secondary aftershocks explicitly.

The way it missed is the point. A statistical forecast is not the statement “day 8 will have 4.3 events.” It is the statement “day 8 is drawn from a distribution with expectation 4.3.” That is usable for a week-scale view of the rate and useless for calling the count on a specific day. When the JMA says “expect aftershocks for about the next week” rather than something sharper, that is not evasion. It is the precise extent of what this mathematics can say.

Multiply the two laws and you get “aftershock probability”

So far Gutenberg-Richter (the magnitude axis) and Omori-Utsu (the time axis) were validated separately. In practice you multiply them. Omori gives the expected total number of aftershocks over, say, the next 3 days; Gutenberg-Richter apportions what fraction of those will be M5 or larger; a Poisson distribution converts the expectation into “the probability of at least one M5+ aftershock in the next 3 days is X%.” This framework is the Reasenberg-Jones model, and it is more or less what sits behind the aftershock probabilities you see in the news.

Running it on the Noto fit: the model expects 19.0 aftershocks (M4+) over days 3 through 6. Apportioning with b = 1.045 gives 1.7 events at M5+. Through the Poisson distribution, the probability of at least one is 82%. So on day 3 after the mainshock you could have announced: “the probability of an M5-class aftershock within the next 3 days is about 80%.” Checking the answer: those 3 days produced 25 events at M4+ and exactly 1 at M5+. The headline claim verified, and the expected count was in a reasonable range.

The raw material was a few hours of aftershock data. From there, two formulas and one probability distribution carry you all the way to a sentence you could broadcast as public safety information. That is the entire pipeline of the forecasting side.

So why is prediction impossible?

If aftershock rates obey formulas this well, it feels like mainshocks should be readable too. Three reasons they are not.

First, the input is unobservable. What determines earthquake initiation is the stress state on a fault plane kilometers to tens of kilometers underground, and no sensor can be placed there. Drilling reaches a few kilometers at best. You cannot predict individual events of a system whose state you cannot see.

Second, rupture growth is sensitive to initial conditions. An earthquake is a small rupture that propagates across the fault surface, and whether it stops at M4 or grows to M7 depends on fine irregularities in friction and stress along the way. Even watching one start does not tell you how big it will get.

Third, precursors do not reproduce. Reports of “X happened before that earthquake” are abundant. Test the forward direction, “when X happens, an earthquake follows,” and you drown in false positives. In hindsight, everything looks like a precursor.

The serious experiments have been run. The famous one is Parkfield, California, where the USGS extrapolated a strikingly regular event history to predict an M6 by 1993 and blanketed the segment with instruments to wait for it. The M6 arrived in 2004, eleven years late. The instruments were still running, and even then no short-term precursor was captured. Japan’s version: the 1978 Large-Scale Earthquake Countermeasures Act was built on the premise that the anticipated Tokai earthquake could be predicted just before it struck. In 2017, about 40 years later, that premise was officially withdrawn. Four decades of dense monitoring and budget produced the conclusion: do not build disaster response on prediction.

Aftershock statistics work for the opposite reason: one mainshock hands you a population of hundreds to thousands of aftershocks. Individuals stay unreadable, but the moment they form a population, the law of large numbers kicks in. “Prediction is impossible, forecasting works” is nothing more or less than the boundary between determinism and probability.

Pointing the same pipeline at the ongoing Kumamoto sequence

As I write this, 0.9 days have passed since yesterday’s M7.1. The same pipeline runs unchanged on a sequence that is still in progress.

$ quake-lens recent --src jma --limit 400 --format json \
    | quake-lens omori - --mainshock 2026-07-28T07:27:15Z
K       = 197.5693
c       = 0.5990
p       = 1.7686
n_used  = 194
window  = [0.0000, 0.9123] days

The fit runs on 194 aftershocks (via the JMA earthquake list), but the resulting p = 1.77 and c = 0.60 are clearly inflated compared to Noto, which is what unstable estimation on less than a day of data looks like. The large c reflects detection saturation right after the mainshock still occupying most of the window. The b-value also comes out low at 0.60, but that one is a catalog completeness problem of the same species we saw at Noto: the JMA list only contains felt earthquakes (JMA seismic intensity 1 or higher), so unfelt M2-class events are missing.

Refit after a few more days of data and the parameters should settle toward stable values. To be clear, publishing probability numbers about an ongoing earthquake is not a job for a personal blog; for actual safety decisions, refer to the JMA and local government announcements. What this post offers is a map of what is being computed behind those announcements.

Summary

Prediction (deterministically calling individual events) is not scientifically established: the input is unobservable, rupture growth is sensitive, and precursors do not reproduce
Forecasting (population statistics) is operational on two century-old empirical laws. On 130 Noto aftershocks, stdlib-only Python reproduces b = 1.04 ± 0.13 and p = 0.89
A statistical forecast returns a distribution, not a point. Forecast 4.3 for day 8, observed 0, then a 4-event burst the next day: that is the demonstration
One Mc setting can break the b-value from 1.04 to 0.27. Before the statistics, interrogate catalog completeness

The CLI used throughout (quake-lens) is a zero-dependency Python tool written with urllib and math. It also runs as an MCP server, so wired into Claude, “get me the b-value for Noto” runs this entire computation from one sentence. The code is on GitHub: kenimo49/quake-lens

Next time a big one hits and the news says “expect aftershocks for about the next week,” these two formulas and a maximum likelihood fit are standing behind that sentence. Say honestly that prediction is impossible, and say only what can be forecast. I keep thinking this line seismology draws is the same one you draw in a postmortem: “we cannot predict the time of recurrence, but here is our estimate of the recurrence rate.”