Merge branch 'main' of https://github.com/DHI/fmskill into main

Author: jsmariegaard

README.md

@@ -39,7 +39,7 @@ Or the development version:
 * [NetCDF_ModelResult.ipynb](https://nbviewer.jupyter.org/github/DHI/fmskill/blob/main/notebooks/NetCDF_ModelResult.ipynb)
 * [Combine_comparers.ipynb](https://nbviewer.jupyter.org/github/DHI/fmskill/blob/main/notebooks/Combine_comparers.ipynb)
 * [DMI_observations.ipynb](https://nbviewer.jupyter.org/github/DHI/fmskill/blob/main/notebooks/DMI_observations.ipynb) (download data from REST api)
-
+* [Altimetry_data.ipynb](https://nbviewer.jupyter.org/github/DHI/fmskill/blob/main/notebooks/Altimetry_data.ipynb) (download data from altimetry api)
 
 ## Workflow
 

docs/api.rst

@@ -111,6 +111,7 @@ Metrics
 	fmskill.metrics.rho
 	fmskill.metrics.r2
 	fmskill.metrics.lin_slope
+	fmskill.metrics.willmott
 	fmskill.metrics.hit_ratio
 
 .. automodule:: fmskill.metrics

fmskill/metrics.py

@@ -45,6 +45,8 @@
 0.637783218973691
 >>> lin_slope(obs, mod)
 0.4724896836313617
+>>> willmott(obs, mod)
+0.7484604452865941
 >>> hit_ratio(obs, mod, a=0.5)
 0.6666666666666666
 """
@@ -59,13 +61,26 @@ def bias(obs, model) -> float:
     .. math::
         bias=\\frac{1}{n}\\sum_{i=1}^n (model_i - obs_i)
 
-    Range: -infinity to infinity; Best: 0.0
+    Range: :math:`(-\\infty, \\infty)`; Best: 0
     """
 
     assert obs.size == model.size
     return np.mean(model.ravel() - obs.ravel())
 
 
+def max_error(obs, model) -> float:
+    """Max (absolute) error
+
+    .. math::
+        max_error = max(|model_i - obs_i|)
+
+    Range: :math:`[0, \\infty)`; Best: 0
+    """
+
+    assert obs.size == model.size
+    return np.max(np.abs(model.ravel() - obs.ravel()))
+
+
 def mae(obs: np.ndarray, model: np.ndarray, weights: np.ndarray = None) -> float:
     """alias for mean_absolute_error"""
     assert obs.size == model.size
@@ -80,7 +95,7 @@ def mean_absolute_error(
     .. math::
         MAE=\\frac{1}{n}\\sum_{i=1}^n|model_i - obs_i|
 
-    Range: 0.0 to infinity; Best: 0.0
+    Range: :math:`[0, \\infty)`; Best: 0
     """
     assert obs.size == model.size
 
@@ -100,7 +115,7 @@ def mean_absolute_percentage_error(obs: np.ndarray, model: np.ndarray) -> float:
     .. math::
         MAPE=\\frac{1}{n}\\sum_{i=1}^n\\frac{|model_i - obs_i|}{obs_i}*100
 
-    Range: 0.0 to infinity; Best: 0.0
+    Range: :math:`[0, \\infty)`; Best: 0
     """
 
     assert obs.size == model.size
@@ -125,7 +140,7 @@ def urmse(obs: np.ndarray, model: np.ndarray, weights: np.ndarray = None) -> flo
 
         uRMSE = \\sqrt{\\frac{1}{n} \\sum_{i=1}^n res_{u,i}^2}
 
-    Range: 0.0 to infinity; Best: 0.0
+    Range: :math:`[0, \\infty)`; Best: 0
 
     See Also
     --------
@@ -165,7 +180,7 @@ def root_mean_squared_error(
 
         uRMSE=\\sqrt{\\frac{1}{n} \\sum_{i=1}^n res_{u,i}^2}
 
-    Range: 0.0 to infinity; Best: 0.0
+    Range: :math:`[0, \\infty)`; Best: 0
 
     """
     assert obs.size == model.size
@@ -191,7 +206,7 @@ def nash_sutcliffe_efficiency(obs: np.ndarray, model: np.ndarray) -> float:
         NSE = 1 - \\frac {\\sum _{i=1}^{n}\\left(model_{i} - obs_{i}\\right)^{2}}
                        {\\sum_{i=1}^{n}\\left(obs_{i} - {\\overline{obs}}\\right)^{2}}
 
-    Range: -infinity to 1.0; Best: 1.0
+    Range: :math:`(-\\infty, 1]`; Best: 1
 
     Note
     ----
@@ -214,16 +229,16 @@ def nash_sutcliffe_efficiency(obs: np.ndarray, model: np.ndarray) -> float:
 
 
 def r2(obs: np.ndarray, model: np.ndarray) -> float:
-    """Coefficient of determination (R2) - pronounced 'R-squared'
+    """Coefficient of determination (R2)
 
-    The proportion of the variation in the dependent variable that is predictable from the independent variable(s), e.g. the proportion of explained variance.
+    Pronounced 'R-squared'; the proportion of the variation in the dependent variable that is predictable from the independent variable(s), i.e. the proportion of explained variance.
 
     .. math::
 
         R^2 = 1 - \\frac{\\sum_{i=1}^n (model_i - obs_i)^2}
                     {\\sum_{i=1}^n (obs_i - \\overline {obs})^2}
 
-    Range: -infinity to 1.0; Best: 1.0
+    Range: :math:`(-\\infty, 1]`; Best: 1
 
     Note
     ----
@@ -262,7 +277,7 @@ def model_efficiency_factor(obs: np.ndarray, model: np.ndarray) -> float:
         MEF = \\frac{RMSE}{STDEV}=\\frac{\\sqrt{\\frac{1}{n} \\sum_{i=1}^n(model_i - obs_i)^2}}
                                         {\\sqrt{\\frac{1}{n} \\sum_{i=1}^n(obs_i - \\overline{obs})^2}}=\\sqrt{1-NSE}
 
-    Range: 0.0 to infinity; Best: 0.0
+    Range: :math:`[0, \\infty)`; Best: 0
 
     See Also
     --------
@@ -288,11 +303,11 @@ def corrcoef(obs, model, weights=None) -> float:
                    {\\sqrt{\\sum_{i=1}^n (model_i - \\overline{model})^2}
                     \\sqrt{\\sum_{i=1}^n (obs_i - \\overline{obs})^2} }
 
-    Range: -1.0 to 1.0; Best: 1.0
+    Range: [-1, 1]; Best: 1
 
     See Also
     --------
-    numpy.corrcoef
+    np.corrcoef
     """
     assert obs.size == model.size
     if len(obs) <= 1:
@@ -321,7 +336,7 @@ def spearmanr(obs: np.ndarray, model: np.ndarray) -> float:
                       {\\sqrt{\\sum_{i=1}^n (rmodel_i - \\overline{rmodel})^2}
                        \\sqrt{\\sum_{i=1}^n (robs_i - \\overline{robs})^2} }
 
-    Range: -1.0 to 1.0; Best: 1.0
+    Range: [-1, 1]; Best: 1
 
     Examples
     --------
@@ -353,7 +368,7 @@ def scatter_index(obs: np.ndarray, model: np.ndarray) -> float:
         \\sqrt {\\frac{\\sum_{i=1}^n \\left( (model_i - \\overline {model}) - (obs_i - \\overline {obs}) \\right)^2}
         {\\sum_{i=1}^n obs_i^2}}
 
-    Range: 0.0 to 100.0; Best: 0.0
+    Range: [0, 100]; Best: 0
     """
     assert obs.size == model.size
     if len(obs) == 0:
@@ -365,19 +380,54 @@ def scatter_index(obs: np.ndarray, model: np.ndarray) -> float:
     )
 
 
+def willmott(obs: np.ndarray, model: np.ndarray) -> float:
+    """willmott's Index of Agreement
+
+    A scaled representation of the predictive accuracy of the model against observations. A value of 1 indicates a perfect match, and 0 indicates no agreement at all.
+
+    .. math::
+
+        willmott = 1 - \\frac{\\frac{1}{n} \\sum_{i=1}^n(model_i - obs_i)^2}
+                           {\\frac{1}{n} \\sum_{i=1}^n(|model_i - \\overline{obs}| + |obs_i - \\overline{obs}|)^2}
+
+    Range: [0, 1]; Best: 1
+
+    Examples
+    --------
+    >>> obs = np.array([1.0, 1.1, 1.2, 1.3, 1.4, 1.4, 1.3])
+    >>> model = np.array([1.02, 1.16, 1.3, 1.38, 1.49, 1.45, 1.32])
+    >>> willmott(obs, model)
+    0.9501403174479723
+
+    References
+    ----------
+    Willmott, C. J. 1981. "On the validation of models". Physical Geography, 2, 184–194.
+    """
+
+    assert obs.size == model.size
+    if len(obs) == 0:
+        return np.nan
+
+    residual = model.ravel() - obs.ravel()
+    nominator = np.sum(residual ** 2)
+    denominator = np.sum((np.abs(model - obs.mean()) + np.abs(obs - obs.mean())) ** 2)
+
+    return 1 - nominator / denominator
+
+
 def hit_ratio(obs: np.ndarray, model: np.ndarray, a=0.1) -> float:
     """Fraction within obs ± acceptable deviation
 
     .. math::
 
         HR = \\frac{1}{n}\\sum_{i=1}^n I_{|(model_i - obs_i)|} < a
 
-    Range: 0.0 to 1.0; Best: 1.0
+    Range: [0, 1]; Best: 1
 
     Examples
     --------
-    >>> obs = np.array([1.0,1.1,1.2,1.3,1.4, 1.4, 1.3])
-    >>> model = np.array([1.02, 1.16, 1.3 , 1.38, 1.49, 1.45, 1.32])
+    >>> obs = np.array([1.0, 1.1, 1.2, 1.3, 1.4, 1.4, 1.3])
+    >>> model = np.array([1.02, 1.16, 1.3, 1.38, 1.49, 1.45, 1.32])
     >>> hit_ratio(obs, model, a=0.05)
     0.2857142857142857
     >>> hit_ratio(obs, model, a=0.1)
@@ -397,7 +447,7 @@ def lin_slope(obs: np.ndarray, model: np.ndarray, reg_method="ols") -> float:
         slope = \\frac{\\sum_{i=1}^n (model_i - \\overline {model})(obs_i - \\overline {obs})}
                       {\\sum_{i=1}^n (obs_i - \\overline {obs})^2}
 
-    Range: -infinity to infinity; Best: 1.0
+    Range: :math:`(-\\infty, \\infty )`; Best: 1
     """
     return _linear_regression(obs, model, reg_method)[0]
 

notebooks/Metrics_widget.ipynb

@@ -3,15 +3,15 @@
   {
    "cell_type": "markdown",
    "source": [
-    "# Metrics widget\r\n",
-    "\r\n",
+    "# Metrics widget\n",
+    "\n",
     "Execute this notebook to test out the different metrics. "
    ],
    "metadata": {}
   },
   {
    "cell_type": "code",
-   "execution_count": 1,
+   "execution_count": null,
    "source": [
     "import numpy as np\n",
     "import matplotlib.pyplot as plt\n",
@@ -23,16 +23,17 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 2,
+   "execution_count": null,
    "source": [
-    "metrics = [mtr.bias, mtr.rmse, mtr.urmse, mtr.mae, mtr.mape, mtr.mef, mtr.si, mtr.cc, mtr.spearmanr, mtr.r2, mtr.lin_slope]"
+    "metrics = [mtr.bias, mtr.max_error, mtr.rmse, mtr.urmse, mtr.mae, mtr.mape, \n",
+    "    mtr.mef, mtr.si, mtr.cc, mtr.spearmanr, mtr.r2, mtr.willmott, mtr.lin_slope]"
    ],
    "outputs": [],
    "metadata": {}
   },
   {
    "cell_type": "code",
-   "execution_count": 3,
+   "execution_count": null,
    "source": [
     "n = 50\n",
     "x = np.linspace(0.0, 6.0, num=n)\n",
@@ -45,9 +46,9 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 23,
+   "execution_count": null,
    "source": [
-    "def plot_metrics(bias, noise_level, fixed_y_axis=False):\n",
+    "def plot_metrics(bias, noise_level, fixed_y_axis=True):\n",
     "    y_mod = y_obs + bias + noise_level*noise_vec\n",
     "    plt.plot(x, y_obs, 'r.-', label=\"obs\")\n",
     "    plt.plot(x, y_mod, 'o-', label=\"model\")\n",
@@ -59,7 +60,7 @@
     "        ymax = 8\n",
     "        ymin = 1\n",
     "    ystep = 1.2*(ymax - ymin)/len(metrics)\n",
-    "    ypos = ymax\n",
+    "    ypos = ymax + 0.5\n",
     "    for m in metrics:\n",
     "        plt.text(6.5, ypos, f\"{m.__name__}:\")\n",
     "        plt.text(8.0, ypos, f\"{m(y_obs,y_mod):.4f}\")\n",
@@ -75,32 +76,10 @@
   },
   {
    "cell_type": "code",
-   "execution_count": 24,
+   "execution_count": null,
    "source": [
     "interact(plot_metrics, bias = (-1,3,0.1), noise_level=(0,2,0.05));"
    ],
-   "outputs": [
-    {
-     "output_type": "display_data",
-     "data": {
-      "application/vnd.jupyter.widget-view+json": {
-       "version_major": 2,
-       "version_minor": 0,
-       "model_id": "0666c7c25f1e449f9198d321ab95a84d"
-      },
-      "text/plain": [
-       "interactive(children=(FloatSlider(value=1.0, description='bias', max=3.0, min=-1.0), FloatSlider(value=1.0, de…"
-      ]
-     },
-     "metadata": {}
-    }
-   ],
-   "metadata": {}
-  },
-  {
-   "cell_type": "code",
-   "execution_count": null,
-   "source": [],
    "outputs": [],
    "metadata": {}
   }
@@ -109,7 +88,7 @@
   "orig_nbformat": 4,
   "language_info": {
    "name": "python",
-   "version": "3.9.6",
+   "version": "3.8.10",
    "mimetype": "text/x-python",
    "codemirror_mode": {
     "name": "ipython",
@@ -121,10 +100,10 @@
   },
   "kernelspec": {
    "name": "python3",
-   "display_name": "Python 3.9.6 64-bit"
+   "display_name": "Python 3.8.10 64-bit ('base': conda)"
   },
   "interpreter": {
-   "hash": "f4041ee05ab07c15354d6207e763f17a216c3f5ccf08906343c2b4fd3fa7a6fb"
+   "hash": "fa576ebcd40e010bdc0ae86b06ce09151f3424f9e9aed6893ff04f39a9299d89"
   }
  },
  "nbformat": 4,