Skip to content
GitLab
Commit 51ad8b59 authored by Thibault Hallouin's avatar Thibault Hallouin
Browse files

Merge branch 'main' into 1-add-probabilistic-scores

parents 0690cded 8bde8cd3
main 11-check-the-computation-of-the-ranks-for-ties add-determinist-to-evalp dev v0.1.2 v0.1.1 v0.1.0 v0.0.1
1 merge request!1implement quantile scores and CRPS
Hide whitespace changes
Inline Side-by-side
Showing
with 76 additions and 36 deletions
+76 -36
include/evalhyd/determinist.hpp
+ 5
− 4
View file @ 51ad8b59
......@@ -2,7 +2,7 @@
#define EVALHYD_DETERMINIST_HPP
#include <unordered_map>
#include <unordered_set>
#include <vector>
#include <xtensor/xexpression.hpp>
#include "utils.hpp"
......@@ -49,8 +49,8 @@ namespace evalhyd
// TODO
// determine required elt/dep to be pre-computed
std::unordered_set<std::string> req_elt = {};
std::unordered_set<std::string> req_dep = {};
std::vector<std::string> req_elt;
std::vector<std::string> req_dep;
eh::utils::find_requirements(metrics, elt, dep, req_elt, req_dep);
......@@ -73,7 +73,8 @@ namespace evalhyd
{
if ( metric == "nse" )
{
if ( req_dep.find(metric) == req_dep.end() )
if (std::find(req_dep.begin(), req_dep.end(), metric)
== req_dep.end())
evaluator.calc_nse();
r.emplace_back(evaluator.nse);
}
......
include/evalhyd/probabilist.hpp
+ 15
− 9
View file @ 51ad8b59
......@@ -3,7 +3,7 @@
#include <utility>
#include <unordered_map>
#include <unordered_set>
#include <vector>
#include <xtensor/xtensor.hpp>
#include <xtensor/xview.hpp>
......@@ -47,16 +47,16 @@ namespace evalhyd
// register potentially recurring computation elements across metrics
elt["bs"] = {"o_k", "y_k"};
elt["bss"] = {"o_k"};
elt["bs_crd"] = {"o_k", "y_k"};
elt["bss"] = {"o_k", "bar_o"};
elt["bs_crd"] = {"o_k", "bar_o", "y_k"};
elt["bs_lbd"] = {"o_k", "y_k"};
// register nested metrics (i.e. metric dependent on another metric)
dep["bss"] = {"bs"};
// determine required elt/dep to be pre-computed
std::unordered_set<std::string> req_elt = {};
std::unordered_set<std::string> req_dep = {};
std::vector<std::string> req_elt;
std::vector<std::string> req_dep;
eh::utils::find_requirements(metrics, elt, dep, req_elt, req_dep);
......@@ -65,6 +65,8 @@ namespace evalhyd
{
if ( element == "o_k" )
evaluator.calc_o_k();
else if ( element == "bar_o" )
evaluator.calc_bar_o();
else if ( element == "y_k" )
evaluator.calc_y_k();
}
......@@ -83,25 +85,29 @@ namespace evalhyd
{
if ( metric == "bs" )
{
if (req_dep.find(metric) == req_dep.end())
if (std::find(req_dep.begin(), req_dep.end(), metric)
== req_dep.end())
evaluator.calc_bs();
r.emplace_back(evaluator.bs);
}
else if ( metric == "bss" )
{
if (req_dep.find(metric) == req_dep.end())
if (std::find(req_dep.begin(), req_dep.end(), metric)
== req_dep.end())
evaluator.calc_bss();
r.emplace_back(evaluator.bss);
}
else if ( metric == "bs_crd" )
{
if (req_dep.find(metric) == req_dep.end())
if (std::find(req_dep.begin(), req_dep.end(), metric)
== req_dep.end())
evaluator.calc_bs_crd();
r.emplace_back(evaluator.bs_crd);
}
else if ( metric == "bs_lbd" )
{
if (req_dep.find(metric) == req_dep.end())
if (std::find(req_dep.begin(), req_dep.end(), metric)
== req_dep.end())
evaluator.calc_bs_lbd();
r.emplace_back(evaluator.bs_lbd);
}
......
include/evalhyd/probabilist/evaluator.h
+ 2
− 0
View file @ 51ad8b59
......@@ -17,6 +17,7 @@ namespace evalhyd
// members for computational elements
xt::xtensor<double, 2> o_k;
xt::xtensor<double, 1> bar_o;
xt::xtensor<double, 2> y_k;
public:
......@@ -37,6 +38,7 @@ namespace evalhyd
// methods to compute elements
void calc_o_k();
void calc_bar_o();
void calc_y_k();
// methods to compute metrics
......
include/evalhyd/probabilist/evaluator_brier.cpp
+ 12
− 15
View file @ 51ad8b59
#include <unordered_map>
#include <xtensor/xmath.hpp>
#include <xtensor/xview.hpp>
#include <xtensor/xoperation.hpp>
......@@ -59,8 +58,6 @@ namespace evalhyd
// declare internal variables
// shape: (bins, thresholds, time)
xt::xtensor<double, 3> msk_bins;
// shape: (thresholds,)
xt::xtensor<double, 1> bar_o;
// shape: (bins, thresholds)
xt::xtensor<double, 2> N_i, bar_o_i;
// shape: (bins,)
......@@ -104,10 +101,6 @@ namespace evalhyd
0.
);
// compute mean "climatology" relative frequency of the event
// $\bar{o} = \frac{1}{n} \sum_{k=1}^{n} o_k$
bar_o = xt::mean(o_k, -1);
// calculate reliability =
// $\frac{1}{n} \sum_{i=1}^{I} N_i (y_i - \bar{o_i})^2$
xt::col(bs_crd, 0) =
......@@ -239,6 +232,9 @@ namespace evalhyd
// \require o_k:
// Observed event outcome.
// shape: (thresholds, time)
// \require bar_o:
// Mean event observed outcome.
// shape: (thresholds,)
// \require bs:
// Brier score(s).
// shape: (thresholds,)
......@@ -247,15 +243,16 @@ namespace evalhyd
// shape: (thresholds, 1)
void Evaluator::calc_bss()
{
// calculate mean observed event outcome
// $\bar{o_k} = \frac{1}{n} \sum_{k=1}^{n} o_k$
xt::xtensor<double, 2> bar_o_k = xt::mean(o_k, -1, xt::keep_dims);
// calculate reference Brier score(s)
// $BS_{ref} = \frac{1}{n} \sum_{k=1}^{n} (o_k - \bar{o_k})^2$
xt::xtensor<double, 2> bs_ref = xt::mean(
xt::square(o_k - bar_o_k), -1, xt::keep_dims
);
// $BS_{ref} = \frac{1}{n} \sum_{k=1}^{n} (o_k - \bar{o})^2$
xt::xtensor<double, 2> bs_ref =
xt::mean(
xt::square(
o_k -
xt::view(bar_o, xt::all(), xt::newaxis())
),
-1, xt::keep_dims
);
// return computed Brier skill score(s)
// $BSS = 1 - \frac{BS}{BS_{ref}}
......
include/evalhyd/probabilist/evaluator_elements.cpp
+ 16
− 1
View file @ 51ad8b59
......@@ -23,6 +23,21 @@ namespace evalhyd
o_k = xt::squeeze(is_above_threshold(q_obs, q_thr));
}
// Determine mean observed realisation of threshold(s) exceedance.
//
// \require o_k:
// 2D array of event observed outcome.
// shape: (thresholds, time)
// \assign bar_o:
// 1D array of mean event observed outcome.
// shape: (thresholds,)
void Evaluator::calc_bar_o()
{
// compute mean "climatology" relative frequency of the event
// $\bar{o} = \frac{1}{n} \sum_{k=1}^{n} o_k$
bar_o = xt::mean(o_k, -1);
}
// Determine forecast probability of threshold(s) exceedance to occur.
//
// \require q_frc:
......@@ -31,7 +46,7 @@ namespace evalhyd
// \require q_thr:
// 1D array of streamflow exceedance threshold(s).
// shape: (thresholds,)
// \assign y_k
// \assign y_k:
// 2D array of event probability forecast.
// shape: (thresholds, time)
void Evaluator::calc_y_k()
......
include/evalhyd/probabilist/evaluator_utils.cpp
+ 4
− 2
View file @ 51ad8b59
......@@ -19,7 +19,8 @@ namespace evalhyd
xt::xtensor<double, 3> is_above_threshold(
const xt::xtensor<double, 2> &q,
const xt::xtensor<double, 1> &thr
) {
)
{
return q >= xt::view(thr, xt::all(), xt::newaxis(), xt::newaxis());
}
......@@ -37,7 +38,8 @@ namespace evalhyd
xt::xtensor<double, 3> is_below_threshold(
const xt::xtensor<double, 2> &q,
const xt::xtensor<double, 1> &thr
) {
)
{
return q < xt::view(thr, xt::all(), xt::newaxis(), xt::newaxis());
}
}
......
include/evalhyd/utils.hpp
+ 22
− 5
View file @ 51ad8b59
......@@ -3,6 +3,7 @@
#include <unordered_map>
#include <unordered_set>
#include <vector>
#include <xtensor/xtensor.hpp>
namespace evalhyd
......@@ -27,24 +28,40 @@ namespace evalhyd
const std::vector<std::string>& metrics,
std::unordered_map<std::string, std::vector<std::string>>& elements,
std::unordered_map<std::string, std::vector<std::string>>& dependencies,
std::unordered_set<std::string>& required_elements,
std::unordered_set<std::string>& required_dependencies
std::vector<std::string>& required_elements,
std::vector<std::string>& required_dependencies
)
{
std::unordered_set<std::string> found_elements;
std::unordered_set<std::string> found_dependencies;
for (const auto& metric : metrics)
{
// add elements to pre-computation set
for (const auto& element : elements[metric])
required_elements.insert(element);
if (found_elements.find(element) == found_elements.end())
{
found_elements.insert(element);
required_elements.push_back(element);
}
// add metric dependencies to pre-computation set
if (dependencies.find(metric) != dependencies.end())
{
for (const auto& dependency : dependencies[metric])
{
required_dependencies.insert(dependency);
if (found_dependencies.find(dependency) == found_dependencies.end())
{
found_dependencies.insert(dependency);
required_dependencies.push_back(dependency);
}
// add dependency elements to pre-computation set
for (const auto& element : elements[dependency])
required_elements.insert(element);
if (found_elements.find(element) == found_elements.end())
{
found_elements.insert(element);
required_elements.push_back(element);
}
}
}
}
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment