The Continuous Ranked Probability Score (CRPS) is given by: \begin{equation} \mathrm{CRPS}(F, x) = \int_{-\infty}^{\infty} \left( F(y) - \mathbb{1}(y - x) \right)^2 \, dy \end{equation}

I am trying to intuitively understand, in steps, how a CRPS of 0.625 is obtained by the properscoring Python library using ensemble values of [1, 2, 3, 4] and an observed value of 3.5. Another python package (CRPS) also gives 0.625, so I am confident that 0.625 is the result I want.

Python code that gives 0.625

import numpy as np
import properscoring as ps

ensemble_values = np.array([1, 2, 3, 4])
observed_value = 3.5

ps_result = ps.crps_ensemble(observed_value, ensemble_values)
print(ps_result)  # 0.625

Step by step calculation that gives 0.875

The CDF of the ensemble is:
F(1) = 0.25, F(2) = 0.5, F(3) = 0.75, F(4) = 1.0

The observed value is 3.5, so the CDF of the observed value H(x) is a step function that is 0 for all values less than 3.5 and 1 for all values greater than or equal to 3.5.

Therefore:
H(1) = H(2) = H(3) = 0 and H(4) = 1

We sum the squared difference at each ensemble point:

For x1: (0.25 - 0)**2 = 0.0625
For x2: (0.5 - 0)**2 = 0.25
For x3: (0.75 - 0)**2 = 0.5625
For x4: (1.0 - 1)**2 = 0

CRPS = 0.0625 + 0.25 + 0.5625 + 0 = 0.875

Any idea what the missing piece of the puzzle is?