Error Metrics

Error analysis and bias quantification for model evaluation.

Error Metrics for Model Evaluation

COE(obs, mod, axis=None)

Center of Mass Error (COE).

The COE measures the displacement between the centroids (centers of mass) of two fields. For spatial data, this represents the shift in the center of a feature (e.g., a storm or a pollutant plume).

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values (typically 2D spatial field). mod : numpy.ndarray or xarray.DataArray Model or predicted values (typically 2D spatial field). axis : int, str, or iterable of such, optional Axis or dimension(s) over which to compute the centroid. If None, computes over all axes.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Center of mass error (Euclidean distance between centroids).

Examples

import numpy as np from monet_stats.error_metrics import COE obs = np.zeros((5, 5)) obs[2, 2] = 1.0 # Peak at center (2, 2) mod = np.zeros((5, 5)) mod[3, 3] = 1.0 # Peak shifted to (3, 3)

Displacement is sqrt(1^2 + 1^2) = sqrt(2) approx 1.414

np.allclose(COE(obs, mod), np.sqrt(2)) True

Source code in src/monet_stats/error_metrics.py

def COE(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Center of Mass Error (COE).

    The COE measures the displacement between the centroids (centers of mass)
    of two fields. For spatial data, this represents the shift in the center
    of a feature (e.g., a storm or a pollutant plume).

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values (typically 2D spatial field).
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values (typically 2D spatial field).
    axis : int, str, or iterable of such, optional
        Axis or dimension(s) over which to compute the centroid.
        If None, computes over all axes.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Center of mass error (Euclidean distance between centroids).

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import COE
    >>> obs = np.zeros((5, 5))
    >>> obs[2, 2] = 1.0  # Peak at center (2, 2)
    >>> mod = np.zeros((5, 5))
    >>> mod[3, 3] = 1.0  # Peak shifted to (3, 3)
    >>> # Displacement is sqrt(1^2 + 1^2) = sqrt(2) approx 1.414
    >>> np.allclose(COE(obs, mod), np.sqrt(2))
    True
    """
    from .utils_stats import _update_history

    def _get_centroid(da: xr.DataArray, dims: Iterable[str]) -> List[xr.DataArray]:
        """Helper to calculate centroid of a DataArray."""
        total = da.sum(dim=dims)
        # Handle zero sum to avoid division by zero
        total_safe = xr.where(total == 0, 1e-10, total)
        coords_list = []
        for d in dims:
            # Check if coord exists and is numeric
            if d in da.coords and np.issubdtype(da.coords[d].dtype, np.number):
                coord = da.coords[d]
            else:
                # Fallback to dimension indices
                coord = xr.DataArray(np.arange(da.sizes[d]), dims=d, name=d)
            # Weighted mean of coordinate
            c = (da * coord).sum(dim=dims) / total_safe
            coords_list.append(c)
        return coords_list

    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        if axis is None:
            dims = list(obs.dims)
        elif isinstance(axis, (int, str)):
            dims = [obs.dims[axis] if isinstance(axis, int) else axis]
        else:
            dims = [obs.dims[d] if isinstance(d, int) else d for d in axis]

        c_obs = _get_centroid(obs, dims)
        c_mod = _get_centroid(mod, dims)

        # Euclidean distance
        dist_sq = sum((cm - co) ** 2 for cm, co in zip(c_mod, c_obs))
        result = dist_sq**0.5

        return _update_history(result, "Center of Mass Error (COE)")

    # Fallback to numpy
    obs_arr = np.asanyarray(obs)
    mod_arr = np.asanyarray(mod)

    if axis is None:
        axes = tuple(range(obs_arr.ndim))
    elif isinstance(axis, int):
        axes = (axis,)
    else:
        axes = tuple(axis)

    def _get_numpy_centroid(arr: np.ndarray, axes_tuple: Tuple[int, ...]) -> List[np.ndarray]:
        """Helper to calculate centroid of a NumPy array."""
        total = np.sum(arr, axis=axes_tuple)
        total_safe = np.where(total == 0, 1e-10, total)
        c_list = []
        for ax in axes_tuple:
            # Create coordinate array for this axis
            shape = [1] * arr.ndim
            shape[ax] = arr.shape[ax]
            coord = np.arange(arr.shape[ax]).reshape(shape)
            c = np.sum(arr * coord, axis=axes_tuple) / total_safe
            c_list.append(c)
        return c_list

    c_obs_np = _get_numpy_centroid(obs_arr, axes)
    c_mod_np = _get_numpy_centroid(mod_arr, axes)

    dist_sq_np = sum((cm - co) ** 2 for cm, co in zip(c_mod_np, c_obs_np))
    return dist_sq_np**0.5

CORR_INDEX(obs, mod, axis=None)

Correlation Index (CORR_INDEX).

Typical Use Cases

• Measuring the linear relationship between observed and modeled values.
• Used as a component in model evaluation.
• Quantifies how well model captures observed patterns.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute correlation index.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Correlation index (unitless, -1 to 1).

Examples

import numpy as np from monet_stats.error_metrics import CORR_INDEX obs = np.array([1, 2, 3, 4]) mod = np.array([2, 4, 6, 8]) CORR_INDEX(obs, mod) 1.0

Source code in src/monet_stats/error_metrics.py

def CORR_INDEX(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Correlation Index (CORR_INDEX).

    Typical Use Cases
    -----------------
    - Measuring the linear relationship between observed and modeled values.
    - Used as a component in model evaluation.
    - Quantifies how well model captures observed patterns.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute correlation index.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Correlation index (unitless, -1 to 1).

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import CORR_INDEX
    >>> obs = np.array([1, 2, 3, 4])
    >>> mod = np.array([2, 4, 6, 8])
    >>> CORR_INDEX(obs, mod)
    1.0
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        # Using xarray's built-in correlation function
        result = xr.corr(obs, mod, dim=dim)
        # Update history
        history = f"CORR_INDEX computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        # Fallback to numpy-compatible logic
        obs = np.asarray(obs)
        mod = np.asarray(mod)
        if axis is None:
            from scipy.stats import pearsonr

            return pearsonr(obs.flatten(), mod.flatten())[0]
        else:
            # Manual vectorized correlation over axis for robustness across scipy versions
            obs_mean = np.mean(obs, axis=axis, keepdims=True)
            mod_mean = np.mean(mod, axis=axis, keepdims=True)
            obs_std = obs - obs_mean
            mod_std = mod - mod_mean
            num = np.sum(obs_std * mod_std, axis=axis)
            den = np.sqrt(np.sum(obs_std**2, axis=axis) * np.sum(mod_std**2, axis=axis))
            return num / den

CRMSE(obs, mod, axis=None)

Centered Root Mean Square Error (CRMSE).

Typical Use Cases

• Quantifying the error between anomalies (deviations from mean) of model and observations.
• Used in Taylor diagrams, model evaluation, and forecast verification.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute CRMSE.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Centered root mean square error.

Examples

import numpy as np from monet_stats.error_metrics import CRMSE obs = np.array([1, 2, 3]) mod = np.array([2, 2, 4]) CRMSE(obs, mod) 0.4714045207910317

Source code in src/monet_stats/error_metrics.py

def CRMSE(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Centered Root Mean Square Error (CRMSE).

    Typical Use Cases
    -----------------
    - Quantifying the error between anomalies (deviations from mean) of model
      and observations.
    - Used in Taylor diagrams, model evaluation, and forecast verification.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute CRMSE.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Centered root mean square error.

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import CRMSE
    >>> obs = np.array([1, 2, 3])
    >>> mod = np.array([2, 2, 4])
    >>> CRMSE(obs, mod)
    0.4714045207910317
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        o_ = obs - obs.mean(dim=dim)
        m_ = mod - mod.mean(dim=dim)
        result = ((m_ - o_) ** 2).mean(dim=dim, keep_attrs=True) ** 0.5
        # Update history
        history = f"CRMSE computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        o_ = np.subtract(obs, np.mean(obs, axis=axis, keepdims=True))
        m_ = np.subtract(mod, np.mean(mod, axis=axis, keepdims=True))
        return (np.ma.abs(m_ - o_) ** 2).mean(axis=axis) ** 0.5

IOA(obs, mod, axis=None)

Index of Agreement (IOA).

Typical Use Cases

• Quantifying the agreement between model and observations, normalized by total deviation.
• Used in model evaluation for skill assessment.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute IOA.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Index of agreement (unitless, 0-1).

Examples

import numpy as np from monet_stats.error_metrics import IOA obs = np.array([1, 2, 3]) mod = np.array([2, 2, 4]) IOA(obs, mod) 0.8

Source code in src/monet_stats/error_metrics.py

def IOA(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Index of Agreement (IOA).

    Typical Use Cases
    -----------------
    - Quantifying the agreement between model and observations, normalized by
      total deviation.
    - Used in model evaluation for skill assessment.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute IOA.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Index of agreement (unitless, 0-1).

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import IOA
    >>> obs = np.array([1, 2, 3])
    >>> mod = np.array([2, 2, 4])
    >>> IOA(obs, mod)
    0.8
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        obs_mean = obs.mean(dim=dim)
        num = ((obs - mod) ** 2).sum(dim=dim)
        denom = ((abs(mod - obs_mean) + abs(obs - obs_mean)) ** 2).sum(dim=dim)
        result = 1.0 - (num / denom)
        # Update history
        history = f"IOA computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        obs_mean = np.mean(obs, axis=axis, keepdims=True)
        num = np.sum((obs - mod) ** 2, axis=axis)
        denom = np.sum((np.abs(mod - obs_mean) + np.abs(obs - obs_mean)) ** 2, axis=axis)
        return 1.0 - (num / denom)

IOA_m(obs, mod, axis=None)

Index of Agreement (IOA) - robust to masked arrays.

Typical Use Cases

• Quantifying the agreement between model and observations, normalized by total deviation, robust to missing data.
• Used in model evaluation for skill assessment with incomplete datasets.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute IOA.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Index of agreement (unitless, 0-1).

Examples

import numpy as np from monet_stats.error_metrics import IOA_m obs = np.array([1, 2, 3]) mod = np.array([2, 2, 4]) IOA_m(obs, mod) 0.8

Source code in src/monet_stats/error_metrics.py

def IOA_m(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Index of Agreement (IOA) - robust to masked arrays.

    Typical Use Cases
    -----------------
    - Quantifying the agreement between model and observations, normalized by
      total deviation, robust to missing data.
    - Used in model evaluation for skill assessment with incomplete datasets.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute IOA.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Index of agreement (unitless, 0-1).

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import IOA_m
    >>> obs = np.array([1, 2, 3])
    >>> mod = np.array([2, 2, 4])
    >>> IOA_m(obs, mod)
    0.8
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        # IOA implementation for xarray already handles NaNs
        return IOA(obs, mod, axis=axis)
    else:
        obs_mean = np.ma.mean(obs, axis=axis, keepdims=True)
        num = np.ma.sum((obs - mod) ** 2, axis=axis)
        denom = np.ma.sum((np.ma.abs(mod - obs_mean) + np.ma.abs(obs - obs_mean)) ** 2, axis=axis)
        return 1.0 - (num / denom)

LOG_ERROR(obs, mod, axis=None)

Logarithmic Error Metric.

Typical Use Cases

• Quantifying errors for variables that span several orders of magnitude.
• Used in atmospheric sciences for concentration data (e.g., pollutants).
• Helpful when relative rather than absolute errors are important.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values (should be positive). mod : numpy.ndarray or xarray.DataArray Model or predicted values (should be positive). axis : int, str, or iterable of such, optional Axis or dimension along which to compute log error.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Logarithmic error metric.

Examples

import numpy as np from monet_stats.error_metrics import LOG_ERROR obs = np.array([1, 100]) mod = np.array([2, 200]) LOG_ERROR(obs, mod) 0.34657359027997264

Source code in src/monet_stats/error_metrics.py

def LOG_ERROR(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Logarithmic Error Metric.

    Typical Use Cases
    -----------------
    - Quantifying errors for variables that span several orders of magnitude.
    - Used in atmospheric sciences for concentration data (e.g., pollutants).
    - Helpful when relative rather than absolute errors are important.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values (should be positive).
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values (should be positive).
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute log error.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Logarithmic error metric.

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import LOG_ERROR
    >>> obs = np.array([1, 100])
    >>> mod = np.array([2, 200])
    >>> LOG_ERROR(obs, mod)
    0.34657359027997264
    """
    # Add small epsilon to avoid log(0) and handle negative values
    epsilon = 1e-10

    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        # Use abs to handle potential negative values, then add epsilon
        obs_safe = abs(obs) + epsilon
        mod_safe = abs(mod) + epsilon
        obs_log = np.log(obs_safe)
        mod_log = np.log(mod_safe)
        result = ((mod_log - obs_log) ** 2).mean(dim=dim, keep_attrs=True) ** 0.5
        # Update history
        history = f"LOG_ERROR computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        # Use abs to handle potential negative values, then add epsilon
        obs_safe = np.abs(obs) + epsilon
        mod_safe = np.abs(mod) + epsilon
        obs_log = np.log(obs_safe)
        mod_log = np.log(mod_safe)

        result = np.sqrt(np.mean((mod_log - obs_log) ** 2, axis=axis))
        # Return 0 for perfect agreement
        if np.array_equal(obs, mod):
            return 0.0
        return result

MAE(obs, mod, axis=None)

Mean Absolute Error (MAE).

Typical Use Cases

• Quantifying the average magnitude of errors between model and observations, regardless of direction.
• Used in model evaluation, forecast verification, and regression analysis.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute MAE.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Mean absolute error.

Examples

import numpy as np from monet_stats.error_metrics import MAE obs = np.array([1, 2, 3]) mod = np.array([2, 2, 4]) MAE(obs, mod) 0.6666666666666666

Source code in src/monet_stats/error_metrics.py

def MAE(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Mean Absolute Error (MAE).

    Typical Use Cases
    -----------------
    - Quantifying the average magnitude of errors between model and observations,
      regardless of direction.
    - Used in model evaluation, forecast verification, and regression analysis.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute MAE.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Mean absolute error.

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import MAE
    >>> obs = np.array([1, 2, 3])
    >>> mod = np.array([2, 2, 4])
    >>> MAE(obs, mod)
    0.6666666666666666
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        result = abs(mod - obs).mean(dim=dim, keep_attrs=True)
        # Update history
        history = f"MAE computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        return np.ma.abs(np.subtract(mod, obs)).mean(axis=axis)

MAE_m(obs, mod, axis=None)

Mean Absolute Error (MAE) - robust to masked arrays.

Typical Use Cases

• Quantifying the average magnitude of errors between model and observations, regardless of direction, robust to missing data.
• Used in model evaluation, forecast verification, and regression analysis with incomplete datasets.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute MAE.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Mean absolute error.

Examples

import numpy as np from monet_stats.error_metrics import MAE_m obs = np.array([1, 2, 3]) mod = np.array([2, 2, 4]) MAE_m(obs, mod) 0.6666666666666666

Source code in src/monet_stats/error_metrics.py

def MAE_m(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Mean Absolute Error (MAE) - robust to masked arrays.

    Typical Use Cases
    -----------------
    - Quantifying the average magnitude of errors between model and
      observations, regardless of direction, robust to missing data.
    - Used in model evaluation, forecast verification, and regression analysis
      with incomplete datasets.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute MAE.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Mean absolute error.

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import MAE_m
    >>> obs = np.array([1, 2, 3])
    >>> mod = np.array([2, 2, 4])
    >>> MAE_m(obs, mod)
    0.6666666666666666
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        # MAE implementation for xarray already handles NaNs
        return MAE(obs, mod, axis=axis)
    else:
        return np.ma.mean(np.ma.abs(np.subtract(mod, obs)), axis=axis)

MAE_norm(obs, mod, axis=None)

Normalized Mean Absolute Error (MAE_norm).

Normalizes MAE by the range of observations.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute normalized MAE.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Normalized mean absolute error (unitless).

Source code in src/monet_stats/error_metrics.py

def MAE_norm(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Normalized Mean Absolute Error (MAE_norm).

    Normalizes MAE by the range of observations.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute normalized MAE.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Normalized mean absolute error (unitless).
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        mae = abs(mod - obs).mean(dim=dim, keep_attrs=True)
        obs_min = obs.min(dim=dim)
        obs_max = obs.max(dim=dim)
        obs_range = obs_max - obs_min
        # Avoid division by zero
        result = xr.where(obs_range == 0, mae, mae / obs_range)
        # Update history
        history = f"MAE_norm computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        mae = np.mean(np.abs(np.subtract(mod, obs)), axis=axis)
        obs_min = np.min(obs, axis=axis)
        obs_max = np.max(obs, axis=axis)
        obs_range = obs_max - obs_min
        # Avoid division by zero
        result = np.where(obs_range == 0, mae, mae / obs_range)
        return result.item() if np.ndim(result) == 0 else result

MAPE(obs, mod, axis=None)

Mean Absolute Percentage Error (MAPE).

Typical Use Cases

• Quantifying the average relative error between model and observations as a percentage.
• Used in time series forecasting, regression, and model evaluation for percentage-based error assessment.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute MAPE.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Mean absolute percentage error (in percent).

Examples

import numpy as np from monet_stats.error_metrics import MAPE obs = np.array([1, 2, 3]) mod = np.array([2, 2, 4]) MAPE(obs, mod) 50.0

Source code in src/monet_stats/error_metrics.py

def MAPE(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Mean Absolute Percentage Error (MAPE).

    Typical Use Cases
    -----------------
    - Quantifying the average relative error between model and observations
      as a percentage.
    - Used in time series forecasting, regression, and model evaluation for
      percentage-based error assessment.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute MAPE.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Mean absolute percentage error (in percent).

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import MAPE
    >>> obs = np.array([1, 2, 3])
    >>> mod = np.array([2, 2, 4])
    >>> MAPE(obs, mod)
    50.0
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        result = (100 * abs(mod - obs) / abs(obs)).mean(dim=dim, keep_attrs=True)
        # Update history
        history = f"MAPE computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        return (100 * np.ma.abs(np.subtract(mod, obs)) / np.ma.abs(obs)).mean(axis=axis)

MAPE_mod(obs, mod, axis=None)

Modified Mean Absolute Percentage Error (MAPE).

This version handles cases where observations might be zero or near zero by using a small epsilon to avoid division by zero.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute MAPE.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Mean absolute percentage error (in percent).

Source code in src/monet_stats/error_metrics.py

def MAPE_mod(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Modified Mean Absolute Percentage Error (MAPE).

    This version handles cases where observations might be zero or near zero
    by using a small epsilon to avoid division by zero.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute MAPE.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Mean absolute percentage error (in percent).
    """
    # Small epsilon to avoid division by zero
    epsilon = 1e-8

    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        # Add epsilon to avoid division by zero
        obs_safe = xr.where(abs(obs) < epsilon, epsilon, obs)
        result = (100 * abs(mod - obs) / abs(obs_safe)).mean(dim=dim, keep_attrs=True)
        # Update history
        history = f"MAPE_mod computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        # Add epsilon to avoid division by zero
        obs_safe = np.where(np.abs(obs) < epsilon, epsilon, obs)
        return (100 * np.abs(np.subtract(mod, obs)) / np.abs(obs_safe)).mean(axis=axis)

MAPEm(obs, mod, axis=None)

Mean Absolute Percentage Error (MAPE) - robust to masked arrays.

Typical Use Cases

• Quantifying the average relative error between model and observations as a percentage, robust to missing data.
• Used in time series forecasting, regression, and model evaluation for percentage-based error assessment.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute MAPE.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Mean absolute percentage error (in percent).

Examples

import numpy as np from monet_stats.error_metrics import MAPEm obs = np.array([1, 2, 3]) mod = np.array([2, 2, 4]) MAPEm(obs, mod) 50.0

Source code in src/monet_stats/error_metrics.py

def MAPEm(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Mean Absolute Percentage Error (MAPE) - robust to masked arrays.

    Typical Use Cases
    -----------------
    - Quantifying the average relative error between model and observations as
      a percentage, robust to missing data.
    - Used in time series forecasting, regression, and model evaluation for
      percentage-based error assessment.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute MAPE.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Mean absolute percentage error (in percent).

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import MAPEm
    >>> obs = np.array([1, 2, 3])
    >>> mod = np.array([2, 2, 4])
    >>> MAPEm(obs, mod)
    50.0
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        # MAPE implementation for xarray already handles NaNs
        return MAPE(obs, mod, axis=axis)
    else:
        return 100 * np.ma.mean(np.ma.abs((mod - obs) / obs), axis=axis)

MASE(obs, mod, axis=None)

Mean Absolute Scaled Error (MASE).

Typical Use Cases

• Quantifying model error relative to the error of a simple baseline model (e.g., naive forecast).
• Used in time series forecasting and model evaluation.
• Provides scale-independent comparison across different datasets.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute MASE.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Mean absolute scaled error (unitless).

Examples

import numpy as np from monet_stats.error_metrics import MASE obs = np.array([1, 2, 3, 4]) mod = np.array([1.1, 2.1, 3.1, 4.1]) MASE(obs, mod) 0.1

Source code in src/monet_stats/error_metrics.py

def MASE(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Mean Absolute Scaled Error (MASE).

    Typical Use Cases
    -----------------
    - Quantifying model error relative to the error of a simple baseline model
      (e.g., naive forecast).
    - Used in time series forecasting and model evaluation.
    - Provides scale-independent comparison across different datasets.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute MASE.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Mean absolute scaled error (unitless).

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import MASE
    >>> obs = np.array([1, 2, 3, 4])
    >>> mod = np.array([1.1, 2.1, 3.1, 4.1])
    >>> MASE(obs, mod)
    0.1
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis

        # Calculate naive forecast error (using previous observation)
        # Assuming 'time' is a dimension for shift
        # If 'time' is not present, we use the provided dimension or axis
        if "time" in obs.dims:
            naive_error = abs(obs - obs.shift(time=1)).mean(dim=dim, skipna=True)
        else:
            # Fallback if time is not named 'time'
            naive_error = abs(obs - obs.shift({obs.dims[0]: 1})).mean(dim=dim, skipna=True)

        model_error = abs(mod - obs).mean(dim=dim, keep_attrs=True)
        result = model_error / naive_error
        # Update history
        history = f"MASE computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        # Calculate naive forecast error (using previous observation)
        if axis is not None:
            naive_diff = np.diff(obs, axis=axis)
            naive_error = np.mean(np.abs(naive_diff), axis=axis)
        else:
            naive_diff = np.diff(obs)
            naive_error = np.mean(np.abs(naive_diff))
        model_error = np.mean(np.abs(np.subtract(mod, obs)), axis=axis)
        return model_error / naive_error

MASE_mod(obs, mod, axis=None)

Modified Mean Absolute Scaled Error (MASE).

This version handles cases where the naive forecast error is zero by using a small epsilon to avoid division by zero.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute MASE.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Mean absolute scaled error (unitless).

Source code in src/monet_stats/error_metrics.py

def MASE_mod(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Modified Mean Absolute Scaled Error (MASE).

    This version handles cases where the naive forecast error is zero
    by using a small epsilon to avoid division by zero.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute MASE.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Mean absolute scaled error (unitless).
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        # Calculate naive forecast error (using previous observation)
        if "time" in obs.dims:
            naive_error = abs(obs - obs.shift(time=1)).mean(dim=dim, skipna=True)
        else:
            naive_error = abs(obs - obs.shift({obs.dims[0]: 1})).mean(dim=dim, skipna=True)

        model_error = abs(mod - obs).mean(dim=dim, keep_attrs=True)
        # Avoid division by zero
        result = xr.where(naive_error == 0, model_error, model_error / naive_error)
        # Update history
        history = f"MASE_mod computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        # Calculate naive forecast error (using previous observation)
        if axis is not None:
            naive_diff = np.diff(obs, axis=axis)
            naive_error = np.mean(np.abs(naive_diff), axis=axis)
        else:
            naive_diff = np.diff(obs)
            naive_error = np.mean(np.abs(naive_diff))
        model_error = np.mean(np.abs(np.subtract(mod, obs)), axis=axis)
        # Avoid division by zero
        result = np.where(naive_error == 0, model_error, model_error / naive_error)
        return result.item() if np.ndim(result) == 0 else result

MASEm(obs, mod, axis=None)

Mean Absolute Scaled Error (MASE) - robust to masked arrays.

Typical Use Cases

• Quantifying model error relative to the error of a simple baseline model (e.g., naive forecast), robust to masked arrays.
• Used in time series forecasting and model evaluation with missing data.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute MASE.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Mean absolute scaled error (unitless).

Examples

import numpy as np from monet_stats.error_metrics import MASEm obs = np.array([1, 2, 3, 4]) mod = np.array([1.1, 2.1, 3.1, 4.1]) MASEm(obs, mod) 0.1

Source code in src/monet_stats/error_metrics.py

def MASEm(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Mean Absolute Scaled Error (MASE) - robust to masked arrays.

    Typical Use Cases
    -----------------
    - Quantifying model error relative to the error of a simple baseline model
      (e.g., naive forecast), robust to masked arrays.
    - Used in time series forecasting and model evaluation with missing data.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute MASE.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Mean absolute scaled error (unitless).

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import MASEm
    >>> obs = np.array([1, 2, 3, 4])
    >>> mod = np.array([1.1, 2.1, 3.1, 4.1])
    >>> MASEm(obs, mod)
    0.1
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        # MASE implementation for xarray already handles NaNs with skipna=True
        return MASE(obs, mod, axis=axis)
    else:
        # Calculate naive forecast error (using previous observation) with masked arrays
        if axis is not None:
            # Use numpy's gradient-like approach for masked arrays
            naive_diff = np.ma.diff(obs, axis=axis)
            naive_error = np.ma.mean(np.ma.abs(naive_diff), axis=axis)
        else:
            naive_diff = np.ma.diff(obs)
            naive_error = np.ma.mean(np.ma.abs(naive_diff))
        model_error = np.ma.mean(np.ma.abs(np.subtract(mod, obs)), axis=axis)
        return model_error / naive_error

MB(obs, mod, axis=None)

Mean Bias (MB).

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute the mean bias.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Mean bias value(s) = mean(model - observation). Positive values indicate model overestimation.

Source code in src/monet_stats/error_metrics.py

def MB(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Mean Bias (MB).

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute the mean bias.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Mean bias value(s) = mean(model - observation).
        Positive values indicate model overestimation.
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        result = (mod - obs).mean(dim=dim, keep_attrs=True)
        # Update history
        history = f"MB computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        return np.ma.mean(np.subtract(mod, obs), axis=axis)

MNB(obs, mod, axis=None)

Mean Normalized Bias (%).

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute the bias.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Mean normalized bias (percent).

Source code in src/monet_stats/error_metrics.py

def MNB(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Mean Normalized Bias (%).

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute the bias.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Mean normalized bias (percent).
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        result = ((mod - obs) / obs).mean(dim=dim, keep_attrs=True) * 100.0
        # Update history
        history = f"MNB computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        return np.ma.masked_invalid((mod - obs) / obs).mean(axis=axis) * 100.0

MNE(obs, mod, axis=None)

Mean Normalized Gross Error (%).

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute the error.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Mean normalized gross error (percent).

Source code in src/monet_stats/error_metrics.py

def MNE(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Mean Normalized Gross Error (%).

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute the error.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Mean normalized gross error (percent).
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        result = (abs(mod - obs) / obs).mean(dim=dim, keep_attrs=True) * 100.0
        # Update history
        history = f"MNE computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        return np.ma.masked_invalid(np.ma.abs(mod - obs) / obs).mean(axis=axis) * 100.0

MO(obs, mod, axis=None)

Mean Error (MO) - Mean of (model - observation).

Typical Use Cases

• Quantifying the average bias between model predictions and observations.
• Used in model evaluation to assess systematic errors.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute the mean error.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Mean error (model - observation) in observation units. Returns 0.0 for perfect agreement.

Examples

import numpy as np from monet_stats.error_metrics import MO obs = np.array([1, 2, 3, 4, 5]) mod = np.array([1.1, 2.1, 3.1, 4.1, 5.1]) MO(obs, mod) 0.1

Source code in src/monet_stats/error_metrics.py

def MO(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Mean Error (MO) - Mean of (model - observation).

    Typical Use Cases
    -----------------
    - Quantifying the average bias between model predictions and observations.
    - Used in model evaluation to assess systematic errors.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute the mean error.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Mean error (model - observation) in observation units.
        Returns 0.0 for perfect agreement.

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import MO
    >>> obs = np.array([1, 2, 3, 4, 5])
    >>> mod = np.array([1.1, 2.1, 3.1, 4.1, 5.1])
    >>> MO(obs, mod)
    0.1
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        result = (mod - obs).mean(dim=dim, keep_attrs=True)
        # Update history
        history = f"MO computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        return np.mean(np.subtract(mod, obs), axis=axis)

MP(obs=None, mod=None, axis=None)

Mean Predictions (model unit).

Typical Use Cases

• Calculating the average value of model predictions for baseline or climatological reference.
• Used in normalization, anomaly calculation, and summary statistics for model output.

Parameters

obs : numpy.ndarray or xarray.DataArray, optional Observed values (not used for MP but included for signature matching). mod : numpy.ndarray or xarray.DataArray Model predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute the mean.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Mean of predictions.

Source code in src/monet_stats/error_metrics.py

def MP(
    obs: Optional[Union[np.ndarray, xr.DataArray]] = None,
    mod: Union[np.ndarray, xr.DataArray] = None,
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Mean Predictions (model unit).

    Typical Use Cases
    -----------------
    - Calculating the average value of model predictions for baseline or
      climatological reference.
    - Used in normalization, anomaly calculation, and summary statistics for
      model output.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray, optional
        Observed values (not used for MP but included for signature matching).
    mod : numpy.ndarray or xarray.DataArray
        Model predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute the mean.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Mean of predictions.
    """
    if isinstance(mod, xr.DataArray):
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = mod.dims[axis]
        else:
            dim = axis
        result = mod.mean(dim=dim, keep_attrs=True)
        # Update history
        history = f"MP computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        return np.mean(mod, axis=axis)

MdnB(obs, mod, axis=None)

Median Bias (MdnB).

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute the median bias.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Median bias value(s) = median(model - observation). Positive values indicate model overestimation.

Source code in src/monet_stats/error_metrics.py

def MdnB(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Median Bias (MdnB).

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute the median bias.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Median bias value(s) = median(model - observation).
        Positive values indicate model overestimation.
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        result = (mod - obs).median(dim=dim, keep_attrs=True)
        # Update history
        history = f"MdnB computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        return np.ma.median(np.subtract(mod, obs), axis=axis)

MdnNB(obs, mod, axis=None)

Median Normalized Bias (%).

Typical Use Cases

• Assessing the central tendency of model bias relative to observations, less sensitive to outliers than mean.
• Useful for robust model evaluation in the presence of skewed or non-normal error distributions.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute the bias.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Median normalized bias (percent).

Source code in src/monet_stats/error_metrics.py

def MdnNB(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Median Normalized Bias (%).

    Typical Use Cases
    -----------------
    - Assessing the central tendency of model bias relative to observations,
      less sensitive to outliers than mean.
    - Useful for robust model evaluation in the presence of skewed or non-normal
      error distributions.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute the bias.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Median normalized bias (percent).
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        result = ((mod - obs) / obs).median(dim=dim, keep_attrs=True) * 100.0
        # Update history
        history = f"MdnNB computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        return np.ma.median(np.ma.masked_invalid((mod - obs) / obs), axis=axis) * 100.0

MdnNE(obs, mod, axis=None)

Median Normalized Gross Error (%).

Typical Use Cases

• Evaluating the typical magnitude of model errors relative to observations, robust to outliers.
• Useful for summarizing error magnitude in non-Gaussian or heavy-tailed error distributions.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute the error.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Median normalized gross error (percent).

Source code in src/monet_stats/error_metrics.py

def MdnNE(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Median Normalized Gross Error (%).

    Typical Use Cases
    -----------------
    - Evaluating the typical magnitude of model errors relative to observations,
      robust to outliers.
    - Useful for summarizing error magnitude in non-Gaussian or heavy-tailed
      error distributions.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute the error.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Median normalized gross error (percent).
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        result = (abs(mod - obs) / obs).median(dim=dim, keep_attrs=True) * 100.0
        # Update history
        history = f"MdnNE computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        return np.ma.median(np.ma.masked_invalid(np.ma.abs(mod - obs) / obs), axis=axis) * 100.0

MdnO(obs, mod, axis=None)

Median Error (MdnO) - Median of (model - observation).

Typical Use Cases

• Quantifying the typical bias between model predictions and observations, robust to outliers.
• Used in robust model evaluation for non-parametric error assessment.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute the median error.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Median error (model - observation) in observation units. Returns 0.0 for perfect agreement.

Examples

import numpy as np from monet_stats.error_metrics import MdnO obs = np.array([1, 2, 3, 4, 5]) mod = np.array([1.1, 2.1, 3.1, 4.1, 5.1]) MdnO(obs, mod) 0.1

Source code in src/monet_stats/error_metrics.py

def MdnO(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Median Error (MdnO) - Median of (model - observation).

    Typical Use Cases
    -----------------
    - Quantifying the typical bias between model predictions and observations,
      robust to outliers.
    - Used in robust model evaluation for non-parametric error assessment.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute the median error.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Median error (model - observation) in observation units.
        Returns 0.0 for perfect agreement.

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import MdnO
    >>> obs = np.array([1, 2, 3, 4, 5])
    >>> mod = np.array([1.1, 2.1, 3.1, 4.1, 5.1])
    >>> MdnO(obs, mod)
    0.1
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        result = (mod - obs).median(dim=dim, keep_attrs=True)
        # Update history
        history = f"MdnO computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        return np.median(np.subtract(mod, obs), axis=axis)

MdnP(obs, mod, axis=None)

Median Error (MdnP) - Median of (model - observation).

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute the median error.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Median error (model - observation) in model units. Returns 0.0 for perfect agreement.

Source code in src/monet_stats/error_metrics.py

def MdnP(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Median Error (MdnP) - Median of (model - observation).

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute the median error.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Median error (model - observation) in model units.
        Returns 0.0 for perfect agreement.
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        result = (mod - obs).median(dim=dim, keep_attrs=True)
        # Update history
        history = f"MdnP computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        return np.median(np.subtract(mod, obs), axis=axis)

MedAE(obs, mod, axis=None)

Median Absolute Error (MedAE).

Typical Use Cases

• Evaluating the typical magnitude of errors, robust to outliers and non-normal error distributions.
• Used in robust regression, model evaluation, and forecast verification.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute MedAE.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Median absolute error.

Examples

import numpy as np from monet_stats.error_metrics import MedAE obs = np.array([1, 2, 3]) mod = np.array([2, 2, 4]) MedAE(obs, mod) 1.0

Source code in src/monet_stats/error_metrics.py

def MedAE(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Median Absolute Error (MedAE).

    Typical Use Cases
    -----------------
    - Evaluating the typical magnitude of errors, robust to outliers and
      non-normal error distributions.
    - Used in robust regression, model evaluation, and forecast verification.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute MedAE.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Median absolute error.

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import MedAE
    >>> obs = np.array([1, 2, 3])
    >>> mod = np.array([2, 2, 4])
    >>> MedAE(obs, mod)
    1.0
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        result = abs(mod - obs).median(dim=dim, keep_attrs=True)
        # Update history
        history = f"MedAE computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        return np.ma.median(np.ma.abs(np.subtract(mod, obs)), axis=axis)

MedAE_m(obs, mod, axis=None)

Median Absolute Error (MedAE) - robust to masked arrays and outliers.

Typical Use Cases

• Evaluating the typical magnitude of errors, robust to outliers and non-normal error distributions with missing data.
• Used in robust regression, model evaluation, and forecast verification with incomplete datasets.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute MedAE.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Median absolute error.

Examples

import numpy as np from monet_stats.error_metrics import MedAE_m obs = np.array([1, 2, 3]) mod = np.array([2, 2, 4]) MedAE_m(obs, mod) 1.0

Source code in src/monet_stats/error_metrics.py

def MedAE_m(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Median Absolute Error (MedAE) - robust to masked arrays and outliers.

    Typical Use Cases
    -----------------
    - Evaluating the typical magnitude of errors, robust to outliers and
      non-normal error distributions with missing data.
    - Used in robust regression, model evaluation, and forecast verification
      with incomplete datasets.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute MedAE.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Median absolute error.

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import MedAE_m
    >>> obs = np.array([1, 2, 3])
    >>> mod = np.array([2, 2, 4])
    >>> MedAE_m(obs, mod)
    1.0
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        # MedAE implementation for xarray already handles NaNs
        return MedAE(obs, mod, axis=axis)
    else:
        return np.ma.median(np.ma.abs(np.subtract(mod, obs)), axis=axis)

NMSE(obs, mod, axis=None)

Normalized Mean Square Error (NMSE).

Typical Use Cases

• Quantifying the normalized squared error between model and observations.
• Used in model evaluation to compare performance across different variables or sites with different scales.
• Provides dimensionless error metric for cross-comparison.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute NMSE.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Normalized mean square error (unitless).

Examples

import numpy as np from monet_stats.error_metrics import NMSE obs = np.array([1, 2, 3, 4]) mod = np.array([2, 2, 2, 2]) NMSE(obs, mod) 0.25

Source code in src/monet_stats/error_metrics.py

def NMSE(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Normalized Mean Square Error (NMSE).

    Typical Use Cases
    -----------------
    - Quantifying the normalized squared error between model and observations.
    - Used in model evaluation to compare performance across different variables
      or sites with different scales.
    - Provides dimensionless error metric for cross-comparison.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute NMSE.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Normalized mean square error (unitless).

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import NMSE
    >>> obs = np.array([1, 2, 3, 4])
    >>> mod = np.array([2, 2, 2, 2])
    >>> NMSE(obs, mod)
    0.25
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        mse = ((mod - obs) ** 2).mean(dim=dim, keep_attrs=True)
        obs_var = obs.var(dim=dim)
        # Handle case where variance is 0 (perfect agreement)
        result = xr.where(obs_var == 0, 0, mse / obs_var)
        # Update history
        history = f"NMSE computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        mse = np.mean((np.subtract(mod, obs)) ** 2, axis=axis)
        obs_var = np.var(obs, axis=axis)
        # Handle case where variance is 0 (perfect agreement)
        result = np.where(obs_var == 0, 0, mse / obs_var)
        return result.item() if np.ndim(result) == 0 else result

NMdnGE(obs, mod, axis=None)

Normalized Median Gross Error (%).

Typical Use Cases

• Comparing the typical (median) error magnitude, normalized by the mean observation, for robust model evaluation.
• Useful for inter-comparison of model performance across sites or variables with different scales.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute the error.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Normalized median gross error (percent).

Examples

import numpy as np from monet_stats.error_metrics import NMdnGE obs = np.array([1, 2, 3, 4, 100]) mod = np.array([1.1, 2.1, 3.1, 4.1, 105]) NMdnGE(obs, mod) 0.45454545454545453

Source code in src/monet_stats/error_metrics.py

def NMdnGE(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Normalized Median Gross Error (%).

    Typical Use Cases
    -----------------
    - Comparing the typical (median) error magnitude, normalized by the mean
      observation, for robust model evaluation.
    - Useful for inter-comparison of model performance across sites or variables
      with different scales.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute the error.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Normalized median gross error (percent).

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import NMdnGE
    >>> obs = np.array([1, 2, 3, 4, 100])
    >>> mod = np.array([1.1, 2.1, 3.1, 4.1, 105])
    >>> NMdnGE(obs, mod)
    0.45454545454545453
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        result = (abs(mod - obs).median(dim=dim) / obs.mean(dim=dim)) * 100.0
        # Update history
        history = f"NMdnGE computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        return np.ma.masked_invalid(np.ma.median(np.ma.abs(mod - obs), axis=axis) / np.ma.mean(obs, axis=axis)) * 100.0

NO(obs, mod=None, axis=None)

N Observations (#).

Typical Use Cases

• Counting the number of valid (non-masked) observations in a dataset.
• Used to report sample size for statistical summaries and model evaluation.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray, optional Model predicted values (not used for NO but included for signature matching). axis : int, str, or iterable of such, optional Axis or dimension along which to count.

Returns

int, numpy.ndarray, or xarray.DataArray Number of valid observations.

Source code in src/monet_stats/error_metrics.py

def NO(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Optional[Union[np.ndarray, xr.DataArray]] = None,
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[int, np.ndarray, xr.DataArray]:
    """
    N Observations (#).

    Typical Use Cases
    -----------------
    - Counting the number of valid (non-masked) observations in a dataset.
    - Used to report sample size for statistical summaries and model evaluation.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray, optional
        Model predicted values (not used for NO but included for signature matching).
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to count.

    Returns
    -------
    int, numpy.ndarray, or xarray.DataArray
        Number of valid observations.
    """
    if isinstance(obs, xr.DataArray):
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        return obs.count(dim=dim)
    else:
        return (~np.ma.getmaskarray(obs)).sum(axis=axis)

NOP(obs, mod, axis=None)

N Observations/Prediction Pairs (#).

Typical Use Cases

• Counting the number of valid observation-prediction pairs for paired statistical analysis.
• Used to ensure sample size consistency in paired model evaluation metrics.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to count.

Returns

int, numpy.ndarray, or xarray.DataArray Number of valid pairs.

Source code in src/monet_stats/error_metrics.py

def NOP(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[int, np.ndarray, xr.DataArray]:
    """
    N Observations/Prediction Pairs (#).

    Typical Use Cases
    -----------------
    - Counting the number of valid observation-prediction pairs for paired
      statistical analysis.
    - Used to ensure sample size consistency in paired model evaluation metrics.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to count.

    Returns
    -------
    int, numpy.ndarray, or xarray.DataArray
        Number of valid pairs.
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        # count() on aligned xarray handles NaNs in both by alignment if NaNs are coords,
        # but if NaNs are in data, we need to mask.
        # However, count() counts non-NaN values.
        # To get pairs where BOTH are not NaN:
        mask = obs.notnull() & mod.notnull()
        return mask.sum(dim=dim)
    else:
        obsc, modc = matchmasks(obs, mod)
        return (~np.ma.getmaskarray(obsc)).sum(axis=axis)

NP(obs=None, mod=None, axis=None)

N Predictions (#).

Typical Use Cases

• Counting the number of valid (non-masked) model predictions in a dataset.
• Used to report sample size for model output and for filtering invalid predictions.

Parameters

obs : numpy.ndarray or xarray.DataArray, optional Observed values (not used for NP but included for signature matching). mod : numpy.ndarray or xarray.DataArray Model predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to count.

Returns

int, numpy.ndarray, or xarray.DataArray Number of valid predictions.

Source code in src/monet_stats/error_metrics.py

def NP(
    obs: Optional[Union[np.ndarray, xr.DataArray]] = None,
    mod: Union[np.ndarray, xr.DataArray] = None,
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[int, np.ndarray, xr.DataArray]:
    """
    N Predictions (#).

    Typical Use Cases
    -----------------
    - Counting the number of valid (non-masked) model predictions in a dataset.
    - Used to report sample size for model output and for filtering invalid
      predictions.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray, optional
        Observed values (not used for NP but included for signature matching).
    mod : numpy.ndarray or xarray.DataArray
        Model predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to count.

    Returns
    -------
    int, numpy.ndarray, or xarray.DataArray
        Number of valid predictions.
    """
    if isinstance(mod, xr.DataArray):
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = mod.dims[axis]
        else:
            dim = axis
        return mod.count(dim=dim)
    else:
        return (~np.ma.getmaskarray(mod)).sum(axis=axis)

NRMSE(obs, mod, axis=None)

Normalized Root Mean Square Error (NRMSE).

Typical Use Cases

• Quantifying the relative error between model and observations, normalized by the range of observations.
• Used in model evaluation to compare performance across different variables or sites with different scales.
• Provides dimensionless error metric for cross-comparison.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute NRMSE.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Normalized root mean square error (unitless).

Examples

import numpy as np from monet_stats.error_metrics import NRMSE obs = np.array([1, 2, 3, 4]) mod = np.array([2, 2, 2, 2]) NRMSE(obs, mod) 0.4714045207910317

Source code in src/monet_stats/error_metrics.py

def NRMSE(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Normalized Root Mean Square Error (NRMSE).

    Typical Use Cases
    -----------------
    - Quantifying the relative error between model and observations, normalized
      by the range of observations.
    - Used in model evaluation to compare performance across different variables
      or sites with different scales.
    - Provides dimensionless error metric for cross-comparison.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute NRMSE.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Normalized root mean square error (unitless).

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import NRMSE
    >>> obs = np.array([1, 2, 3, 4])
    >>> mod = np.array([2, 2, 2, 2])
    >>> NRMSE(obs, mod)
    0.4714045207910317
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        rmse = ((mod - obs) ** 2).mean(dim=dim, keep_attrs=True) ** 0.5
        obs_range = obs.max(dim=dim) - obs.min(dim=dim)
        result = xr.where(obs_range == 0, 0, rmse / obs_range)
        # Update history
        history = f"NRMSE computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        rmse = np.ma.sqrt(np.ma.mean((np.subtract(mod, obs)) ** 2, axis=axis))
        obs_range = np.ma.max(obs, axis=axis) - np.ma.min(obs, axis=axis)
        with np.errstate(divide="ignore", invalid="ignore"):
            result = np.where(obs_range == 0, 0, rmse / obs_range)
            return result.item() if np.ndim(result) == 0 else result

NSC(obs, mod, axis=None)

Nash-Sutcliffe Coefficient (NSC) - Alternative to NSE.

Typical Use Cases

• Quantifying the predictive power of hydrological models relative to the mean of observations.
• Used in hydrology, meteorology, and environmental model evaluation.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute NSC.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Nash-Sutcliffe coefficient (unitless).

Examples

import numpy as np from monet_stats.error_metrics import NSC obs = np.array([1, 2, 3, 4]) mod = np.array([2, 2, 2, 2]) NSC(obs, mod) -0.33333333333333326

Source code in src/monet_stats/error_metrics.py

def NSC(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Nash-Sutcliffe Coefficient (NSC) - Alternative to NSE.

    Typical Use Cases
    -----------------
    - Quantifying the predictive power of hydrological models relative to
      the mean of observations.
    - Used in hydrology, meteorology, and environmental model evaluation.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute NSC.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Nash-Sutcliffe coefficient (unitless).

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import NSC
    >>> obs = np.array([1, 2, 3, 4])
    >>> mod = np.array([2, 2, 2, 2])
    >>> NSC(obs, mod)
    -0.33333333333333326
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        obs_mean = obs.mean(dim=dim)
        numerator = ((obs - mod) ** 2).sum(dim=dim)
        denominator = ((obs - obs_mean) ** 2).sum(dim=dim)
        result = 1.0 - (numerator / denominator)
        # Update history
        history = f"NSC computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        obs_mean = np.mean(obs, axis=axis, keepdims=True)
        numerator = np.sum((obs - mod) ** 2, axis=axis)
        denominator = np.sum((obs - obs_mean) ** 2, axis=axis)
        return 1.0 - (numerator / denominator)

NSE_alpha(obs, mod, axis=None)

NSE Alpha - Decomposed NSE component measuring ratio of standard deviations.

Typical Use Cases

• Quantifying the model's ability to capture the variability of observations.
• Used in model evaluation to assess how well model represents observed variability.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute NSE_alpha.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray NSE alpha component (unitless).

Examples

import numpy as np from monet_stats.error_metrics import NSE_alpha obs = np.array([1, 2, 3, 4]) mod = np.array([2, 2, 2, 2]) NSE_alpha(obs, mod) 0.0

Source code in src/monet_stats/error_metrics.py

def NSE_alpha(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    NSE Alpha - Decomposed NSE component measuring ratio of standard deviations.

    Typical Use Cases
    -----------------
    - Quantifying the model's ability to capture the variability of observations.
    - Used in model evaluation to assess how well model represents observed
      variability.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute NSE_alpha.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        NSE alpha component (unitless).

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import NSE_alpha
    >>> obs = np.array([1, 2, 3, 4])
    >>> mod = np.array([2, 2, 2, 2])
    >>> NSE_alpha(obs, mod)
    0.0
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        result = mod.std(dim=dim) / obs.std(dim=dim)
        # Update history
        history = f"NSE_alpha computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        return np.std(mod, axis=axis) / np.std(obs, axis=axis)

NSE_beta(obs, mod, axis=None)

NSE Beta - Decomposed NSE component measuring bias.

Typical Use Cases

• Quantifying the systematic bias between model and observations.
• Used in model evaluation to assess mean differences between model and observations.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute NSE_beta.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray NSE beta component (unitless).

Examples

import numpy as np from monet_stats.error_metrics import NSE_beta obs = np.array([1, 2, 3, 4]) mod = np.array([2, 2, 2, 2]) NSE_beta(obs, mod) 0.5

Source code in src/monet_stats/error_metrics.py

def NSE_beta(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    NSE Beta - Decomposed NSE component measuring bias.

    Typical Use Cases
    -----------------
    - Quantifying the systematic bias between model and observations.
    - Used in model evaluation to assess mean differences between model and
      observations.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute NSE_beta.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        NSE beta component (unitless).

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import NSE_beta
    >>> obs = np.array([1, 2, 3, 4])
    >>> mod = np.array([2, 2, 2, 2])
    >>> NSE_beta(obs, mod)
    0.5
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        result = mod.mean(dim=dim) / obs.mean(dim=dim)
        # Update history
        history = f"NSE_beta computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        return np.mean(mod, axis=axis) / np.mean(obs, axis=axis)

RM(obs, mod, axis=None)

Root Mean Error (RM) - Root of mean squared error.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute the error.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Root of mean squared error (observation units). Returns 0.0 for perfect agreement.

Source code in src/monet_stats/error_metrics.py

def RM(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Root Mean Error (RM) - Root of mean squared error.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute the error.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Root of mean squared error (observation units).
        Returns 0.0 for perfect agreement.
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        result = np.sqrt(((obs - mod) ** 2).mean(dim=dim, keep_attrs=True))
        # Update history
        history = f"RM computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        return np.sqrt(np.mean((np.subtract(obs, mod)) ** 2, axis=axis))

RMSE(obs, mod, axis=None)

Root Mean Square Error (RMSE).

Typical Use Cases

• Quantifying the average magnitude of errors between model and observations, accounting for large errors more heavily than MAE.
• Used in model evaluation, forecast verification, and regression analysis.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute RMSE.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Root mean square error.

Examples

import numpy as np from monet_stats.error_metrics import RMSE obs = np.array([1, 2, 3]) mod = np.array([2, 2, 4]) RMSE(obs, mod) 0.816496580927726

Source code in src/monet_stats/error_metrics.py

def RMSE(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Root Mean Square Error (RMSE).

    Typical Use Cases
    -----------------
    - Quantifying the average magnitude of errors between model and observations,
      accounting for large errors more heavily than MAE.
    - Used in model evaluation, forecast verification, and regression analysis.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute RMSE.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Root mean square error.

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import RMSE
    >>> obs = np.array([1, 2, 3])
    >>> mod = np.array([2, 2, 4])
    >>> RMSE(obs, mod)
    0.816496580927726
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        result = ((mod - obs) ** 2).mean(dim=dim, keep_attrs=True) ** 0.5
        # Update history
        history = f"RMSE computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        return np.sqrt(np.mean((np.subtract(mod, obs)) ** 2, axis=axis))

RMSE_m(obs, mod, axis=None)

Root Mean Square Error (RMSE) - robust to masked arrays.

Typical Use Cases

• Quantifying the average magnitude of errors between model and observations, accounting for large errors more heavily than MAE, robust to missing data.
• Used in model evaluation, forecast verification, and regression analysis with incomplete datasets.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute RMSE.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Root mean square error.

Examples

import numpy as np from monet_stats.error_metrics import RMSE_m obs = np.array([1, 2, 3]) mod = np.array([2, 2, 4]) RMSE_m(obs, mod) 0.816496580927726

Source code in src/monet_stats/error_metrics.py

def RMSE_m(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Root Mean Square Error (RMSE) - robust to masked arrays.

    Typical Use Cases
    -----------------
    - Quantifying the average magnitude of errors between model and
      observations, accounting for large errors more heavily than MAE,
      robust to missing data.
    - Used in model evaluation, forecast verification, and regression analysis
      with incomplete datasets.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute RMSE.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Root mean square error.

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import RMSE_m
    >>> obs = np.array([1, 2, 3])
    >>> mod = np.array([2, 2, 4])
    >>> RMSE_m(obs, mod)
    0.816496580927726
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        # RMSE implementation for xarray already handles NaNs
        return RMSE(obs, mod, axis=axis)
    else:
        return np.ma.sqrt(np.ma.mean((np.subtract(mod, obs)) ** 2, axis=axis))

RMSE_norm(obs, mod, axis=None)

Normalized Root Mean Square Error (RMSE_norm).

Normalizes RMSE by the range of observations.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute normalized RMSE.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Normalized root mean square error (unitless).

Source code in src/monet_stats/error_metrics.py

def RMSE_norm(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Normalized Root Mean Square Error (RMSE_norm).

    Normalizes RMSE by the range of observations.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute normalized RMSE.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Normalized root mean square error (unitless).
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        rmse = ((mod - obs) ** 2).mean(dim=dim, keep_attrs=True) ** 0.5
        obs_min = obs.min(dim=dim)
        obs_max = obs.max(dim=dim)
        obs_range = obs_max - obs_min
        # Avoid division by zero
        result = xr.where(obs_range == 0, rmse, rmse / obs_range)
        # Update history
        history = f"RMSE_norm computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        rmse = np.sqrt(np.mean((np.subtract(mod, obs)) ** 2, axis=axis))
        obs_min = np.min(obs, axis=axis)
        obs_max = np.max(obs, axis=axis)
        obs_range = obs_max - obs_min
        # Avoid division by zero
        result = np.where(obs_range == 0, rmse, rmse / obs_range)
        return result.item() if np.ndim(result) == 0 else result

RMSPE(obs, mod, axis=None)

Root Mean Square Percentage Error (RMSPE).

Typical Use Cases

• Quantifying the average relative error between model and observations as a percentage, emphasizing larger errors.
• Used in time series forecasting, regression, and model evaluation for percentage-based error assessment.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute RMSPE.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Root mean square percentage error (in percent).

Examples

import numpy as np from monet_stats.error_metrics import RMSPE obs = np.array([1, 2, 3]) mod = np.array([2, 2, 4]) RMSPE(obs, mod) 50.0

Source code in src/monet_stats/error_metrics.py

def RMSPE(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Root Mean Square Percentage Error (RMSPE).

    Typical Use Cases
    -----------------
    - Quantifying the average relative error between model and observations as
      a percentage, emphasizing larger errors.
    - Used in time series forecasting, regression, and model evaluation for
      percentage-based error assessment.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute RMSPE.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Root mean square percentage error (in percent).

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import RMSPE
    >>> obs = np.array([1, 2, 3])
    >>> mod = np.array([2, 2, 4])
    >>> RMSPE(obs, mod)
    50.0
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        result = (100 * ((mod - obs) / obs) ** 2).mean(dim=dim, keep_attrs=True) ** 0.5
        # Update history
        history = f"RMSPE computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        return 100 * np.ma.sqrt(np.ma.mean(((mod - obs) / obs) ** 2, axis=axis))

RMdn(obs, mod, axis=None)

Root Median Error (RMdn) - Root of median squared error.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute the error.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Root of median squared error (observation units). Returns 0.0 for perfect agreement.

Source code in src/monet_stats/error_metrics.py

def RMdn(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Root Median Error (RMdn) - Root of median squared error.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute the error.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Root of median squared error (observation units).
        Returns 0.0 for perfect agreement.
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        result = np.sqrt(((obs - mod) ** 2).median(dim=dim, keep_attrs=True))
        # Update history
        history = f"RMdn computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        squared_errors = (np.subtract(obs, mod)) ** 2
        return np.sqrt(np.median(squared_errors, axis=axis))

STDO(obs, mod, axis=None)

Standard deviation of Observation Errors (obs - mod).

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute the standard deviation. If None, computes over all axes.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Standard deviation of (observation - model) errors. Returns 0.0 for perfect agreement.

Examples

import numpy as np from monet_stats.error_metrics import STDO obs = np.array([1.0, 2.0, 3.0]) mod = np.array([1.1, 1.9, 3.2]) STDO(obs, mod) 0.1247219128924647

Source code in src/monet_stats/error_metrics.py

def STDO(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Standard deviation of Observation Errors (obs - mod).

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute the standard deviation.
        If None, computes over all axes.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Standard deviation of (observation - model) errors.
        Returns 0.0 for perfect agreement.

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import STDO
    >>> obs = np.array([1.0, 2.0, 3.0])
    >>> mod = np.array([1.1, 1.9, 3.2])
    >>> STDO(obs, mod)
    0.1247219128924647
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        errors = obs - mod
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        result = errors.std(dim=dim, keep_attrs=True)
        # Update history
        history = f"STDO computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result

    # Fallback to numpy-compatible logic
    errors = np.subtract(obs, mod)
    return np.std(errors, axis=axis)

STDP(obs, mod, axis=None)

Standard deviation of Prediction Errors (mod - obs).

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute the standard deviation. If None, computes over all axes.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Standard deviation of (model - observation) errors. Returns 0.0 for perfect agreement.

Examples

import numpy as np from monet_stats.error_metrics import STDP obs = np.array([1.0, 2.0, 3.0]) mod = np.array([1.1, 1.9, 3.2]) STDP(obs, mod) 0.1247219128924647

Source code in src/monet_stats/error_metrics.py

def STDP(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Standard deviation of Prediction Errors (mod - obs).

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute the standard deviation.
        If None, computes over all axes.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Standard deviation of (model - observation) errors.
        Returns 0.0 for perfect agreement.

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import STDP
    >>> obs = np.array([1.0, 2.0, 3.0])
    >>> mod = np.array([1.1, 1.9, 3.2])
    >>> STDP(obs, mod)
    0.1247219128924647
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        errors = mod - obs
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        result = errors.std(dim=dim, keep_attrs=True)
        # Update history
        history = f"STDP computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result

    # Fallback to numpy-compatible logic
    errors = np.subtract(mod, obs)
    return np.std(errors, axis=axis)

VOLUMETRIC_ERROR(obs, mod, axis=None)

Volumetric Error Metric.

Typical Use Cases

• Quantifying the volume difference between observed and modeled features.
• Used in hydrology for flood extent verification.
• Applied in meteorology for precipitation volume verification.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute volumetric error.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Volumetric error metric.

Examples

import numpy as np from monet_stats.error_metrics import VOLUMETRIC_ERROR obs = np.array([1, 2, 3]) mod = np.array([2, 2, 4]) VOLUMETRIC_ERROR(obs, mod) 0.2

Source code in src/monet_stats/error_metrics.py

def VOLUMETRIC_ERROR(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Volumetric Error Metric.

    Typical Use Cases
    -----------------
    - Quantifying the volume difference between observed and modeled features.
    - Used in hydrology for flood extent verification.
    - Applied in meteorology for precipitation volume verification.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute volumetric error.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Volumetric error metric.

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import VOLUMETRIC_ERROR
    >>> obs = np.array([1, 2, 3])
    >>> mod = np.array([2, 2, 4])
    >>> VOLUMETRIC_ERROR(obs, mod)
    0.2
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        obs_sum = obs.sum(dim=dim)
        mod_sum = mod.sum(dim=dim)
        result = abs(mod_sum - obs_sum) / abs(obs_sum)
        # Update history
        history = f"VOLUMETRIC_ERROR computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        obs_sum = np.sum(obs, axis=axis)
        mod_sum = np.sum(mod, axis=axis)
        return np.abs(mod_sum - obs_sum) / np.abs(obs_sum)

WDMB(obs, mod, axis=None)

Wind Direction Mean Bias (WDMB, standard version).

This version uses circlebias, which is not robust to masked arrays. Use this if your data are dense and do not contain missing values.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed wind direction values (degrees). mod : numpy.ndarray or xarray.DataArray Model predicted wind direction values (degrees). axis : int, str, or iterable of such, optional Axis or dimension along which to compute the mean bias.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Mean wind direction bias (degrees).

Source code in src/monet_stats/error_metrics.py

def WDMB(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Wind Direction Mean Bias (WDMB, standard version).

    This version uses circlebias, which is not robust to masked arrays.
    Use this if your data are dense and do not contain missing values.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed wind direction values (degrees).
    mod : numpy.ndarray or xarray.DataArray
        Model predicted wind direction values (degrees).
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute the mean bias.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Mean wind direction bias (degrees).
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        result = circlebias(mod - obs).mean(dim=dim, keep_attrs=True)
        # Update history
        history = f"WDMB computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        return np.ma.mean(circlebias(np.subtract(mod, obs)), axis=axis)

WDMB_m(obs, mod, axis=None)

Wind Direction Mean Bias (WDMB, robust version for masked arrays).

This version uses circlebias_m, which is robust to masked arrays and missing data. Use this if your data may contain NaNs or masked values.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed wind direction values (degrees). mod : numpy.ndarray or xarray.DataArray Model predicted wind direction values (degrees). axis : int, str, or iterable of such, optional Axis or dimension along which to compute the mean bias.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Mean wind direction bias (degrees).

Source code in src/monet_stats/error_metrics.py

def WDMB_m(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Wind Direction Mean Bias (WDMB, robust version for masked arrays).

    This version uses circlebias_m, which is robust to masked arrays and
    missing data. Use this if your data may contain NaNs or masked values.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed wind direction values (degrees).
    mod : numpy.ndarray or xarray.DataArray
        Model predicted wind direction values (degrees).
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute the mean bias.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Mean wind direction bias (degrees).
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        result = circlebias_m(mod - obs).mean(dim=dim, keep_attrs=True)
        # Update history
        history = f"WDMB_m computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        return np.ma.mean(circlebias_m(np.subtract(mod, obs)), axis=axis)

WDMdnB(obs, mod, axis=None)

Wind Direction Median Bias (WDMdnB).

Parameters

obs : numpy.ndarray or xarray.DataArray Observed wind direction values (degrees). mod : numpy.ndarray or xarray.DataArray Model predicted wind direction values (degrees). axis : int, str, or iterable of such, optional Axis or dimension along which to compute the median bias.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Median wind direction bias (degrees).

Source code in src/monet_stats/error_metrics.py

def WDMdnB(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Wind Direction Median Bias (WDMdnB).

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed wind direction values (degrees).
    mod : numpy.ndarray or xarray.DataArray
        Model predicted wind direction values (degrees).
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute the median bias.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Median wind direction bias (degrees).
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        result = circlebias(mod - obs).median(dim=dim, keep_attrs=True)
        # Update history
        history = f"WDMdnB computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        return np.ma.median(circlebias(np.subtract(mod, obs)), axis=axis)

bias_fraction(obs, mod, axis=None)

Bias Fraction (BF).

Quantifies the fraction of total error that is due to systematic bias.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute bias fraction.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Bias fraction (unitless, 0-1).

Source code in src/monet_stats/error_metrics.py

def bias_fraction(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Bias Fraction (BF).

    Quantifies the fraction of total error that is due to systematic bias.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute bias fraction.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Bias fraction (unitless, 0-1).
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        bias = (mod - obs).mean(dim=dim)
        total_error = np.sqrt(((mod - obs) ** 2).mean(dim=dim, keep_attrs=True))
        # Avoid division by zero
        result = xr.where(total_error == 0, 0, (bias**2) / (total_error**2))
        # Update history
        history = f"bias_fraction computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        bias = np.mean(np.subtract(mod, obs), axis=axis)
        total_error = np.sqrt(np.mean((np.subtract(mod, obs)) ** 2, axis=axis))
        # Avoid division by zero
        result = np.where(total_error == 0, 0, (bias**2) / (total_error**2))
        return result.item() if np.ndim(result) == 0 else result

sMAPE(obs, mod, axis=None)

Symmetric Mean Absolute Percentage Error (sMAPE).

Typical Use Cases

• Quantifying the average relative error between model and observations, normalized by their mean.
• Used in time series forecasting, regression, and model evaluation for percentage-based error assessment.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute sMAPE.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Symmetric mean absolute percentage error (in percent).

Examples

import numpy as np from monet_stats.error_metrics import sMAPE obs = np.array([1, 2, 3]) mod = np.array([2, 2, 4]) sMAPE(obs, mod) 28.57142857142857

Source code in src/monet_stats/error_metrics.py

def sMAPE(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Symmetric Mean Absolute Percentage Error (sMAPE).

    Typical Use Cases
    -----------------
    - Quantifying the average relative error between model and observations,
      normalized by their mean.
    - Used in time series forecasting, regression, and model evaluation for
      percentage-based error assessment.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute sMAPE.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Symmetric mean absolute percentage error (in percent).

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import sMAPE
    >>> obs = np.array([1, 2, 3])
    >>> mod = np.array([2, 2, 4])
    >>> sMAPE(obs, mod)
    28.57142857142857
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        obs, mod = xr.align(obs, mod, join="inner")
        # Handle axis vs dim
        if axis is not None and isinstance(axis, int):
            dim = obs.dims[axis]
        else:
            dim = axis
        result = (200 * abs(mod - obs) / (abs(mod) + abs(obs))).mean(dim=dim, keep_attrs=True)
        # Update history
        history = f"sMAPE computed at {pd.Timestamp.now().isoformat()}"
        result.attrs["history"] = f"{result.attrs.get('history', '')}\n{history}".strip()
        return result
    else:
        return (200 * np.ma.abs(np.subtract(mod, obs)) / (np.ma.abs(mod) + np.ma.abs(obs))).mean(axis=axis)

sMAPEm(obs, mod, axis=None)

Symmetric Mean Absolute Percentage Error (sMAPE) - robust to masked arrays.

Typical Use Cases

• Quantifying the average relative error between model and observations, normalized by their mean, robust to missing data.
• Used in time series forecasting, regression, and model evaluation for percentage-based error assessment.

Parameters

obs : numpy.ndarray or xarray.DataArray Observed values. mod : numpy.ndarray or xarray.DataArray Model or predicted values. axis : int, str, or iterable of such, optional Axis or dimension along which to compute sMAPE.

Returns

numpy.number, numpy.ndarray, or xarray.DataArray Symmetric mean absolute percentage error (in percent).

Examples

import numpy as np from monet_stats.error_metrics import sMAPEm obs = np.array([1, 2, 3]) mod = np.array([2, 2, 4]) sMAPEm(obs, mod) 28.57142857142857

Source code in src/monet_stats/error_metrics.py

def sMAPEm(
    obs: Union[np.ndarray, xr.DataArray],
    mod: Union[np.ndarray, xr.DataArray],
    axis: Optional[Union[int, str, Iterable[Union[int, str]]]] = None,
) -> Union[np.number, np.ndarray, xr.DataArray]:
    """
    Symmetric Mean Absolute Percentage Error (sMAPE) - robust to masked arrays.

    Typical Use Cases
    -----------------
    - Quantifying the average relative error between model and observations,
      normalized by their mean, robust to missing data.
    - Used in time series forecasting, regression, and model evaluation for
      percentage-based error assessment.

    Parameters
    ----------
    obs : numpy.ndarray or xarray.DataArray
        Observed values.
    mod : numpy.ndarray or xarray.DataArray
        Model or predicted values.
    axis : int, str, or iterable of such, optional
        Axis or dimension along which to compute sMAPE.

    Returns
    -------
    numpy.number, numpy.ndarray, or xarray.DataArray
        Symmetric mean absolute percentage error (in percent).

    Examples
    --------
    >>> import numpy as np
    >>> from monet_stats.error_metrics import sMAPEm
    >>> obs = np.array([1, 2, 3])
    >>> mod = np.array([2, 2, 4])
    >>> sMAPEm(obs, mod)
    28.57142857142857
    """
    if isinstance(obs, xr.DataArray) and isinstance(mod, xr.DataArray):
        # sMAPE implementation for xarray already handles NaNs
        return sMAPE(obs, mod, axis=axis)
    else:
        return 200 * np.ma.mean(np.ma.abs(mod - obs) / (np.ma.abs(mod) + np.ma.abs(obs)), axis=axis)