ML weather and climate models

After this lesson you can

describe the forecast problem (a learned map from the state at $t$ to the state at $t + 6$ h, rolled out to 10–15 days) and the data it is trained on (ERA5);
name the main model families (graph networks on a mesh, 3D transformers, spherical Fourier neural operators, hybrid models, generative ensembles) and say where earlier lessons of this track enter each one;
explain rollout training, and why forecast errors grow with lead time and level off;
prove that the MSE-optimal forecast is the conditional mean, and use it to explain blurring and the double penalty;
compute latitude-weighted RMSE and ACC, the fair ensemble CRPS and the spread–skill ratio, and prove that a perfect $M$ -member ensemble has spread/skill $\approx \sqrt{M/(M+1)}$ ;
set up a fair comparison with the WeatherBench 2 protocol, and say what changes for climate emulation.

Before this lesson: lesson 05 (chaos and error growth), lesson 06 (4D-Var, which makes the initial states), lesson 07 (the CRPS, PIT and calibration) and lesson 09 (spherical harmonics and the SFNO). Helpful: lesson 03 (gradients through a solver), lesson 13 (diffusion models), and in the deep-learning track graph neural networks and transformers. Time: about 90 minutes, plus the notebook.

In plain words

A weather forecast starts from today's best estimate of the atmosphere and steps it forward in time. Classic forecasts solve the equations of physics on a supercomputer. An ML weather model learns the step from 40 years of past weather: it maps the state now to the state 6 hours later. Applied to its own output 40 times, it gives a 10-day forecast in about a minute.

These models now run every day next to the physics models at large weather centres. The lesson shows the main designs and how they are trained. It then spends most of its time on evaluation, because it is easy to get a good-looking score for the wrong reason. Two facts explain most surprises. First, a model trained to minimise the squared error learns the average of the possible futures, so its maps become smooth and blurry at long lead times. Second, a score computed on a latitude–longitude grid must weight each point by the area it covers. A fair comparison fixes the test year, the ground truth, the grid and the scores.

Why this matters #

ML runs in operational forecasting. As of 2026, ECMWF runs the Artificial Intelligence Forecasting System next to its physics model: AIFS Single became operational on 25 February 2025, and the 51-member ensemble AIFS ENS on 1 July 2025.
ML models match or beat physics on standard scores. GraphCast (Lam et al. 2023) beat the ECMWF high-resolution forecast on 90 % of 1380 verification targets; GenCast (Price et al. 2024) beat the ECMWF ensemble on more than 97 % of 1320 targets. These claims depend on how the evaluation is done, which is the second half of this lesson.
Climate needs long, stable runs. Hybrid models such as NeuralGCM (Kochkov et al. 2024) and emulators such as ACE (Watt-Meyer et al. 2023) run for decades, where new questions appear: stability, conservation and changing forcing.
The whole track meets here. The initial states come from 4D-Var (lesson 06). Errors grow because the atmosphere is chaotic (lesson 05). SFNO and NeuralGCM are spectral on the sphere (lessons 08, 09). Training through many steps is reverse-mode differentiation (lessons 01, 03). Ensembles come from diffusion models (lesson 13) and are scored with the CRPS (lesson 07).

1. Predict first #

2. The picture #

The forecast problem #

The state of the atmosphere is a stack of fields on the sphere. For a model like GraphCast it is 6 variables (temperature, wind components, humidity, geopotential, vertical velocity) on 37 pressure levels, plus 5 surface variables, at each of the $721 \times 1440$ points of a 0.25° latitude–longitude grid. That is about 235 million numbers. An ML weather model is a learned map

x_{t + 6\text{h}} = \Phi_\theta(x_t, x_{t - 6\text{h}}),

with weights $\theta$ . A forecast rolls out the map: it applies $\Phi_\theta$ to its own output again and again, 40 times for 10 days. The lead time is the time since the start of the forecast (the initialisation).

The training data is ERA5. ERA5 is ECMWF's reanalysis (Hersbach et al. 2020): one fixed version of the physics model, kept close to all available observations with 4D-Var (lesson 06), run over the past. As of 2026 it is hourly, 0.25°, on 37 pressure levels, from 1940 to a few days ago. It is the best estimate of the past atmosphere. It is not observations: in regions with few observations (the upper atmosphere over oceans, early decades) it is mostly the physics model, and its errors become the "truth" an ML model learns.

Model families #

All the families below keep the same outer loop: state on the grid in, next state on the grid out, rolled out. They differ in how they move information across the sphere.

