Traditional Binary Forecast Verification

This example demonstrates how to perform traditional binary forecast verification

Out:

Hits                              15.000000
Misses                             7.000000
False Alarm                        6.000000
Correct Negative                   8.000000
Probability of Detection           0.681818
False Alarm Ratio                  0.285714
Critical Success Index             0.535714
Frequency Bias                     0.954545
Accuracy                           0.638889
Probability of False Detection     0.428571
Heidke Skill Score                 0.250000
Equitable Threat Score             0.142857
dtype: float64

import numpy as np
import pandas as pd
from swirlspy.ver.crosstab import contingency
import swirlspy.ver.metric as mt
import xarray as xr


# Generating sample forecast and observed data
# Dimensions: x, y, time (3 x 3 x 4)
forecast_data = [
    [[7, 4, 6, 5], [4, 7, 9, 2], [3, 1, 6, 5]],
    [[4, 6, 8, 1], [11, 0, 3, 6], [0, 3, 6, 8]],
    [[5, 3, 7, 5], [7, 7, 6, 2], [3, 8, 9, 2]]
]

observed_data = [
    [[8, 7, 11, 4], [6, 7, 6, 0], [3, 3, 7, 3]],
    [[8, 8, 9, 5], [12, 1, 2, 8], [1, 4, 3, 12]],
    [[5, 6, 6, 2], [3, 6, 4, 7], [5, 7, 8, 2]]
]

# Creating xarrays
x = np.arange(3)
y = np.arange(3)
time = pd.date_range('1/1/2011', periods=4, freq='D')
forecast = xr.DataArray(
        forecast_data,
        coords=[('x', x), ('y', y), ('time', time)])
observed = xr.DataArray(
        observed_data,
        coords=[('x', x), ('y', y), ('time', time)])

# Applying the contingency function
cont = contingency(5, forecast, observed)

# Generating skill metrics
pod = mt.pod(cont)
far = mt.far(cont)
csi = mt.csi(cont)
freq_bias = mt.freq_bias(cont)
accuracy = mt.accuracy(cont)
pofd = mt.pofd(cont)
hss = mt.hss(cont)
ets = mt.ets(cont)

# Displaying results as pandas series
d = {
    'Hits': cont[0],
    'Misses': cont[1],
    'False Alarm': cont[2],
    'Correct Negative': cont[3],
    'Probability of Detection': pod,
    'False Alarm Ratio': far,
    'Critical Success Index': csi,
    'Frequency Bias': freq_bias,
    'Accuracy': accuracy,
    'Probability of False Detection': pofd,
    'Heidke Skill Score': hss,
    'Equitable Threat Score': ets
    }

pd_series_d = pd.Series(d)
print(pd_series_d)

Total running time of the script: ( 0 minutes 2.246 seconds)