xskillscore.brier_score¶
-
xskillscore.
brier_score
(observations, forecasts, member_dim='member', fair=False, dim=None, weights=None, keep_attrs=False)¶ Calculate Brier score (BS).
- Parameters
observations (xarray.Dataset or xarray.DataArray) – The observations or set of observations of the event. Data should be boolean or logical (True or 1 for event occurance, False or 0 for non-occurance).
forecasts (xarray.Dataset or xarray.DataArray) – The forecast likelihoods of the event. If
fair==False
, forecasts should be between 0 and 1 without a dimensionmember_dim
or should be boolean (True,False) or binary (0, 1) containing a member dimension (probabilities will be internally calculated by.mean(member_dim))
. Iffair==True
, forecasts must be boolean (True,False) or binary (0, 1) containing dimensionmember_dim
.member_dim (str, optional) – Name of ensemble member dimension. By default, ‘member’.
fair (boolean) – Apply ensemble member-size adjustment for unbiased, fair metric; see Ferro (2013). Defaults to False.
dim (str or list of str, optional) – Dimension over which to compute mean after computing
brier_score
. Defaults to None implying averaging over all dimensions.weights (xr.DataArray with dimensions from dim, optional) – Weights for weighted.mean(dim). Defaults to None, such that no weighting is applied.
keep_attrs (bool) – If True, the attributes (attrs) will be copied from the first input to the new one. If False (default), the new object will be returned without attributes.
- Returns
- Return type
Examples
>>> observations = xr.DataArray(np.random.normal(size=(3, 3)), ... coords=[('x', np.arange(3)), ... ('y', np.arange(3))]) >>> forecasts = xr.DataArray(np.random.normal(size=(3, 3, 3)), ... coords=[('x', np.arange(3)), ... ('y', np.arange(3)), ... ('member', np.arange(3))]) >>> xs.brier_score(observations > .5, ... (forecasts > .5).mean('member'), ... dim="y") <xarray.DataArray (x: 3)> array([0.37037037, 0.14814815, 0.51851852]) Coordinates: * x (x) int64 0 1 2
See also
properscoring.brier_score
References
Gneiting, Tilmann, and Adrian E Raftery. “Strictly Proper Scoring Rules, Prediction, and Estimation.” Journal of the American Statistical Association 102, no. 477 (March 1, 2007): 359–78. https://doi.org/10/c6758w.
Brier, Glenn W. “VERIFICATION OF FORECASTS EXPRESSED IN TERMS OF PROBABILITY.” Monthly Weather Review, 78(1): 1-3 https://journals.ametsoc.org/doi/abs/10.1175/1520-0493%281950%29078%3C0001%3AVOFEIT%3E2.0.CO%3B2
C. A. T. Ferro. Fair scores for ensemble forecasts. Q.R.J. Meteorol. Soc., 140: 1917–1923, 2013. doi: 10.1002/qj.2270.