Graph networks on a mesh (GraphCast; AIFS). An encoder sends the grid values to the nodes of a mesh that covers the sphere evenly (GraphCast: a refined icosahedron with 40,962 nodes). A processor updates the mesh many times, by message passing (GraphCast) or a sliding-window transformer (AIFS, Lang et al. 2024). A decoder sends the result back to the grid. The mesh avoids the crowding of grid points near the poles.
3D transformers (Pangu-Weather, Bi et al. 2023). The state is a 3D cube (pressure level × latitude × longitude) cut into patches. Attention works in local 3D windows with a learned bias that depends on the position on the Earth. Four models with steps of 1, 3, 6 and 24 hours are combined, so a long forecast needs fewer steps.
Spherical Fourier neural operators (SFNO, Bonev et al. 2023; the base of FourCastNet v2 and of the ACE climate emulator). Each layer transforms to spherical harmonics, multiplies each degree $l$ by learned weights, and transforms back (lesson 09). The layer is a convolution on the sphere, so it respects rotations; the paper reports stable rollouts of one year (1460 steps).
Hybrid models (NeuralGCM, Kochkov et al. 2024). A differentiable dynamical core solves the primitive equations of the atmosphere with a spectral method. A neural network, shared by all atmospheric columns, adds the tendencies of the processes that the grid does not resolve (clouds, convection). The model is trained end to end through the solver (lesson 03).
Generative ensembles (GenCast, Price et al. 2024; AIFS ENS). GenCast is a diffusion model (lesson 13): given the two last states, it samples one possible next state (12-hour steps) from noise. New noise gives a new ensemble member. AIFS ENS is instead trained directly with a CRPS loss (Lang et al. 2024, AIFS-CRPS).

Training: weighted MSE and rollouts #

A deterministic model is trained to minimise a weighted mean squared error between the predicted and the true next state. Each grid point is weighted by its area (section 3), each variable and pressure level by a chosen weight. After this one-step training, the model is fine-tuned on rollouts: it runs $k$ steps on its own output, and the loss adds the errors of all $k$ steps,

L_k(\theta) = \sum_{j=1}^{k} \big\| \hat x_{t + j} - x_{t + j} \big\|_W^2 ,

with $\hat x_t = x_t$ and $\hat x_{t + j} = \Phi_\theta(\hat x_{t + j - 1})$ (written with one input state for short). $\|\cdot\|_W^2$ is the weighted sum of squares.

The gradient flows back through all $k$ applications of $\Phi_\theta$ : reverse-mode differentiation through a time-stepper, as in lesson 03. GraphCast is fine-tuned with rollouts that grow from 2 to 12 steps (3 days). Why does this help? A model trained on one step never sees its own errors as input. In a rollout it does, and the loss teaches it to keep these errors small over several steps.

Errors grow, then level off #

The atmosphere is chaotic (lesson 05): a small initial error grows roughly like $e^{\lambda t}$ , with $\lambda$ the largest Lyapunov exponent. A learned model adds its own error at every step. The error cannot grow forever: when the forecast has lost all information, it is like a random state of the climate. Two independent states with climatological standard deviation $\sigma_c$ differ by $\sqrt 2\,\sigma_c$ in RMSE. The climatological mean itself has RMSE $\sigma_c$ , and an ensemble mean tends to it.

A toy atmosphere: Lorenz-96 with 40 variables (one step = 6 hours). Learned surrogates (orange) were trained on noisy data with a one-step loss (dashed) or a $k$ -step rollout loss (solid). Compare with the true equations started from a small initial error: one run (blue) and the mean of 20 runs (violet). Change $k$ : the rollout loss lowers the error over the first days. Move $\sigma_0$ : a better start helps the physics model for a while, then chaos wins. Notice where the curves level off: one run near $\sqrt2\,\sigma_c$ , the ensemble mean near $\sigma_c$ .

Quick check

A model is trained with a one-step loss and then fine-tuned with a 4-step rollout loss. What does the fine-tuning change?

3. Definitions and main results #

Notation. $f_{t,i,j}$ is a forecast of one variable at one lead time, for initialisation $t$ ( $t = 1, \dots, T$ ), latitude row $i$ ( $1, \dots, I$ ) and longitude $j$ ( $1, \dots, J$ ). $o_{t,i,j}$ is the ground truth at the matching valid time, and $c_{i,j}$ is a climatology (the long-term mean for that day of the year and time of day). An ensemble has members $x_1, \dots, x_M$ (here we use latitude $\phi$ , not the colatitude of lesson 09).

Area weights #

On a latitude–longitude grid every row has the same number of points, but a row near a pole covers a tiny area. The area of the band between latitudes $\phi^l$ and $\phi^u$ is $2\pi R^2(\sin\phi^u - \sin\phi^l)$ . WeatherBench 2 (WB2; Rasp et al. 2024, eq. 1) uses the latitude weights

w_i = \frac{\sin\phi_i^u - \sin\phi_i^l}{\frac1I\sum_{k}\big(\sin\phi_k^u - \sin\phi_k^l\big)}, \qquad \text{mean } 1 .

