Circular metrics for directional quantities everywhere; raise on linear ones

src/modelskill/comparison/_collection.py

@@ -418,6 +418,16 @@ def query(self, query: str) -> "ComparerCollection":
 
         return ComparerCollection(cmps_with_data)
 
+    def _is_directional(self) -> bool:
+        directional = {cmp.quantity.is_directional for cmp in self}
+        if len(directional) > 1:
+            raise ValueError(
+                "The collection mixes directional and non-directional quantities, "
+                "which need different metrics. Select one quantity first, "
+                "e.g. cc.sel(quantity=...)."
+            )
+        return directional.pop() if directional else False
+
     def skill(
         self,
         by: str | Iterable[str] | None = None,
@@ -482,7 +492,7 @@ def skill(
         """
         cc = self
 
-        pmetrics = _parse_metric(metrics)
+        pmetrics = _parse_metric(metrics, directional=self._is_directional())
 
         agg_cols = _parse_groupby(by, n_mod=cc.n_models, n_qnt=cc.n_quantities)
         agg_cols, attrs_keys = self._attrs_keys_in_by(agg_cols)
@@ -634,7 +644,7 @@ def gridded_skill(
         """
         cmp = self
 
-        metrics = _parse_metric(metrics)
+        metrics = _parse_metric(metrics, directional=self._is_directional())
 
         df = cmp._to_long_dataframe()
         df = _add_spatial_grid_to_df(df=df, bins=bins, binsize=binsize)
@@ -724,7 +734,7 @@ def mean_skill(
             case _:
                 raise ValueError("Invalid weights specification")
 
-        pmetrics = _parse_metric(metrics)
+        pmetrics = _parse_metric(metrics, directional=self._is_directional())
         skilldf = (
             self.skill(metrics=pmetrics)
             .to_dataframe()
@@ -761,7 +771,7 @@ def _mean_skill_by(self, skilldf: pd.DataFrame) -> list[str]:
 
     def score(
         self,
-        metric: str | Callable = mtr.rmse,
+        metric: str | Callable | None = None,
         weights: str | list[float] | dict[str, float] | None = None,
     ) -> Dict[str, float]:
         """Weighted mean score of model(s) over all observations
@@ -779,8 +789,9 @@ def score(
             - "points": giving all points equal weight,
             - list of weights e.g. [0.3, 0.3, 0.4] per observation,
             - dictionary of observations with special weigths, others will be set to 1.0
-        metric : list, optional
+        metric : str or callable, optional
             a single metric from modelskill.metrics, by default rmse
+            (c_rmse for a directional quantity)
 
         Returns
         -------
@@ -809,7 +820,10 @@ def score(
         {'mod': 8.414442957854142}
         """
 
-        metric = _parse_metric(metric)[0]
+        directional = self._is_directional()
+        if metric is None:
+            metric = mtr.c_rmse if directional else mtr.rmse
+        metric = _parse_metric(metric, directional=directional)[0]
 
         if not (callable(metric) or isinstance(metric, str)):
             raise ValueError("metric must be a string or a function")

src/modelskill/comparison/_collection_plotter.py

@@ -500,6 +500,12 @@ def taylor(
         Copin, Y. (2018). https://gist.github.com/ycopin/3342888, Yannick Copin <yannick.copin@laposte.net>
         """
 
+        if self.cc._is_directional():
+            raise ValueError(
+                "A Taylor diagram uses linear standard deviation and correlation, "
+                "which are not defined for a directional quantity."
+            )
+
         if (not aggregate_observations) and (not normalize_std):
             raise ValueError(
                 "aggregate_observations=False is only possible if normalize_std=True!"

src/modelskill/comparison/_comparer_plotter.py

@@ -723,6 +723,11 @@ def taylor(
         -----
         Copin, Y. (2018). https://gist.github.com/ycopin/3342888, Yannick Copin <yannick.copin@laposte.net>
         """
+        if self.is_directional:
+            raise ValueError(
+                "A Taylor diagram uses linear standard deviation and correlation, "
+                "which are not defined for a directional quantity."
+            )
         cmp = self.comparer
 
         # TODO consider if this round-trip  via mtr is necessary to get the std:s

src/modelskill/comparison/_comparison.py

@@ -1173,14 +1173,15 @@ def _add_as_col_if_not_in_index(
 
     def score(
         self,
-        metric: str | Callable = mtr.rmse,
+        metric: str | Callable | None = None,
     ) -> Dict[str, float]:
         """Model skill score
 
         Parameters
         ----------
-        metric : list, optional
+        metric : str or callable, optional
             a single metric from modelskill.metrics, by default rmse
+            (c_rmse for a directional quantity)
 
         Returns
         -------
@@ -1202,7 +1203,10 @@ def score(
         >>> cmp.score(metric="mape")
         {'mod': 11.567399646108198}
         """
-        metric = _parse_metric(metric)[0]
+        directional = self.quantity.is_directional
+        if metric is None:
+            metric = mtr.c_rmse if directional else mtr.rmse
+        metric = _parse_metric(metric, directional=directional)[0]
         if not (callable(metric) or isinstance(metric, str)):
             raise ValueError("metric must be a string or a function")
 
@@ -1283,7 +1287,7 @@ def gridded_skill(
 
         cmp = self
 
-        metrics = _parse_metric(metrics)
+        metrics = _parse_metric(metrics, directional=self.quantity.is_directional)
         if cmp.n_points == 0:
             raise ValueError("No data to compare")
 

src/modelskill/metrics.py

@@ -1201,9 +1201,30 @@ def _parse_metric(
         else:
             raise TypeError(f"metric {m} must be a string or callable")
 
+    if directional:
+        for m in parsed_metrics:
+            _check_circular(m)
+
     return parsed_metrics
 
 
+def _check_circular(metric: Callable) -> None:
+    # A user-defined metric may well handle angles; only built-ins are known not to.
+    name = metric.__name__
+    if metric.__module__ != __name__ or name.startswith("c_"):
+        return
+    circular = f"c_{name}"
+    hint = (
+        f"Use '{circular}' instead."
+        if hasattr(sys.modules[__name__], circular)
+        else f"Use one of the circular metrics: {', '.join(m.__name__ for m in default_circular_metrics)}."
+    )
+    raise ValueError(
+        f"'{name}' treats values as linear, but the quantity is directional, "
+        f"so 359° vs 1° would count as 358° apart. {hint}"
+    )
+
+
 def is_best(metric: str, expected: str | int) -> bool:
     try:
         func = get_metric(metric)

tests/test_directional_skill.py

@@ -0,0 +1,97 @@
+import pandas as pd
+import pytest
+
+import modelskill as ms
+
+DIRECTION = ms.Quantity("Wave direction", "degree", is_directional=True)
+HEIGHT = ms.Quantity("Wave height", "m")
+
+
+def _cmp(quantity: ms.Quantity, name: str) -> ms.Comparer:
+    # 1° observed vs 359° modelled is a 2° error, not 358°
+    df = pd.DataFrame(
+        {"obs": [1.0, 2.0, 3.0], "mod": [359.0, 0.0, 1.0]},
+        index=pd.date_range("2020", periods=3, freq="h"),
+    )
+    return ms.from_matched(df, obs_item="obs", quantity=quantity, name=name)
+
+
+@pytest.fixture
+def cmp() -> ms.Comparer:
+    return _cmp(DIRECTION, "buoy")
+
+
+@pytest.fixture
+def cc(cmp) -> ms.ComparerCollection:
+    return ms.ComparerCollection([cmp])
+
+
+def test_comparer_skill_default_is_circular(cmp):
+    assert cmp.skill().to_dataframe().loc["buoy", "c_bias"] == pytest.approx(-2.0)
+
+
+def test_comparer_score_default_is_circular(cmp):
+    assert cmp.score()["mod"] == pytest.approx(2.0)
+
+
+def test_collection_skill_default_is_circular(cc):
+    assert cc.skill().to_dataframe().loc["buoy", "c_bias"] == pytest.approx(-2.0)
+
+
+def test_collection_score_default_is_circular(cc):
+    assert cc.score()["mod"] == pytest.approx(2.0)
+
+
+def test_collection_mean_skill_default_is_circular(cc):
+    assert cc.mean_skill().to_dataframe().loc["mod", "c_rmse"] == pytest.approx(2.0)
+
+
+@pytest.mark.parametrize("metric", ["bias", "rmse", "mae", "cc"])
+def test_comparer_skill_with_linear_metric_raises(cmp, metric):
+    with pytest.raises(ValueError, match="directional"):
+        cmp.skill(metrics=[metric])
+
+
+def test_error_names_the_circular_metric(cmp):
+    with pytest.raises(ValueError, match="c_rmse"):
+        cmp.skill(metrics=["rmse"])
+
+
+def test_collection_skill_with_linear_metric_raises(cc):
+    with pytest.raises(ValueError, match="directional"):
+        cc.skill(metrics=["rmse"])
+
+
+def test_score_with_linear_metric_raises(cmp, cc):
+    with pytest.raises(ValueError, match="directional"):
+        cmp.score(metric="rmse")
+    with pytest.raises(ValueError, match="directional"):
+        cc.score(metric="rmse")
+
+
+def test_custom_metric_is_passed_through(cmp):
+    def my_metric(obs, model):
+        return 42.0
+
+    assert cmp.skill(metrics=[my_metric]).to_dataframe().loc[
+        "buoy", "my_metric"
+    ] == pytest.approx(42.0)
+
+
+def test_collection_mixing_directional_and_scalar_raises():
+    cc = ms.ComparerCollection([_cmp(DIRECTION, "dir"), _cmp(HEIGHT, "hm0")])
+    with pytest.raises(ValueError, match="quantity"):
+        cc.skill()
+
+
+def test_collection_mixing_directional_and_scalar_works_per_quantity():
+    cc = ms.ComparerCollection([_cmp(DIRECTION, "dir"), _cmp(HEIGHT, "hm0")])
+    sk = cc.sel(quantity="Wave direction").skill()
+    assert sk.to_dataframe().loc["dir", "c_bias"] == pytest.approx(-2.0)
+
+
+def test_taylor_diagram_raises_for_directional(cmp, cc):
+    with pytest.raises(ValueError, match="Taylor diagram"):
+        cmp.plot.taylor()
+    with pytest.raises(ValueError, match="Taylor diagram"):
+        cc.plot.taylor()