weighted_interval_score adds w_median * median instead of w_median * |obs - median| on all non-numba backends

Abierto
  • ¿Algo no está bien?
#140 0 comentarios 0 reacciones 0 asignados Ver en GitHub

Los mantenedores suelen responder en 1 día

@lucaluo925 ya está trabajando en esto.

Desde el 29/9/2026.

  • #141 de @lucaluo925 — abierto

Evaluación

Dificultad
2/5
Tiempo estimado
1-3 horas
Aptitud para principiantes
35/100
Tipo de issue
Error
Claridad
Bien especificado
Estado de actividad
Estancado
Stack tecnológico
numpy, python
Área
data

Línea de trabajo

The issue identifies scoringrules/core/interval/_score.py and tests/test_interval.py. Read the weighted interval score implementation and its existing test, then use a median different from the observation in a regression case. Done means the non-numba backends match the stated formula and the regression test catches the previous result; a linked open pull request is already addressing this issue.

Escrito por el modelo de indexación a partir del texto del issue.

Descripción

bug
What happened?

weighted_interval_score returns incorrect values under the numpy backend. The numba backend produces correct results, meaning the same function yields different outputs across backends.

Line 39 in scoringrules/core/interval/_score.py reads WIS += w_median * median. This adds the median itself instead of the absolute deviation between the observation and the median. Per equation (3) in Bracher et al. (2021), this term should be w_median * |obs - median|. The numba implementation in core/interval/_gufunc.py correctly implements weight_median * np.abs(obs - median).

Since _score.py is the shared implementation for all non-numba backends, numpy, jax and torch are all affected. I have verified numpy and jax: both return 0.458, while numba returns 0.338.

The issue is more than a numerical discrepancy: WIS is a weighted average of non-negative terms and mathematically cannot be negative. However, the numpy backend returns negative WIS values when the median is negative. The third example in my reproduction script returns -19.82.

Why existing tests did not catch this: test_weighted_interval_score in tests/test_interval.py calls sr.weighted_interval_score(obs, obs, lower, upper, alpha), passing obs as both the observation and the median. In this case median = obs and |obs - median| = 0, so both formulas evaluate to zero and the faulty term vanishes. As a result, CI passes for all four backends. This also explains why the issue persisted even after #81 and #82 improved backend test coverage.

Fix: change line 39 to WIS += w_median * B.abs(obs - median) (using the absolute value function of the active backend) to match the numba implementation. For tests, it is recommended to use a median value different from obs; otherwise the bug will remain untestable after the fix.

There is also a mismatch between documentation and code: the docstring states the default weights are w_k = 2/alpha_k, but the code uses alpha/2. The code is correct (consistent with the paper), so the docstring needs updating.

I also compared all 97 exported functions one-by-one between numpy and numba backends, and only this function shows inconsistency. I independently verified all closed-form CRPS variants, energy score, variogram score, interval score and pinball score via numerical integration and reference paper formulas, and all are correct. This is therefore an isolated bug.

If this assessment holds, I can open a PR to fix this single line and add a regression test covering cases where median != obs.

AI note: This bug was discovered with assistance from an AI assistant. The reproduction script and validations above were run by the assistant, and I have manually verified every line prior to submission.

Minimal reproducer
import numpy as np, scoringrules as sr

# Bracher et al. (2021), Eq. (3):
#   WIS = (1/(K+1/2)) * ( w0*|y - m| + sum_k w_k * IS_k ),  w0 = 1/2, w_k = alpha_k/2
def wis_reference(y, m, lower, upper, alpha):
    IS = (upper - lower) \
         + 2/alpha*(lower - y)*(y < lower) \
         + 2/alpha*(y - upper)*(y > upper)
    return (0.5*abs(y - m) + np.sum(alpha/2 * IS)) / (len(alpha) + 0.5)

cases = {
    "median = obs (the case the test suite uses)":
        (0.0, 0.0, np.array([-1.2, -0.7]), np.array([1.2, 0.7])),
    "median != obs":
        (1.4, 1.0, np.array([0.5, 0.8]), np.array([1.2, 1.1])),
    "negative location":
        (-100.0, -100.0, np.array([-101.0, -100.5]), np.array([-99.0, -99.5])),
}
alpha = np.array([0.2, 0.5])
for name, (y, m, lo, hi) in cases.items():
    ref = wis_reference(y, m, lo, hi, alpha)
    npv = float(sr.weighted_interval_score(y, m, lo, hi, alpha, backend="numpy"))
    nbv = float(sr.weighted_interval_score(y, m, lo, hi, alpha, backend="numba"))
    print(f"{name}\n  reference={ref:12.6f}  numpy={npv:12.6f}  numba={nbv:12.6f}")

# Output on scoringrules 0.10.0 / numpy 2.4.6 / numba 0.67.0 / Python 3.11:
#
# median = obs (the case the test suite uses)
#   reference=    0.236000  numpy=    0.236000  numba=    0.236000
# median != obs
#   reference=    0.338000  numpy=    0.458000  numba=    0.338000
# negative location
#   reference=    0.180000  numpy=  -19.820000  numba=    0.180000
Backend

numpy, jax, torch

Python version

3.11

scoringrules version

0.10.0

Lenguaje dominante
Python
Estrellas
108
Forks
15
Merge medio
1 h 42 min
PR fusionados (30 d)
9

Preparar el entorno

Primeros pasos

  1. Lee el issue completo y luego la guía de contribución del proyecto.
  2. Comenta en el issue que vas a ocuparte — evita que dos personas hagan lo mismo.
  3. Haz un fork del repositorio y trabaja en una rama.
  4. Abre un pull request que haga referencia al número del issue.

Más de frazane/scoringrules

Todos los issues de frazane/scoringrules

Issues similares

Más issues de Python

Recibe los nuevos issues en tu correo

Un resumen breve de issues de GitHub para principiantes.