For a row of width $\Delta$ around $\phi_i$ ,

\sin(\phi_i + \tfrac\Delta2) - \sin(\phi_i - \tfrac\Delta2) = 2\sin(\tfrac\Delta2)\cos\phi_i ,

so $w_i \propto \cos\phi_i$ . Without weights, the many small polar cells dominate every mean.

Quick check

On the 0.25° grid, 242 of the 721 latitude rows (33.6 %) have $|\text{lat}| \ge 60°$ . What fraction of the Earth's area lies poleward of 60°?

Deterministic scores #

These are the definitions of WB2 §4.3.

RMSE (eq. 2):

\operatorname{RMSE} = \sqrt{\frac{1}{TIJ}\sum_{t,i,j} w_i\,(f_{t,i,j} - o_{t,i,j})^2}.

The time mean is inside the square root. (WeatherBench 1 took it outside; WB2 found the difference to be below 2 % on average.)

ACC, the anomaly correlation coefficient (eqs. 4–5): with anomalies $f' = f - c$ and $o' = o - c$ ,

\operatorname{ACC} = \frac1T\sum_t \frac{\sum_{i,j} w_i f'_{t,i,j}\, o'_{t,i,j}}{\sqrt{\sum_{i,j} w_i f'^2_{t,i,j}\ \sum_{i,j} w_i o'^2_{t,i,j}}}.

It asks whether the pattern of the anomalies is right, whatever its amplitude. The climatology is the ERA5 mean over 1990–2019 for each day of the year and time of day, smoothed over 61 days. ECMWF's rule of thumb (quoted in WB2): below an ACC of 0.6, the positions of large weather systems are no longer useful.

Bias maps (the time mean of $f - o$ ), and SEEPS for precipitation (a score on three categories: dry, light, heavy).
Energy spectra (WB2 §4.5): the power of the Fourier modes along latitude circles. A blurry forecast has too little power at high wavenumbers (small scales).

Theorem 1 (the MSE-optimal forecast is the conditional mean) #

Let $\mathcal I$ be all the information available at the initialisation, let the future truth $o$ be random given $\mathcal I$ , and let $\mu = \mathbb E[o \mid \mathcal I]$ . For every forecast $f$ that is a function of $\mathcal I$ :

\mathbb E\big[(f - o)^2 \mid \mathcal I\big] = (f - \mu)^2 + \operatorname{Var}(o \mid \mathcal I).

So $f = \mu$ minimises the expected squared error, and no forecast can go below the floor $\operatorname{Var}(o \mid \mathcal I)$ .

Proof.

Write $f - o = (f - \mu) + (\mu - o)$ .
Square and take $\mathbb E[\,\cdot \mid \mathcal I]$ . Since $f - \mu$ is known given $\mathcal I$ , the result is $(f - \mu)^2 + 2(f - \mu)\,\mathbb E[\mu - o \mid \mathcal I]$ $+\ \mathbb E[(\mu - o)^2 \mid \mathcal I]$ .
$\mathbb E[\mu - o \mid \mathcal I] = \mu - \mu = 0$ , so the middle term vanishes. The last term is $\operatorname{Var}(o \mid \mathcal I)$ by definition. $\blacksquare$

Consequences.

Blurring. As the lead time grows, $\operatorname{Var}(o \mid \mathcal I)$ grows, and $\mu$ is an average over many possible weather states. A front that may be here or 300 km further east is smeared over both places. A model trained with MSE learns this average: low RMSE, high ACC, smooth fields, too little small-scale energy. WB2 shows this blurring for the deterministic ML models.
Double penalty. A sharp, realistic forecast $o'$ , drawn from the same conditional distribution as $o$ but independent of it, has $\mathbb E[(o' - o)^2 \mid \mathcal I] = 2\operatorname{Var}(o \mid \mathcal I)$ : twice the floor. It is penalised once for missing the feature where it happened, and once for putting it where it did not happen.
So RMSE and ACC alone cannot rank a smooth model against a sharp one. Compare spectra, and compare ensembles with probabilistic scores.

Move the shift $d$ of the true feature (blue). Forecast A (orange, dashed) puts a sharp feature at the expected place. Forecast B (violet) is the conditional mean: the feature averaged over the possible positions, with uncertainty $s$ . Notice three things. For a clear miss, A has a larger RMSE than even "no feature". Averaged over the possible shifts, B always wins. In the spectra (second chart), A has exactly the spectrum of the truth, while B has too little power at high wavenumbers.

Quick check

The truth at a point is $o = \mu + \varepsilon$ , where $\mu$ is known at the start and $\varepsilon \sim \mathcal N(0, 1)$ is unpredictable. Which forecast has the smallest expected squared error?

Probabilistic scores #

