diff --git a/mlair/helpers/statistics.py b/mlair/helpers/statistics.py
index 5c424aba0d4656441a5f66745ad6d09874bb4e2b..57143c9aed0e730c81adbed33f7ba62fea39b298 100644
--- a/mlair/helpers/statistics.py
+++ b/mlair/helpers/statistics.py
@@ -569,9 +569,10 @@ def create_n_bootstrap_realizations(data: xr.DataArray, dim_name_time: str, dim_
"""
res_dims = [dim_name_boots]
dims = list(data.dims)
- coords = {dim_name_boots: range(n_boots), dim_name_model: data.coords[dim_name_model] }
+ other_dims = [v for v in dims if v in set(dims).difference([dim_name_time])]
+ coords = {dim_name_boots: range(n_boots), **{dim_name: data.coords[dim_name] for dim_name in other_dims}}
if len(dims) > 1:
- res_dims = res_dims + dims[1:]
+ res_dims = res_dims + other_dims
res = xr.DataArray(np.nan, dims=res_dims, coords=coords)
for boot in range(n_boots):
res[boot] = (calculate_average(
diff --git a/mlair/plotting/postprocessing_plotting.py b/mlair/plotting/postprocessing_plotting.py
index 0cdfada5f201fa0cd7b7eb69c6bc135d53138f86..bd2012c3fe9f53e7e07bfb4bfc4cde096c2dc891 100644
--- a/mlair/plotting/postprocessing_plotting.py
+++ b/mlair/plotting/postprocessing_plotting.py
@@ -6,7 +6,8 @@ import logging
import math
import os
import warnings
-from typing import Dict, List, Tuple
+from typing import Dict, List, Tuple, Union
+import itertools
import matplotlib
import matplotlib.pyplot as plt
@@ -538,9 +539,10 @@ class PlotCompetitiveSkillScore(AbstractPlotClass): # pragma: no cover
def _plot(self, single_model_comparison=False):
"""Plot skill scores of the comparisons."""
- size = max([len(np.unique(self._data.comparison)), 6])
- fig, ax = plt.subplots(figsize=(size, size * 0.8))
data = self._filter_comparisons(self._data) if single_model_comparison is True else self._data
+ max_label_size = len(max(np.unique(data.comparison).tolist(), key=len))
+ size = max([len(np.unique(data.comparison)), 6])
+ fig, ax = plt.subplots(figsize=(size, 5 * max(0.8, max_label_size/20)))
order = self._create_pseudo_order(data)
sns.boxplot(x="comparison", y="data", hue="ahead", data=data, whis=1.5, ax=ax, palette="Blues_d",
showmeans=True, meanprops={"markersize": 3, "markeredgecolor": "k"}, flierprops={"marker": "."},
@@ -554,8 +556,10 @@ class PlotCompetitiveSkillScore(AbstractPlotClass): # pragma: no cover
def _plot_vertical(self, single_model_comparison=False):
"""Plot skill scores of the comparisons, but vertically aligned."""
- fig, ax = plt.subplots()
data = self._filter_comparisons(self._data) if single_model_comparison is True else self._data
+ max_label_size = len(max(np.unique(data.comparison).tolist(), key=len))
+ size = max([len(np.unique(data.comparison)), 6])
+ fig, ax = plt.subplots(figsize=(5 * max(0.8, max_label_size/20), size))
order = self._create_pseudo_order(data)
sns.boxplot(y="comparison", x="data", hue="ahead", data=data, whis=1.5, ax=ax, palette="Blues_d",
showmeans=True, meanprops={"markersize": 3, "markeredgecolor": "k"}, flierprops={"marker": "."},
@@ -1090,34 +1094,36 @@ class PlotSampleUncertaintyFromBootstrap(AbstractPlotClass): # pragma: no cover
def __init__(self, data: xr.DataArray, plot_folder: str = ".", model_type_dim: str = "type",
error_measure: str = "mse", error_unit: str = None, dim_name_boots: str = 'boots',
- block_length: str = None, model_name: str = "NN", model_indicator: str = "nn"):
+ block_length: str = None, model_name: str = "NN", model_indicator: str = "nn",
+ ahead_dim: str = "ahead", sampling: Union[str, Tuple[str]] = ""):
super().__init__(plot_folder, "sample_uncertainty_from_bootstrap")
- default_name = self.plot_name
+ self.default_plot_name = self.plot_name
self.model_type_dim = model_type_dim
+ self.ahead_dim = ahead_dim
self.error_measure = error_measure
self.dim_name_boots = dim_name_boots
self.error_unit = error_unit
self.block_length = block_length
self.model_name = model_name
+ self.sampling = {"daily": "d", "hourly": "H"}.get(sampling[1] if isinstance(sampling, tuple) else sampling, "")
data = self.rename_model_indicator(data, model_name, model_indicator)
self.prepare_data(data)
- variants = [("_vertical", "v", False), ("_vertical_u_test", "v", True),
- ("_horizontal", "h", False), ("_horizontal_u_test", "h", True)]
+ # create all combinations to plot (h/v, utest/notest, single/multi)
+ variants = list(itertools.product(*[["v", "h"], [True, False], ["single", "multi"]]))
- for name_tag, orientation, utest in variants:
- self.plot_name = default_name + name_tag
- self._plot(orientation=orientation, apply_u_test=utest)
+ # plot raw metric (mse)
+ for orientation, utest, agg_type in variants:
+ self._plot(orientation=orientation, apply_u_test=utest, agg_type=agg_type)
+ # plot root of metric (rmse)
self._apply_root()
- variants = [(tag+"_sqrt", ori, ut) for tag, ori, ut in variants]
-
- for name_tag, orientation, utest in variants:
- self.plot_name = default_name + name_tag
- self._plot(orientation=orientation, apply_u_test=utest)
+ for orientation, utest, agg_type in variants:
+ self._plot(orientation=orientation, apply_u_test=utest, agg_type=agg_type, tag="_sqrt")
self._data_table = None
self._n_boots = None
+ self._factor = None
@property
def get_asteriks_from_mann_whitney_u_result(self):
@@ -1131,20 +1137,29 @@ class PlotSampleUncertaintyFromBootstrap(AbstractPlotClass): # pragma: no cover
return data
def prepare_data(self, data: xr.DataArray):
- data_table = data.to_pandas()
- self._data_table = data_table[data_table.mean().sort_values().index]
- self._n_boots = self._data_table.shape[0]
+ data_table = data.to_dataframe(self.model_type_dim).unstack()
+ factor = len(data.coords[self.ahead_dim]) if self.ahead_dim in data.dims else 1
+ self._data_table = data_table[data_table.mean().sort_values().index].droplevel(0, axis=1)
+ self._n_boots = int(self._data_table.shape[0] / factor)
+ self._factor = factor
def _apply_root(self):
self._data_table = np.sqrt(self._data_table)
self.error_measure = f"root {self.error_measure}"
self.error_unit = self.error_unit.replace("$^2$", "")
- def _plot(self, orientation: str = "v", apply_u_test: bool = False):
+ def _plot(self, orientation: str = "v", apply_u_test: bool = False, agg_type="single", tag=""):
+ self.plot_name = self.default_plot_name + {"v": "_vertical", "h": "_horizontal"}[orientation] + \
+ {True: "_u_test", False: ""}[apply_u_test] + "_" + agg_type + tag
+ if apply_u_test is True and agg_type == "multi":
+ return # not implemented
data_table = self._data_table
+ if self.ahead_dim not in data_table.index.names and agg_type == "multi":
+ return # nothing to do
n_boots = self._n_boots
size = len(np.unique(data_table.columns))
- asteriks = self.get_asteriks_from_mann_whitney_u_result
+ asteriks = self.get_asteriks_from_mann_whitney_u_result if apply_u_test is True else None
+ color_palette = sns.color_palette("Blues_d", self._factor).as_hex()
if orientation == "v":
figsize, width = (size, 5), 0.4
elif orientation == "h":
@@ -1152,11 +1167,24 @@ class PlotSampleUncertaintyFromBootstrap(AbstractPlotClass): # pragma: no cover
else:
raise ValueError(f"orientation must be `v' or `h' but is: {orientation}")
fig, ax = plt.subplots(figsize=figsize)
- sns.boxplot(data=data_table, ax=ax, whis=1.5, color="white",
- showmeans=True, meanprops={"markersize": 6, "markeredgecolor": "k"},
- flierprops={"marker": "o", "markerfacecolor": "black", "markeredgecolor": "none", "markersize": 3},
- boxprops={'facecolor': 'none', 'edgecolor': 'k'},
- width=width, orient=orientation)
+ if agg_type == "single":
+ if self.ahead_dim in data_table.index.names:
+ data_table = data_table.groupby(level=0).mean()
+ sns.boxplot(data=data_table, ax=ax, whis=1.5, color="white",
+ showmeans=True, meanprops={"markersize": 6, "markeredgecolor": "k"},
+ flierprops={"marker": "o", "markerfacecolor": "black", "markeredgecolor": "none", "markersize": 3},
+ boxprops={'facecolor': 'none', 'edgecolor': 'k'}, width=width, orient=orientation)
+ else:
+ xy = {"x": self.model_type_dim, "y": 0} if orientation == "v" else {"x": 0, "y": self.model_type_dim}
+ sns.boxplot(data=data_table.stack(self.model_type_dim).reset_index(), ax=ax, whis=1.5, palette=color_palette,
+ showmeans=True, meanprops={"markersize": 6, "markeredgecolor": "k", "markerfacecolor": "white"},
+ flierprops={"marker": "o", "markerfacecolor": "black", "markeredgecolor": "none", "markersize": 3},
+ boxprops={'edgecolor': 'k'}, width=.9, orient=orientation, **xy, hue=self.ahead_dim)
+
+ _labels = [str(i) + self.sampling for i in data_table.index.levels[1].values]
+ handles, _ = ax.get_legend_handles_labels()
+ ax.legend(handles, _labels)
+
if orientation == "v":
if apply_u_test:
ax = self.set_significance_bars(asteriks, ax, data_table, orientation)
diff --git a/mlair/run_modules/post_processing.py b/mlair/run_modules/post_processing.py
index 8b69034fa5d3b69b39ba74aa8a6bbab873c06cc8..8f6bf05d29b2534e4918d72aa59f91aace0ec982 100644
--- a/mlair/run_modules/post_processing.py
+++ b/mlair/run_modules/post_processing.py
@@ -111,7 +111,7 @@ class PostProcessing(RunEnvironment):
# sample uncertainty
if self.data_store.get("do_uncertainty_estimate", "postprocessing"):
- self.estimate_sample_uncertainty()
+ self.estimate_sample_uncertainty(separate_ahead=True)
# feature importance bootstraps
if self.data_store.get("evaluate_feature_importance", "postprocessing"):
@@ -170,12 +170,27 @@ class PostProcessing(RunEnvironment):
# store statistics
if percentiles is None:
percentiles = [.05, .1, .25, .5, .75, .9, .95]
- df_descr = self.uncertainty_estimate.to_pandas().describe(percentiles=percentiles).astype("float32")
- column_format = tables.create_column_format_for_tex(df_descr)
- file_name = os.path.join(report_path, "uncertainty_estimate_statistics.%s")
- tables.save_to_tex(report_path, file_name % "tex", column_format=column_format, df=df_descr)
- tables.save_to_md(report_path, file_name % "md", df=df_descr)
- df_descr.to_csv(file_name % "csv", sep=";")
+
+ for ahead_steps in ["single", "multi"]:
+ if ahead_steps == "single":
+ try:
+ df_descr = self.uncertainty_estimate.to_pandas().describe(percentiles=percentiles).astype("float32")
+ except ValueError:
+ df_descr = self.uncertainty_estimate.mean(self.ahead_dim).to_pandas().describe(percentiles=percentiles).astype("float32")
+ else:
+ if self.ahead_dim not in self.uncertainty_estimate.dims:
+ continue
+ df_descr = self.uncertainty_estimate.to_dataframe(self.model_type_dim).unstack().groupby(level=self.ahead_dim).describe(
+ percentiles=percentiles).astype("float32")
+ df_descr = df_descr.stack(-1)
+ df_descr = df_descr.reorder_levels(df_descr.index.names[::-1])
+ df_sorter = ["count", "mean", "std", "min", *[f"{round(p * 100)}%" for p in percentiles], "max"]
+ df_descr = df_descr.loc[df_sorter]
+ column_format = tables.create_column_format_for_tex(df_descr)
+ file_name = os.path.join(report_path, f"uncertainty_estimate_statistics_{ahead_steps}.%s")
+ tables.save_to_tex(report_path, file_name % "tex", column_format=column_format, df=df_descr)
+ tables.save_to_md(report_path, file_name % "md", df=df_descr)
+ df_descr.to_csv(file_name % "csv", sep=";")
def calculate_block_mse(self, evaluate_competitors=True, separate_ahead=False, block_length="1m"):
"""
@@ -566,6 +581,18 @@ class PostProcessing(RunEnvironment):
logging.error(f"Could not create plot PlotTimeSeries due to the following error:\n{sys.exc_info()[0]}\n"
f"{sys.exc_info()[1]}\n{sys.exc_info()[2]}\n{traceback.format_exc()}")
+ try:
+ if "PlotSampleUncertaintyFromBootstrap" in plot_list and self.uncertainty_estimate is not None:
+ block_length = self.data_store.get_default("block_length", default="1m", scope="uncertainty_estimate")
+ PlotSampleUncertaintyFromBootstrap(
+ data=self.uncertainty_estimate, plot_folder=self.plot_path, model_type_dim=self.model_type_dim,
+ dim_name_boots=self.uncertainty_estimate_boot_dim, error_measure="mean squared error",
+ error_unit=r"ppb$^2$", block_length=block_length, model_name=self.model_display_name,
+ model_indicator=self.forecast_indicator, sampling=self._sampling)
+ except Exception as e:
+ logging.error(f"Could not create plot PlotSampleUncertaintyFromBootstrap due to the following error: {e}"
+ f"\n{sys.exc_info()[0]}\n{sys.exc_info()[1]}\n{sys.exc_info()[2]}")
+
try:
if "PlotStationMap" in plot_list:
if self.data_store.get("hostname")[:2] in self.data_store.get("hpc_hosts") or self.data_store.get(
@@ -602,15 +629,6 @@ class PostProcessing(RunEnvironment):
logging.error(f"Could not create plot PlotAvailabilityHistogram due to the following error: {e}"
f"\n{sys.exc_info()[0]}\n{sys.exc_info()[1]}\n{sys.exc_info()[2]}")
- try:
- if "PlotPeriodogram" in plot_list:
- PlotPeriodogram(self.train_data, plot_folder=self.plot_path, time_dim=time_dim,
- variables_dim=target_dim, sampling=self._sampling,
- use_multiprocessing=use_multiprocessing)
- except Exception as e:
- logging.error(f"Could not create plot PlotPeriodogram due to the following error: {e}"
- f"\n{sys.exc_info()[0]}\n{sys.exc_info()[1]}\n{sys.exc_info()[2]}")
-
try:
if "PlotDataHistogram" in plot_list:
upsampling = self.data_store.get_default("upsampling", scope="train", default=False)
@@ -622,15 +640,12 @@ class PostProcessing(RunEnvironment):
f"\n{sys.exc_info()[0]}\n{sys.exc_info()[1]}\n{sys.exc_info()[2]}")
try:
- if "PlotSampleUncertaintyFromBootstrap" in plot_list and self.uncertainty_estimate is not None:
- block_length = self.data_store.get_default("block_length", default="1m", scope="uncertainty_estimate")
- PlotSampleUncertaintyFromBootstrap(
- data=self.uncertainty_estimate, plot_folder=self.plot_path, model_type_dim=self.model_type_dim,
- dim_name_boots=self.uncertainty_estimate_boot_dim, error_measure="mean squared error",
- error_unit=r"ppb$^2$", block_length=block_length, model_name=self.model_display_name,
- model_indicator=self.forecast_indicator)
+ if "PlotPeriodogram" in plot_list:
+ PlotPeriodogram(self.train_data, plot_folder=self.plot_path, time_dim=time_dim,
+ variables_dim=target_dim, sampling=self._sampling,
+ use_multiprocessing=use_multiprocessing)
except Exception as e:
- logging.error(f"Could not create plot PlotSampleUncertaintyFromBootstrap due to the following error: {e}"
+ logging.error(f"Could not create plot PlotPeriodogram due to the following error: {e}"
f"\n{sys.exc_info()[0]}\n{sys.exc_info()[1]}\n{sys.exc_info()[2]}")
@TimeTrackingWrapper