How to measure a rain forecast

Hits, misses and false alarms

Pick a rain threshold, say 1 mm/h. Lay the forecast over what really happened. Every place and time falls in one of four boxes:

• hit: rain forecast, and it rained;
• miss: it rained, but no rain was forecast;
• false alarm: rain forecast, and it stayed dry;
• correct dry: no rain forecast, none fell.

Correct dry cases are left out of most rain scores on purpose: on most days most of Europe is dry, and a forecast of “no rain anywhere” would look excellent.

The scores built from them

Read the hit score and the bias together. A forecast can raise its detection by spreading rain everywhere (more hits, more false alarms), or look safe by forecasting rain rarely. A bias far from 1 tells you which trick is at work.

One picture for all of them

The performance diagram puts these scores on one chart. Across: the success ratio. Up: detection. A perfect forecast sits in the top right corner. The curved lines are equal hit scores, the dashed lines through the origin equal bias. Since hit score and bias fix the other two, any published pair of them can be placed on it.

In [1]:

# from CSI and bias to detection and success ratio
pod = csi * (1 + bias) / (1 + csi)
sr = pod / bias
plt.plot(sr, pod, marker="o")  # one point per time ahead

Out[1]:

The three forecasts slide toward the bottom left as they look further ahead. Nimbo and the still radar image stay on the bias 1 line: they forecast rain about as often as it falls. pySTEPS sits on the bias 2 line: it forecasts rain of 1 mm/h twice as often as it falls, so it catches most of the rain (detection 0.82 at 15 minutes) at the price of about three false alarms for every two hits. Its hit score alone would not tell you that.

Near misses: the neighbourhood score

Rain forecast 4 km east of where it fell counts twice against a forecast square by square: once as a miss, once as a false alarm. A forecast that blurs its rain into a wide grey area suffers less from this than a sharp one, which is not what users want. The fractions skill score (FSS) answers it: in a window around each square, compare the share of rainy squares forecast with the share observed.

FSS = 1 − Σ(F − O)² / (ΣF² + ΣO²), where F and O are the forecast and observed shares in each window. With a window of one square, it is exactly the F-score. A common rule of thumb calls a forecast useful at the scale where FSS passes 0.5 + f/2, with f the share of the map that is rainy.

In [2]:

def fss(fcst, obs, thr, size):
    f = uniform_filter((fcst >= thr) * 1.0, size)
    o = uniform_filter((obs >= thr) * 1.0, size)
    return 1 - ((f - o)**2).sum() / ((f**2).sum() + (o**2).sum())

Out[2]:

Read across the dashed line: one hour ahead, this forecast is useful at about 6 km; two hours ahead, at about 26 km; four hours ahead, not even at 50 km. The scale at which a forecast can be trusted grows with time ahead, which is a more honest summary than any single score. Nimbo’s scoreboard reports the near-miss score at 10 km.

Score against a truth the forecast never saw

There are three common truths, and they are not equal:

• Radar covers every square, but it is also the forecast’s input. A forecast scored on the radar it was built from is partly marking its own homework.
• Rain gauges measure the water that actually fell, hour by hour, at thousands of points. No radar nowcast uses them as input, which makes them the fairest check.
• Airport reports say whether it is raining every half hour or so at hundreds of airports. Coarse, but independent and worldwide.

Nimbo’s scoreboard uses all three: radar, 5,369 gauges in 10 national networks, and 757 airports.

Every day is different

A score pooled over many days hides how much it moves from one day to the next. Figure 3 takes the pooled comparison apart: one dot per rainy day, Nimbo’s hit score minus pySTEPS’s, one hour ahead, on each truth.

In [3]:

for truth in ["radar", "gauges"]:
    days = pd.DataFrame(bench["by_day"][truth])
    diff = days["product"] - days["sprog_raw"]
    print(truth, (diff > 0).sum(), "of", len(diff))

Out[3]:

Against radar the gap is steady. Against gauges, a far smaller sample of points, it swings from day to day, and on 9 of 39 days pySTEPS did better. This is why a difference between two forecasts needs an interval: resample the days many times, recompute the gap each time, and report the range it moves in. The scoreboard does this with whole blocks of neighbouring days, since days a day or two apart often share the same weather system.

The traps

• Days the model was trained or tuned on. A model scored on days it learned from will look better than it is. Only days set aside before training count.
• Different events for different forecasts. Two forecasts must be scored on exactly the same places and times, or the comparison measures the weather, not the forecasts.
• Dry days. Scores on days with little rain swing wildly. Choose rainy days, and say how they were chosen.
• One number. A forecast can be good one hour ahead and poor at three, good for drizzle and poor for downpours. Look at skill against time ahead and by rain threshold.
• No error bars. With a few dozen days, small differences are noise. Resample the days to see how much a score could move.
• Blurring to win. Smoothing a forecast can raise some scores while making the map useless. A gain should hold on independent truths, not only on the score it was tuned for.

A worked example

On the scoreboard’s rainy days, Nimbo’s hit score for a wet hour (1 mm or more) in rain gauges is about 0.47 one hour ahead, 0.38 two hours ahead and 0.31 three hours ahead. pySTEPS, scored on the same gauges and hours, reaches about 0.40, 0.31 and 0.25. Both fall with time; the gap between them is the measure that matters, and it comes with an interval from resampling the days.

Next: Radar coverage in Europe