CRPS of an ensemble. From lesson 07, $\operatorname{CRPS}(F, y) = \mathbb E|X - y| - \tfrac12\mathbb E|X - X'|$ with $X, X'$ independent draws from $F$ . For an ensemble $x_1, \dots, x_M$ , WB2 (eq. 8) uses the fair estimator (Zamo and Naveau 2018) at each grid point, then the weighted mean over space and time:

\operatorname{CRPS}_{\text{fair}} = \tfrac1M\sum_{m} |x_m - y| - \tfrac{1}{2M(M-1)}\sum_{m,n} |x_m - x_n| .

The standard estimator divides the second sum by $2M^2$ . It is the CRPS of the empirical distribution of the members (a CDF with $M$ steps).

Why "fair". The diagonal terms $m = n$ of the double sum are 0. Each of the other $M(M - 1)$ pairs has expectation $\mathbb E|X - X'|$ . So the fair estimator is unbiased for the CRPS of the distribution the members come from: a small ensemble is not punished for being small. The standard estimator keeps the zero diagonal in a sum divided by $2M^2$ , so it is too large on average by $\mathbb E|X - X'|/(2M)$ .
A degeneracy. For $M = 2$ the fair CRPS is 0 whenever $y$ lies between the two members, however far apart they are. AIFS-CRPS (Lang et al. 2024) therefore trains with an "almost fair" CRPS, which removes most of the bias and avoids this degeneracy.
Spread–skill ratio (WB2 eqs. 9–10): the spread is $\sqrt{\frac{1}{TIJ}\sum w_i \operatorname{var}_m(x)}$ , with the ensemble variance at each point (with $1/(M - 1)$ ). The skill is the RMSE of the ensemble mean, and the ratio is spread/skill. A ratio below 1 means the ensemble is under-dispersive (overconfident); above 1, over-dispersive. The rank histogram (the rank of $y$ among the members, $M + 1$ bins) is the ensemble version of the PIT histogram of lesson 07.

Theorem 2 (spread and skill of a perfect ensemble) #

Let the truth $y$ and the members $x_1, \dots, x_M$ be independent draws, given $\mathcal I$ , from the same distribution with mean $\mu$ and variance $\sigma^2$ . Let $\bar x$ be the ensemble mean. Then

\mathbb E\big[(\bar x - y)^2 \mid \mathcal I\big] = \Big(1 + \frac1M\Big)\sigma^2,

while the member variance with $1/(M - 1)$ has expectation $\sigma^2$ .

Proof.

$\bar x - y = (\bar x - \mu) - (y - \mu)$ . The two parts are independent and have mean 0.
$\operatorname{Var}(\bar x) = \sigma^2/M$ (a mean of $M$ independent draws) and $\operatorname{Var}(y) = \sigma^2$ . Variances of independent parts add: $\mathbb E(\bar x - y)^2 = \sigma^2/M + \sigma^2$ .
The sample variance with $1/(M - 1)$ is unbiased: its expectation is $\sigma^2$ . $\blacksquare$

So skill² $= \frac{M + 1}{M}$ spread² on average, and even a perfect ensemble has a spread–skill ratio of about $\sqrt{M/(M + 1)}$ : $0.953$ for 10 members, $0.990$ for 50. WB2 eq. 9 has no correction for this; Fortin et al. (2014) multiply the spread by $\sqrt{(M + 1)/M}$ . A small shortfall below 1 is expected; a large one is a real problem.

Change the ensemble size $M$ and the spread factor $s$ . Each of 3000 cases has a predictable part and unit noise; the members are drawn with spread $s$ ( $s = 1$ is the perfect ensemble). Left: the rank histogram. Right: the mean CRPS against $M$ , with the standard and the fair estimator (lines: expected values; dots: this sample). Try "perfect, M = 2": the fair CRPS stays at $1/\sqrt\pi$ , the standard one is much higher, and spread/skill is about $0.82$ . Then try "too narrow".

Quick check

Two members $(-1, 1)$ , outcome $y = 2$ . What are the fair and the standard CRPS?

Quick check

A perfect 10-member ensemble is evaluated on many cases. What spread–skill ratio (WB2 eq. 9, no correction) do you expect?

A fair comparison: the WeatherBench 2 protocol #

Choice	WB2 rule	Why
Test period	all 00 and 12 UTC initialisations in 2020 (732 forecasts)	held out from the training data of the ML models in the benchmark
Ground truth	ERA5 for models started from ERA5; for the ECMWF physics forecasts, the operational analysis (in practice the forecast at lead time 0)	a model verified against a different analysis than the one it starts from has an error even at lead time 0
Resolution	forecasts and truth regridded first-order conservatively to 1.5° ( $240 \times 121$ )	scores change with resolution and with the regridding method
Time convention	by initialisation time	"2020" means forecasts started in 2020
Metrics	RMSE, ACC, bias, SEEPS; CRPS and spread–skill for ensembles; spectra	each one sees a different failure

