Skip to content

Commit 8684990

Browse files
authored
Merge pull request #671 from DHI/docs/metrics-user-guide
Add metrics page to user guide
2 parents d1916ee + 9eed954 commit 8684990

4 files changed

Lines changed: 230 additions & 6 deletions

File tree

docs/_quarto.yml

Lines changed: 2 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -48,6 +48,7 @@ website:
4848
- user-guide/selecting-data.qmd
4949
- user-guide/plotting.qmd
5050
- user-guide/statistics.qmd
51+
- user-guide/metrics.qmd
5152
- section: "Extensions"
5253
contents:
5354
- user-guide/network.qmd
@@ -286,6 +287,7 @@ quartodoc:
286287
- c_mae
287288
- c_mean_absolute_error
288289
- mape
290+
- mean_absolute_percentage_error
289291
- cc
290292
- corrcoef
291293
- rho

docs/user-guide/index.qmd

Lines changed: 1 addition & 1 deletion
Original file line numberDiff line numberDiff line change
@@ -8,7 +8,7 @@ format-links: false
88
ModelSkill compares model results with observations. The workflow can be split in two phases:
99

1010
1. [Matching](matching.qmd) - making sure that observations and model results are in the same space and time
11-
2. Analysis - [plots](plotting.qmd) and [statistics](statistics.qmd) of the matched data
11+
2. Analysis - [plots](plotting.qmd) and [statistics](statistics.qmd) of the matched data (see the [metrics guide](metrics.qmd) for choosing skill metrics)
1212

1313
If the observations and model results are already matched (i.e. are stored in the same data source),
1414
the `from_matched()` function can be used to go directly to the analysis phase.

docs/user-guide/metrics.qmd

Lines changed: 220 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,220 @@
1+
# Metrics
2+
3+
ModelSkill comes with a comprehensive set of skill metrics. This page is a *when-to-use*
4+
guide: how the metrics differ, which one answers which question, and the common pitfalls.
5+
Each metric name links to its full mathematical definition in the
6+
[metrics API reference](../api/metrics.qmd).
7+
8+
Metrics are passed to [](`~modelskill.Comparer.skill`) (and
9+
[](`~modelskill.ComparerCollection.mean_skill`), [](`~modelskill.Comparer.score`)) as
10+
lower-case strings, or as callables:
11+
12+
```{python}
13+
#| code-fold: true
14+
#| code-summary: "Construct a comparer from observation and model data"
15+
import modelskill as ms
16+
o1 = ms.observation("../data/SW/HKNA_Hm0.dfs0", item=0,
17+
x=4.2420, y=52.6887, name="HKNA")
18+
o2 = ms.observation("../data/SW/eur_Hm0.dfs0", item=0,
19+
x=3.2760, y=51.9990, name="EPL")
20+
mr = ms.model_result("../data/SW/HKZN_local_2017_DutchCoast.dfsu",
21+
item="Sign. Wave Height", name="m1")
22+
cc = ms.match([o1, o2], mr)
23+
```
24+
25+
```{python}
26+
cc.skill(metrics=["bias", "rmse", "kge"])
27+
```
28+
29+
In the tables below the first column gives the metric string ModelSkill accepts; the
30+
**alias** (a longer synonym that also works) is listed underneath it. Both link to the API.
31+
32+
When you pass no `metrics=` argument, [](`~modelskill.Comparer.skill`) reports the default
33+
set — the point count `n` plus `bias`, `rmse`, `urmse`, `mae`, `cc`, `si`, `r2`. You can
34+
change the default via `ms.options.metrics.list = [...]`.
35+
36+
```{=html}
37+
<style>
38+
.metrics-ref table { table-layout: fixed; width: 100%; }
39+
.metrics-ref table th:nth-child(1), .metrics-ref table td:nth-child(1) { width: 22%; }
40+
.metrics-ref table th:nth-child(2), .metrics-ref table td:nth-child(2) { width: 12%; }
41+
.metrics-ref table th:nth-child(3), .metrics-ref table td:nth-child(3) { width: 6%; }
42+
.metrics-ref table td { vertical-align: top; }
43+
.metrics-ref table code { overflow-wrap: anywhere; }
44+
</style>
45+
```
46+
47+
---
48+
49+
## Error magnitude — same units as the data (scale-dependent)
50+
51+
These are in the data's units (m, m³/s, …), so they are **not** comparable across stations
52+
of different magnitude — a 5,000 m³/s river will always out-RMSE a 50 m³/s tributary.
53+
54+
::: {.metrics-ref}
55+
| Metric | Range | Best | What it is |
56+
|---|---|---|---|
57+
| [](`~modelskill.metrics.bias`) | (−∞, ∞) | 0 | Mean error, `mean(model − obs)`. **Sign matters**: + = model runs high, − = runs low. |
58+
| [](`~modelskill.metrics.mae`)<br>[](`~modelskill.metrics.mean_absolute_error`) | [0, ∞) | 0 | Mean absolute error. Typical error magnitude, **robust to outliers**. |
59+
| [](`~modelskill.metrics.rmse`)<br>[](`~modelskill.metrics.root_mean_squared_error`) | [0, ∞) | 0 | Root-mean-square error. Like MAE but **penalises large misses** disproportionately. |
60+
| [](`~modelskill.metrics.urmse`) | [0, ∞) | 0 | Unbiased RMSE — the **random/scatter** part after removing bias. `rmse² = bias² + urmse²`. |
61+
| [](`~modelskill.metrics.max_error`) | [0, ∞) | 0 | Largest single absolute error — the **worst-case** miss. |
62+
:::
63+
64+
**Use when:** you want the error in physical units. `bias` for systematic offset, `rmse` as
65+
the default, `mae` when a few outliers shouldn't dominate, `urmse` to isolate random error
66+
from bias, `max_error` for safety-/threshold-critical checks.
67+
68+
---
69+
70+
## Relative / normalised error — dimensionless
71+
72+
Divide the error by a scale, so they *can* be compared across stations — but they divide by
73+
observed values, so they are **unstable when observations pass through zero** (watch out for
74+
water level around a datum).
75+
76+
::: {.metrics-ref}
77+
| Metric | Range | Best | What it is |
78+
|---|---|---|---|
79+
| [](`~modelskill.metrics.si`)<br>[](`~modelskill.metrics.scatter_index`) | [0, ∞) | 0 | Scatter index = `urmse / mean(obs)`. Random error as a fraction of the mean. |
80+
| [](`~modelskill.metrics.mape`)<br>[](`~modelskill.metrics.mean_absolute_percentage_error`) | [0, ∞) | 0 | Mean absolute % error, `mean(\|model−obs\| / \|obs\|)`. Intuitive %, but blows up near obs = 0. |
81+
:::
82+
83+
**Use when:** comparing error across stations of very different magnitude, **and** the
84+
observed values stay comfortably away from zero. For water level, prefer the efficiency
85+
scores below over `mape`/`si`.
86+
87+
---
88+
89+
## Efficiency / skill scores — dimensionless, comparable across stations
90+
91+
The "one number for overall skill" family. All reference the observed mean or variance, so a
92+
score of 0 (for NSE/KGE) means "no better than predicting the mean."
93+
94+
::: {.metrics-ref}
95+
| Metric | Range | Best | What it is |
96+
|---|---|---|---|
97+
| [](`~modelskill.metrics.nse`)<br>[](`~modelskill.metrics.nash_sutcliffe_efficiency`) | (−∞, 1] | 1 | Nash–Sutcliffe efficiency. 0 = no better than the obs mean; < 0 = worse. |
98+
| [](`~modelskill.metrics.r2`) | (−∞, 1] | 1 | Coefficient of determination. Identical to NSE (see below). |
99+
| [](`~modelskill.metrics.kge`)<br>[](`~modelskill.metrics.kling_gupta_efficiency`) | (−∞, 1] | 1 | Kling–Gupta — composite of **correlation + bias ratio + variability ratio**. |
100+
| [](`~modelskill.metrics.willmott`) | [0, 1] | 1 | Willmott's Index of Agreement — bounded [0,1], less harsh on bias than NSE. |
101+
| [](`~modelskill.metrics.ev`)<br>[](`~modelskill.metrics.explained_variance`) | (−∞, 1] | 1 | Proportion of variance explained (differs from NSE when the model is biased). |
102+
| [](`~modelskill.metrics.mef`)<br>[](`~modelskill.metrics.model_efficiency_factor`) | [0, ∞) | 0 | `RMSE / std(obs) = √(1 − NSE)`. Same information as NSE, expressed as an **error** (lower is better). |
103+
:::
104+
105+
**Use when:** you need a single dimensionless score to rank models or compare stations.
106+
Reach for **`kge`** when you want to *diagnose why* a model fails (it separates correlation,
107+
bias, and variance); **`nse`** is the hydrology standard for overall predictive power;
108+
`willmott` if you want a strictly bounded [0,1] score.
109+
110+
### `cc`, `ev` and `r2`: a nested hierarchy
111+
112+
A common question is how `cc`, `ev` and `r2` relate. They answer increasingly strict
113+
versions of *"how much of the observed variation does the model capture?"*, each penalising
114+
one more kind of error:
115+
116+
| Score | Penalises bias (offset)? | Penalises wrong amplitude? | Scored against |
117+
|---|---|---|---|
118+
| `cc²` (square of `cc`) | no | no | the best-fit line (free slope + intercept) |
119+
| `ev` | no | yes | obs variance, ignoring a constant offset |
120+
| `r2` (= `nse`) | yes | yes | the 1:1 line |
121+
122+
For the same data this gives the ordering **`cc²``ev``r2`**: `cc²` forgives both a
123+
constant offset and a wrong amplitude, `ev` forgives only the offset, and `r2`/`nse` forgives
124+
nothing (it scores against the 1:1 line). The gap `ev − r2` is exactly the squared,
125+
normalised bias.
126+
127+
In practice, report `cc` (timing/phase) **plus** one efficiency score (`nse` or `kge`)
128+
**plus** `bias` separately — that trio localises a failure to phase, offset, or amplitude,
129+
which no single number can. `cc²` rarely adds anything once you already report `cc`, and `ev`
130+
(scikit-learn's `explained_variance_score`) sits between the two and is seldom reported on its
131+
own in water modelling.
132+
133+
::: {.callout-warning}
134+
## `r2` is the coefficient of determination, not squared correlation
135+
ModelSkill's `r2` equals **NSE** — they are the same number under two names — *not* squared
136+
Pearson correlation. The two generally differ (`r2` penalises bias, `cc²` does not). If you
137+
want squared correlation, compute `cc` and square it yourself.
138+
:::
139+
140+
---
141+
142+
## Correlation & amplitude — dimensionless
143+
144+
These ignore systematic bias, so always read them **alongside** `bias`.
145+
146+
::: {.metrics-ref}
147+
| Metric | Range | Best | What it is |
148+
|---|---|---|---|
149+
| [](`~modelskill.metrics.cc`)<br>[](`~modelskill.metrics.corrcoef`) | [−1, 1] | 1 | Pearson correlation — **linear co-variation / timing**. Blind to bias and amplitude. |
150+
| [](`~modelskill.metrics.rho`)<br>[](`~modelskill.metrics.spearmanr`) | [−1, 1] | 1 | Spearman rank correlation — monotonic, **robust** to outliers and non-linearity. |
151+
| [](`~modelskill.metrics.lin_slope`) | (−∞, ∞) | 1 | Slope of the model-vs-obs regression. < 1 = model **under-responds** in amplitude. |
152+
:::
153+
154+
**Use when:** the question is about **phase/timing** (`cc`), a monotonic-but-non-linear
155+
relationship (`rho`), or **amplitude** of the response (`lin_slope`).
156+
157+
---
158+
159+
## Event & distribution
160+
161+
::: {.metrics-ref}
162+
| Metric | Range | Best | What it is |
163+
|---|---|---|---|
164+
| [](`~modelskill.metrics.peak_ratio`)<br>[](`~modelskill.metrics.pr`) | [0, ∞) | 1 | Ratio of modelled to observed **peaks** (mean over matched peak events). < 1 = peaks under-predicted. |
165+
| [](`~modelskill.metrics.hit_ratio`) | [0, 1] | 1 | Fraction of points within an acceptable deviation `a` of the observation. Takes an `a=` argument. |
166+
:::
167+
168+
**Use when:** `peak_ratio` for storm-surge / flood-peak capture. `hit_ratio` for
169+
acceptance-criterion reporting — "X % of points within ±0.1 m" (set the tolerance with `a=`).
170+
171+
---
172+
173+
## Directional / circular — for direction variables only
174+
175+
For wind / wave / current **direction**, where 359° and 1° are 2° apart, not 358°. These
176+
require the quantity to be flagged directional
177+
(`Quantity(..., is_directional=True)`); they handle the 0–360° wrap-around. See the
178+
[directional data example](../examples/Directional_data_comparison.qmd).
179+
180+
::: {.metrics-ref}
181+
| Metric | Range | Best | What it is |
182+
|---|---|---|---|
183+
| [](`~modelskill.metrics.c_bias`) | [−180, 180] | 0 | Circular bias (mean angular error). |
184+
| [](`~modelskill.metrics.c_mae`)<br>[](`~modelskill.metrics.c_mean_absolute_error`) | [0, 180] | 0 | Circular mean absolute error. |
185+
| [](`~modelskill.metrics.c_rmse`)<br>[](`~modelskill.metrics.c_root_mean_squared_error`) | [0, 180] | 0 | Circular RMSE. |
186+
| [](`~modelskill.metrics.c_urmse`)<br>[](`~modelskill.metrics.c_unbiased_root_mean_squared_error`) | [0, 180] | 0 | Circular unbiased RMSE. |
187+
| [](`~modelskill.metrics.c_max_error`) | [0, 180] | 0 | Largest circular (angular) error. |
188+
:::
189+
190+
**Use when:** the variable is an angle. Don't use the scalar metrics above on direction data.
191+
192+
---
193+
194+
## Quick decision guide
195+
196+
| Your question | Reach for |
197+
|---|---|
198+
| Is the model systematically high or low? | `bias` |
199+
| Does it under-/over-respond in amplitude? | `lin_slope` |
200+
| How big is a typical error, in my units? | `rmse` (penalise big misses) or `mae` (robust) |
201+
| Split random vs systematic error? | `urmse` vs `bias` (`rmse² = bias² + urmse²`) |
202+
| What's the worst single miss? | `max_error` |
203+
| One dimensionless score, comparable across stations? | `nse` or `kge` |
204+
| *Why* does the model fail (corr / bias / variance)? | `kge` |
205+
| Compare error across very different-sized stations? | `si` or `nse`/`kge` (avoid `rmse`) |
206+
| Timing / phase agreement? | `cc` |
207+
| Do we capture the storm peaks? | `peak_ratio` |
208+
| What fraction meets an acceptance tolerance? | `hit_ratio` (set `a=`) |
209+
| Working with direction (wind/wave/current)? | the `c_*` family |
210+
211+
## Cross-cutting reminders
212+
213+
- **Scale-dependent vs dimensionless.** The error-magnitude metrics are in data units — never
214+
compare them across stations of different magnitude. The rest are dimensionless and comparable.
215+
- **Always pair correlation with bias.** `cc`/`rho` are blind to offset; a model can correlate
216+
perfectly while sitting 1 m too high.
217+
- **Division-by-obs metrics** (`mape`, `si`) are fragile near zero observations.
218+
- **Custom metrics**: any `f(obs, model) -> float` callable drops straight into
219+
`metrics=[...]`; the column takes the function name. See the
220+
[custom metric example](../examples/Metrics_custom_metric.qmd).

src/modelskill/metrics.py

Lines changed: 7 additions & 5 deletions
Original file line numberDiff line numberDiff line change
@@ -483,7 +483,7 @@ def scatter_index2(obs: ArrayLike, model: ArrayLike) -> Any:
483483
{\sum_{i=1}^n obs_i^2}}
484484
$$
485485
486-
Range: [0, 100]; Best: 0
486+
Range: $[0, \infty)$; Best: 0
487487
"""
488488
assert obs.size == model.size
489489
if len(obs) == 0:
@@ -506,8 +506,10 @@ def explained_variance(obs: ArrayLike, model: ArrayLike) -> Any:
506506
r"""EV: Explained variance
507507
508508
EV is the explained variance and measures the proportion
509-
[0 - 1] to which the model accounts for the variation
510-
(dispersion) of the observations.
509+
to which the model accounts for the variation (dispersion)
510+
of the observations. A perfect match gives 1; a model that
511+
explains less variance than the observed mean gives a
512+
negative value.
511513
512514
In cases with no bias, EV is equal to r2
513515
@@ -518,7 +520,7 @@ def explained_variance(obs: ArrayLike, model: ArrayLike) -> Any:
518520
(obs_i - \overline{obs})^2}
519521
$$
520522
521-
Range: [0, 1]; Best: 1
523+
Range: $(-\infty, 1]$; Best: 1
522524
523525
See Also
524526
--------
@@ -958,7 +960,7 @@ def c_max_error(obs: ArrayLike, model: ArrayLike) -> Any:
958960
959961
Notes
960962
-----
961-
Range: $[0, \\infty)$; Best: 0
963+
Range: $[0, 180]$; Best: 0
962964
963965
Returns
964966
-------

0 commit comments

Comments
 (0)