From 2c07a21a7607ba4f5200bd1ae5b06aa2509b810b Mon Sep 17 00:00:00 2001
From: Thibault Hallouin <thibault.hallouin@inrae.fr>
Date: Fri, 27 Jan 2023 09:49:52 +0100
Subject: [PATCH] add unit tests for intervals-based metrics
---
tests/test_probabilist.cpp | 90 ++++++++++++++++++++++++++++++++------
1 file changed, 76 insertions(+), 14 deletions(-)
diff --git a/tests/test_probabilist.cpp b/tests/test_probabilist.cpp
index e25233e..1e15bf9 100644
--- a/tests/test_probabilist.cpp
+++ b/tests/test_probabilist.cpp
@@ -29,7 +29,8 @@ std::vector<std::string> all_metrics_p = {
"BS", "BSS", "BS_CRD", "BS_LBD",
"QS", "CRPS",
"POD", "POFD", "FAR", "CSI", "ROCSS",
- "RANK_DIAG", "DS", "AS"
+ "RANK_DIAG", "DS", "AS",
+ "CR", "AW", "AWN", "AWI", "WS", "WSS"
};
std::tuple<xt::xtensor<double, 1>, xt::xtensor<double, 2>> load_data_p()
@@ -406,7 +407,6 @@ TEST(ProbabilistTests, TestContingency)
);
}
-
TEST(ProbabilistTests, TestRanks)
{
// read in data
@@ -424,6 +424,7 @@ TEST(ProbabilistTests, TestRanks)
{"RANK_DIAG", "DS", "AS"},
xt::xtensor<double, 2>({}),
"high", // events
+ {}, // c_lvl
xt::xtensor<bool, 4>({}), // t_msk
{}, // m_cdt
xtl::missing<const std::unordered_map<std::string, int>>(), // bootstrap
@@ -497,6 +498,52 @@ TEST(ProbabilistTests, TestRanks)
);
}
+TEST(ProbabilistTests, TestIntervals)
+{
+ // read in data
+ xt::xtensor<double, 1> observed;
+ xt::xtensor<double, 2> predicted;
+ std::tie(observed, predicted) = load_data_p();
+
+ // compute scores
+ std::vector<xt::xarray<double>> metrics =
+ evalhyd::evalp(
+ // shape: (sites [1], time [t])
+ xt::eval(xt::view(observed, xt::newaxis(), xt::all())),
+ // shape: (sites [1], lead times [1], members [m], time [t])
+ xt::eval(xt::view(predicted, xt::newaxis(), xt::newaxis(), xt::all(), xt::all())),
+ {"CR", "AW", "AWN", "AWI", "WS", "WSS"},
+ xt::xtensor<double, 2>({}),
+ "", // events
+ {30., 80.} // c_lvl
+ );
+
+ // check results
+ // coverage ratios
+ xt::xtensor<double, 5> cr = {{{{{ 0.006431, 0.03537 }}}}};
+ EXPECT_TRUE(xt::all(xt::isclose(metrics[0], cr, 1e-05, 1e-06, true)));
+
+ // average widths
+ xt::xtensor<double, 5> aw = {{{{{ 9.27492, 31.321543}}}}};
+ EXPECT_TRUE(xt::all(xt::isclose(metrics[1], aw, 1e-05, 1e-06, true)));
+
+ // average widths normalised
+ xt::xtensor<double, 5> awn = {{{{{ 0.007383, 0.024931}}}}};
+ EXPECT_TRUE(xt::all(xt::isclose(metrics[2], awn, 1e-05, 1e-06, true)));
+
+ // average widths indices
+ xt::xtensor<double, 5> awi = {{{{{ 0.982112, 0.988095}}}}};
+ EXPECT_TRUE(xt::all(xt::isclose(metrics[3], awi, 1e-05, 1e-06, true)));
+
+ // Winkler scores
+ xt::xtensor<double, 5> ws = {{{{{ 764.447175, 2578.138264}}}}};
+ EXPECT_TRUE(xt::all(xt::isclose(metrics[4], ws, 1e-05, 1e-06, true)));
+
+ // Winkler skill scores
+ xt::xtensor<double, 5> wss = {{{{{ 0.662189, 0.436039}}}}};
+ EXPECT_TRUE(xt::all(xt::isclose(metrics[5], wss, 1e-05, 1e-06, true)));
+}
+
TEST(ProbabilistTests, TestMasks)
{
// read in data
@@ -512,6 +559,7 @@ TEST(ProbabilistTests, TestMasks)
// compute scores using masks to subset whole record
xt::xtensor<double, 2> thresholds = {{690, 534, 445}};
+ std::vector<double> confidence_levels = {30., 80.};
std::vector<xt::xarray<double>> metrics_masked =
evalhyd::evalp(
@@ -522,6 +570,7 @@ TEST(ProbabilistTests, TestMasks)
all_metrics_p,
thresholds,
"high",
+ confidence_levels,
// shape: (sites [1], lead times [1], subsets [1], time [t])
masks
);
@@ -535,7 +584,8 @@ TEST(ProbabilistTests, TestMasks)
xt::eval(xt::view(predicted, xt::newaxis(), xt::newaxis(), xt::all(), xt::range(20, _))),
all_metrics_p,
thresholds,
- "high"
+ "high",
+ confidence_levels
);
// check results are identical
@@ -563,6 +613,7 @@ TEST(ProbabilistTests, TestMasks)
TEST(ProbabilistTests, TestMaskingConditions)
{
xt::xtensor<double, 2> thresholds = {{690, 534, 445}};
+ std::vector<double> confidence_levels = {30., 80.};
// read in data
xt::xtensor<double, 1> observed_;
@@ -589,6 +640,7 @@ TEST(ProbabilistTests, TestMaskingConditions)
all_metrics_p,
thresholds,
"high",
+ confidence_levels,
masks,
q_conditions
);
@@ -600,7 +652,8 @@ TEST(ProbabilistTests, TestMaskingConditions)
xt::eval(xt::view(predicted, xt::newaxis(), xt::newaxis(), xt::all(), xt::all())),
all_metrics_p,
thresholds,
- "high"
+ "high",
+ confidence_levels
);
// check results are identical
@@ -644,6 +697,7 @@ TEST(ProbabilistTests, TestMaskingConditions)
all_metrics_p,
thresholds,
"high",
+ confidence_levels,
masks,
q_conditions_
);
@@ -655,7 +709,8 @@ TEST(ProbabilistTests, TestMaskingConditions)
xt::eval(xt::view(predicted, xt::newaxis(), xt::newaxis(), xt::all(), xt::all())),
all_metrics_p,
thresholds,
- "high"
+ "high",
+ confidence_levels
);
// check results are identical
@@ -696,6 +751,7 @@ TEST(ProbabilistTests, TestMaskingConditions)
all_metrics_p,
thresholds,
"high",
+ confidence_levels,
masks,
t_conditions
);
@@ -707,7 +763,8 @@ TEST(ProbabilistTests, TestMaskingConditions)
xt::eval(xt::view(predicted, xt::newaxis(), xt::newaxis(), xt::all(), xt::range(0, 100))),
all_metrics_p,
thresholds,
- "high"
+ "high",
+ confidence_levels
);
// check results are identical
@@ -737,8 +794,8 @@ TEST(ProbabilistTests, TestMaskingConditions)
TEST(ProbabilistTests, TestMissingData)
{
- xt::xtensor<double, 2> thresholds
- {{ 4., 5. }};
+ xt::xtensor<double, 2> thresholds = {{ 4., 5. }};
+ std::vector<double> confidence_levels = {30., 80.};
// compute metrics on series with NaN
xt::xtensor<double, 4> forecast_nan {{
@@ -761,7 +818,8 @@ TEST(ProbabilistTests, TestMissingData)
forecast_nan,
all_metrics_p,
thresholds,
- "high"
+ "high",
+ confidence_levels
);
// compute metrics on manually subset series (one leadtime at a time)
@@ -781,7 +839,8 @@ TEST(ProbabilistTests, TestMissingData)
forecast_pp1,
all_metrics_p,
thresholds,
- "high"
+ "high",
+ confidence_levels
);
xt::xtensor<double, 4> forecast_pp2 {{
@@ -800,7 +859,8 @@ TEST(ProbabilistTests, TestMissingData)
forecast_pp2,
all_metrics_p,
thresholds,
- "high"
+ "high",
+ confidence_levels
);
// check that numerical results are identical
@@ -824,8 +884,8 @@ TEST(ProbabilistTests, TestMissingData)
TEST(ProbabilistTests, TestBootstrap)
{
- xt::xtensor<double, 2> thresholds
- {{ 33.87, 55.67 }};
+ xt::xtensor<double, 2> thresholds = {{ 33.87, 55.67 }};
+ std::vector<double> confidence_levels = {30., 80.};
// read in data
std::ifstream ifs;
@@ -854,6 +914,7 @@ TEST(ProbabilistTests, TestBootstrap)
all_metrics_p,
thresholds,
"high", // events
+ confidence_levels,
xt::xtensor<bool, 4>({}), // t_msk
{}, // m_cdt
bootstrap,
@@ -877,7 +938,8 @@ TEST(ProbabilistTests, TestBootstrap)
xt::eval(xt::view(predicted_x3, xt::newaxis(), xt::newaxis(), xt::all(), xt::all())),
all_metrics_p,
thresholds,
- "high"
+ "high",
+ confidence_levels
);
// check results are identical
--
GitLab