Rules for your own results:

Never train or tune on the test year. If your model has seen 2020, evaluate it on a later year, and evaluate the baselines on the same year.
Same ground truth, grid and regridding for every model in one plot.
Scores against lead time, not one number.
Test differences with paired scores. Score both models on the same initialisations and study the difference per initialisation. Forecasts started on nearby days are correlated, so resample whole blocks of days (a block bootstrap; see the bootstrap lesson), not single days.
Look beyond RMSE: spectra, bias maps, the CRPS for ensembles, case studies (a cyclone, a heat wave).
Watch the units: WB2 stores geopotential in m² s⁻². Geopotential height (m) times $g = 9.80665$ m s⁻² gives geopotential.

Weather versus climate #

A weather model is judged on single forecasts up to two weeks. A climate model is judged on statistics over decades: long-term means, variability, extremes and trends. This changes the problem.

Stability. A rollout of 100 years is about 146,000 steps of 6 hours. Small errors that do not matter at day 10 can make the model drift or blow up. Spherical layers (SFNO) and hybrid models help; ACE, an SFNO-based emulator of a physics atmosphere model, is stable for 100 years (Watt-Meyer et al. 2023), and NeuralGCM tracks global mean temperature over multiple decades with prescribed sea-surface temperatures.
Conservation. Over long runs, small leaks of mass, moisture or energy add up. ACE nearly conserves column moisture without an explicit constraint; other models add constraints.
Changing forcing. A climate projection asks what happens with more greenhouse gases than in any training year. That is out of distribution for a model trained on the past. Hybrid models keep the physics of the dynamical core; pure ML models must be tested on held-out forcing (for example a warmer or a historical sea-surface temperature data set they have not seen).

Beyond the atmosphere: ML in seismology #

The same ideas appear in seismology. Deep-learning phase pickers such as PhaseNet (Zhu and Beroza 2019) and EQTransformer (Mousavi et al. 2020) read three-component seismograms and output the probability of a P- or S-wave arrival at each time. They find many more small earthquakes than manual picking, and these denser catalogs show fault structures in more detail (Mousavi and Beroza 2022). The evaluation questions are the same as here: a strict hold-out of events and regions, proper scores for probabilistic outputs, and spectra to catch blurred waveforms. For imaging the Earth's interior, see lesson 04.

Answers to section 1.

RMSE cannot decide it. By Theorem 1, the smooth model is closer to the conditional mean, which minimises RMSE; the realistic model pays the double penalty. Compare spectra, and compare ensembles with the CRPS and the spread–skill ratio.
Rows with $|\text{lat}| \ge 60°$ are $121$ per hemisphere, $242$ of $721$ , so 33.6 % of the points. They cover $1 - \sin 60° = 13.4\,\%$ of the area. An unweighted mean counts the polar caps about 2.5 times too much.
Smaller. By Theorem 2, the RMSE of the mean is $\sqrt{1 + 1/10}\,\sigma$ and the spread is $\sigma$ , so spread/skill $\approx \sqrt{10/11} = 0.953$ .

4. Worked example #

(a) Area weights. Three bands of 60°: $[-90°, -30°]$ , $[-30°, 30°]$ , $[30°, 90°]$ .

$\sin\phi^u - \sin\phi^l$ : $-0.5 - (-1) = 0.5$ ; $0.5 - (-0.5) = 1$ ; $1 - 0.5 = 0.5$ .
Mean: $2/3$ . Weights: $w = (0.5, 1, 0.5)/(2/3) = (0.75, 1.5, 0.75)$ .
Suppose the squared errors are 4, 1 and 4 K² in the three bands (4 points each).
Unweighted: $\sqrt{(16 + 4 + 16)/12} = \sqrt3 = 1.732$ K.
Weighted: $\sqrt{(0.75 \cdot 4 + 1.5 \cdot 1 + 0.75 \cdot 4)/3} = \sqrt{7.5/3} = \sqrt{2.5} = 1.581$ K.

The unweighted score is pulled up by the polar bands, which cover only half of the area.

(b) Double penalty. A 1D strip of four cells; the truth has a feature in cell 2: $o = (0, 1, 0, 0)$ .

A sharp forecast one cell off, $(0, 0, 1, 0)$ : two cells are wrong by 1, MSE $= 2/4 = 0.5$ .
No feature at all, $(0, 0, 0, 0)$ : MSE $= 1/4 = 0.25$ .
A smeared forecast, $(0, 0.5, 0.5, 0)$ : two cells are wrong by $0.5$ , MSE $= 0.5/4 = 0.125$ .

If the feature is equally likely in cells 2 and 3, the smeared forecast is the conditional mean. Its MSE is $0.125$ in both cases. The sharp forecast at cell 2 has MSE 0 or $0.5$ , so $0.25$ on average: twice as much. The smeared forecast wins (Theorem 1), although it never looks like a real outcome.

(c) CRPS of a small ensemble. Members $(0, 2, 4)$ , outcome $y = 1$ .

First term: $\frac13(|0 - 1| + |2 - 1| + |4 - 1|)$ $= \frac13(1 + 1 + 3) = 5/3$ .
Pairs: $|0 - 2| + |0 - 4| + |2 - 4| = 8$ ; each pair counts twice, so $\sum_{m,n}|x_m - x_n| = 16$ .
Fair: $5/3 - 16/(2 \cdot 3 \cdot 2)$ $= 5/3 - 4/3 = 1/3$ .
Standard: $5/3 - 16/(2 \cdot 9)$ $= 15/9 - 8/9 = 7/9 \approx 0.778$ . This equals $\int (F_M(x) - \mathbb 1\{x \ge 1\})^2\, dx$ for the step CDF $F_M$ of the members: $\int_0^1 (\tfrac13)^2 dx + \int_1^2 (\tfrac23)^2 dx + \int_2^4 (\tfrac13)^2 dx$ $= \tfrac19 + \tfrac49 + \tfrac29 = \tfrac79$ .
The ensemble mean $2$ as a deterministic forecast: CRPS $=$ absolute error $= 1$ .

(d) Spread–skill of a perfect ensemble. $M = 10$ : $\sqrt{10/11} = 0.953$ . $M = 50$ : $\sqrt{50/51} = 0.990$ .

(e) Units. A Z500 RMSE of 35 m of geopotential height is $35 \times 9.80665 = 343$ m² s⁻² of geopotential.

5. In code #

NumPy checks the formulas on a toy world where everything is known. The truth is a predictable signal plus unpredictable noise of variance 1, so the conditional mean is the signal, and a perfect ensemble draws its members exactly like the truth.

import numpy as np

# WeatherBench 2 latitude weights (eq. 1) on the 1.5-degree grid (240 x 121, poles included)
lat = np.linspace(-90, 90, 121)
edges = np.concatenate([[-90], (lat[:-1] + lat[1:]) / 2, [90]])
w = np.diff(np.sin(np.deg2rad(edges)))
w /= w.mean()
print(w[[0, 30, 60]].round(4))            # pole row, 45 S, equator: [0.0052 1.1199 1.5838]

# A toy world: truth = predictable signal + unpredictable noise (variance 1 each)
rng = np.random.default_rng(0)
T, I, J, M = 8, lat.size, 240, 10
signal = rng.standard_normal((T, I, J))
obs = signal + rng.standard_normal((T, I, J))
ens = signal[:, None] + rng.standard_normal((T, M, I, J))   # a perfect ensemble
wb = w[None, :, None]

def rmse(f, o):                            # WB2 eq. 2: time mean inside the square root
    return np.sqrt(np.mean(wb * (f - o) ** 2))

def acc(f, o, c=0.0):                      # WB2 eqs. 4-5 (the climatology c is 0 here)
    fa, oa = f - c, o - c
    num = (wb * fa * oa).mean(axis=(1, 2))
    den = np.sqrt((wb * fa**2).mean(axis=(1, 2)) * (wb * oa**2).mean(axis=(1, 2)))
    return np.mean(num / den)

def crps(ens, o, fair=True):               # WB2 eq. 8 (fair) or the 1/(2M^2) version
    M = ens.shape[1]
    skill = np.abs(ens - o[:, None]).mean(axis=1)
    pair = sum(np.abs(ens[:, m:m + 1] - ens).sum(axis=1) for m in range(M))
    return np.mean(wb * (skill - pair / (2 * M * ((M - 1) if fair else M))))

print(rmse(signal, obs).round(3), rmse(ens[:, 0], obs).round(3))   # 1.0 and 1.416: the double penalty
print(acc(signal, obs).round(3), acc(ens[:, 0], obs).round(3))     # 0.707 and 0.499
print(crps(ens, obs).round(3), crps(ens, obs, fair=False).round(3),
      np.mean(wb * np.abs(signal - obs)).round(3))                 # 0.564 0.621 0.798
spread = np.sqrt(np.mean(wb * ens.var(axis=1, ddof=1)))
print((spread / rmse(ens.mean(axis=1), obs)).round(3), np.sqrt(M / (M + 1)).round(3))   # 0.953 0.953

The numbers repeat section 3. The blurry conditional mean wins on RMSE and ACC. A single sharp member pays the double penalty ( $\sqrt2$ instead of 1). The ensemble wins on the CRPS: $0.564 = 1/\sqrt\pi$ , the CRPS of the true predictive distribution $\mathcal N(\text{signal}, 1)$ , against $0.798 = \sqrt{2/\pi}$ for the conditional mean. A perfect 10-member ensemble has spread/skill $\sqrt{10/11}$ , not 1.

6. Where it shows up #

Operational centres. As of 2026, AIFS Single and AIFS ENS run at ECMWF, built with the open-source framework Anemoi, which ECMWF develops with national weather services in Europe. GraphCast, GenCast and NeuralGCM come from Google. New models are reported with WB2-style scores on the WeatherBench 2 website.
Climate emulation. ACE and NeuralGCM run multi-decade atmosphere simulations at a small fraction of the cost of a physics model. Climate evaluation looks at long-term means, trends and variability, not at day-5 RMSE.
Regional and short-range models. The same encoder–processor–decoder design works on limited areas with a coarser global model at the boundary.
Seismology. Deep-learning phase pickers (PhaseNet, EQTransformer) build earthquake catalogs. The evaluation logic of this lesson applies to ground-motion emulators too: proper scores (the CRPS of peak ground acceleration), spectra of generated waveforms to catch blurring, and a strict hold-out of events.
Beyond weather. Any learned simulator that is rolled out (ocean models, turbulence, plasma) faces the same blurring, error growth and evaluation questions.

7. Common mistakes #

Unweighted means on a latitude–longitude grid. The polar rows dominate. Use the area weights $w_i$ .
Different resolutions or regridding. Bilinear and conservative regridding give different scores. Regrid every model, and the truth, the same way.
The wrong ground truth. Verify each model against the analysis it was started from (ERA5 for ERA5-started models, the operational analysis for operational forecasts).
Leakage. Training or tuning on the test year, or normalisation statistics computed on it.
RMSE alone for a blurry model. A lower RMSE can mean a smoother forecast, not a better one. Add spectra; for ensembles, the CRPS and the spread–skill ratio.
The standard CRPS for small ensembles. It punishes small $M$ . WB2 uses the fair estimator.
Expecting spread/skill = 1 for a small ensemble. A perfect ensemble gives $\sqrt{M/(M + 1)}$ .
Initialisation time versus valid time. WB2 selects forecasts by initialisation time.
Comparing with WeatherBench 1 numbers. WB1 took the time mean outside the square root and used other grids and years.
Units. Geopotential (m² s⁻²) and geopotential height (m) differ by $g$ .

8. Check yourself #

A grid has 5 latitude bands of 36° each. Give the WB2 weights.
Answer
The edges are $-90°, -54°, -18°, 18°, 54°, 90°$ with sines $-1, -0.809, -0.309, 0.309, 0.809, 1$ . The differences are $0.191, 0.5, 0.618, 0.5, 0.191$ (sum 2, mean 0.4). Normalised: $(0.477, 1.25, 1.545, 1.25, 0.477)$ .
Why does an MSE-trained model blur more at day 10 than at day 2?
Answer
It learns the conditional mean $\mathbb E[o \mid \mathcal I]$ (Theorem 1). The conditional variance grows with lead time, so the mean averages over more possible states. Small features whose position is uncertain are smoothed out.
Members $(1, 3)$ , outcome $y = 2$ . Compute the fair and the standard CRPS. Why is the fair one 0?
Answer
$\frac12(1 + 1) = 1$ and $\sum_{m,n}|x_m - x_n| = 4$ . Fair: $1 - 4/(2 \cdot 2 \cdot 1) = 0$ . Standard: $1 - 4/8 = 0.5$ . The fair CRPS is 0 because $y$ lies between the two members: the degeneracy that the "almost fair" CRPS of AIFS-CRPS avoids.
A 10-member ensemble has a spread–skill ratio of 0.7. What does it mean, and what value is expected for a perfect 10-member ensemble?
Answer
The ensemble is under-dispersive: its spread is too small for its error, so it is overconfident. A perfect ensemble gives about $\sqrt{10/11} = 0.953$ . Remedies: CRPS training, stochastic perturbations, or inflation (lesson 06).
Your model was trained on ERA5 from 1979 to 2020, and you evaluate it on 2020 with WB2. What is wrong, and how do you fix it?
Answer
Test-set leakage: 2020 is in the training data, so the scores are too optimistic and not comparable with the WB2 baselines. Retrain without 2020 (for example train to 2018, validate on 2019), or evaluate on a later held-out year and score the baselines on that year too.
Why does a single forecast level off at an RMSE of about $\sqrt2\,\sigma_c$ at very long lead times, but an ensemble mean at about $\sigma_c$ ?
Answer
At very long lead times the single forecast is like an independent random state of the climate: $\mathbb E(a - b)^2 = 2\sigma_c^2$ for two independent states. The ensemble mean averages many such states and tends to the climatological mean, whose RMSE is $\sigma_c$ .

Key takeaways

An ML weather model learns the map $x_t \to x_{t + 6\text{h}}$ from reanalysis data (ERA5, made with 4D-Var) and rolls it out to 10–15 days; families differ in how they move information on the sphere (mesh GNNs, 3D transformers, SFNO, hybrid physics–ML, diffusion ensembles).
Rollout fine-tuning trains on several autoregressive steps with the gradient through all of them, so the model learns to handle its own errors; it cannot remove chaotic error growth.
Forecast errors grow with lead time and level off: a single forecast near $\sqrt2\,\sigma_c$ , an ensemble mean near the climatological $\sigma_c$ .
Theorem 1: $\mathbb E[(f - o)^2 \mid \mathcal I] = (f - \mu)^2 + \operatorname{Var}(o \mid \mathcal I)$ , so MSE training produces the smooth conditional mean, and a sharp realistic forecast pays the double penalty (twice the floor).
On a latitude–longitude grid, weight each row by its area, $w_i \propto \sin\phi_i^u - \sin\phi_i^l \approx \cos\phi_i$ ; poleward of 60° are 33.6 % of the rows of the 0.25° grid but 13.4 % of the area.
The fair ensemble CRPS $\frac1M\sum|x_m - y| - \frac{1}{2M(M-1)}\sum_{m,n}|x_m - x_n|$ is unbiased for the CRPS of the member distribution; the standard $1/(2M^2)$ version punishes small ensembles.
Theorem 2: for a perfect $M$ -member ensemble, skill² $= \frac{M+1}{M}$ spread², so spread/skill $\approx \sqrt{M/(M+1)}$ (0.953 for $M = 10$ ).
A fair comparison fixes the test year (WB2: 2020), the ground truth, the grid and regridding, the time convention and the units, and tests differences with paired scores and a block bootstrap.

Recall cardsAnswer in your head first, then reveal. Add the cards to your review deck to see them again just before you would forget.

What is the WeatherBench 2 latitude weight of row $i$ ?
State the decomposition behind "the MSE-optimal forecast is the conditional mean".
What is the double penalty?
How do RMSE and ACC differ?
Write the fair CRPS of an ensemble $x_1, \dots, x_M$ for outcome $y$ .
Compute the fair and standard CRPS for members $(0, 2, 4)$ and $y = 1$ .
What spread–skill ratio does a perfect 10-member ensemble have, and why?
Why fine-tune a weather model on multi-step rollouts?
What is ERA5, and what is it not?
Which year does WeatherBench 2 use for testing, and what follows for training?

Open the review →

9. Interview angle #

"How would you evaluate a new global ML weather model?" — The WB2 protocol: a held-out test year, 00 and 12 UTC initialisations, conservative regridding to a common grid, latitude-weighted RMSE and ACC per lead time against the right ground truth, spectra and bias maps, the CRPS and the spread–skill ratio for ensembles, and paired, block-bootstrapped differences.
"Your model beats the physics model on RMSE at day 10. Is it better?" — Not necessarily. Check the verification target and the resolution, and check blurring with spectra. At day 10 the conditional mean is smooth; the fair comparison is an ensemble against the physics ensemble, with the CRPS.
"Why the fair CRPS?" — It is unbiased for the CRPS of the distribution the members come from, so ensembles of different sizes can be compared.
"Hybrid (NeuralGCM) or pure ML (GraphCast, AIFS)?" — Hybrid: physical structure, multi-decade climate runs, interpretable parts, but a more complex model and coarser grids so far. Pure ML: simple to scale and very accurate at medium range at 0.25°, but long-term stability, conservation and behaviour under new forcing are not guaranteed.
"Why do deterministic ML forecasts look smooth?" — MSE training learns the conditional mean (Theorem 1). Remedies: generative or CRPS-trained ensembles, spectral losses, or reporting ensembles instead of single forecasts.

10. Practice #

The notebook for this lesson runs in Colab and in the browser. It has five exercises, each with a check cell and a folded solution:

Compute the WB2 latitude weights from band edges, and the polar over-count of an unweighted mean.
Write latitude-weighted RMSE and ACC; reproduce the toy-world numbers (1 and $\sqrt2$ ; 0.707 and 0.5).
Write the fair and the standard ensemble CRPS, and measure the bias of the standard one for small $M$ .
Compute the spread–skill ratio of a perfect ensemble, with and without the $\sqrt{(M + 1)/M}$ correction.
Fit a Lorenz-96 surrogate by one-step least squares on noisy data, then with a two-step rollout loss, and compare their error growth with lead time.

11. Watch and read #

The history of ML weather prediction and the design of the WeatherBench 2 benchmark: the protocol and scores of section 3.