diff --git a/mlair/helpers/statistics.py b/mlair/helpers/statistics.py index 87f780f9a6133edfcb2f9c71c2956b92f332e915..8b510f704396b35cf022089e7d94e037f4c62a2b 100644 --- a/mlair/helpers/statistics.py +++ b/mlair/helpers/statistics.py @@ -287,7 +287,7 @@ class SkillScores: combination_strings = [f"{first}-{second}" for (first, second) in combinations] return combinations, combination_strings - def skill_scores(self) -> pd.DataFrame: + def skill_scores(self) -> [pd.DataFrame, pd.DataFrame]: """ Calculate skill scores for all combinations of model names. @@ -296,6 +296,7 @@ class SkillScores: ahead_names = list(self.external_data[self.ahead_dim].data) combinations, combination_strings = self.get_model_name_combinations() skill_score = pd.DataFrame(index=combination_strings) + count = pd.DataFrame(index=combination_strings) for iahead in ahead_names: data = self.external_data.sel({self.ahead_dim: iahead}) skill_score[iahead] = [self.general_skill_score(data, @@ -303,7 +304,12 @@ class SkillScores: reference_name=second, observation_name=self.observation_name) for (first, second) in combinations] - return skill_score + count[iahead] = [self.get_count(data, + forecast_name=first, + reference_name=second, + observation_name=self.observation_name) + for (first, second) in combinations] + return skill_score, count def climatological_skill_scores(self, internal_data: Data, forecast_name: str) -> xr.DataArray: """ @@ -376,6 +382,21 @@ class SkillScores: skill_score = 1 - mse(observation, forecast) / mse(observation, reference) return skill_score.values + def get_count(self, data: Data, forecast_name: str, reference_name: str, + observation_name: str = None) -> np.ndarray: + r""" + Calculate general skill score based on mean squared error. + + :param data: internal data containing data for observation, forecast and reference + :param observation_name: name of observation + :param forecast_name: name of forecast + :param reference_name: name of reference + + :return: skill score of forecast + """ + data = data.dropna("index") + return data.count("index").max().values + def skill_score_pre_calculations(self, data: Data, observation_name: str, forecast_name: str) -> Tuple[np.ndarray, np.ndarray, np.ndarray, diff --git a/mlair/plotting/postprocessing_plotting.py b/mlair/plotting/postprocessing_plotting.py index 7ba84656cc00a89a7ec88efcc30bf665e0849e2f..5e8b121ae26ed0d4a967b9c440071d61c1ad703b 100644 --- a/mlair/plotting/postprocessing_plotting.py +++ b/mlair/plotting/postprocessing_plotting.py @@ -491,12 +491,12 @@ class PlotCompetitiveSkillScore(AbstractPlotClass): """ - def __init__(self, data: pd.DataFrame, plot_folder=".", model_setup="NN"): + def __init__(self, data: Dict[str, pd.DataFrame], plot_folder=".", model_setup="NN"): """Initialise.""" super().__init__(plot_folder, f"skill_score_competitive_{model_setup}") self._model_setup = model_setup self._labels = None - self._data = self._prepare_data(data) + self._data = self._prepare_data(helpers.remove_items(data, "total")) default_plot_name = self.plot_name # draw full detail plot self.plot_name = default_plot_name + "_full_detail" diff --git a/mlair/run_modules/post_processing.py b/mlair/run_modules/post_processing.py index 62cc85edefbaede065d6a571531a13da2c04aac3..3ce53c3c90921aa6bed742c5702104a5349812b8 100644 --- a/mlair/run_modules/post_processing.py +++ b/mlair/run_modules/post_processing.py @@ -117,7 +117,7 @@ class PostProcessing(RunEnvironment): # skill scores and error metrics with TimeTracking(name="calculate skill scores"): - skill_score_competitive, skill_score_climatological, errors = self.calculate_error_metrics() + skill_score_competitive, _, skill_score_climatological, errors = self.calculate_error_metrics() self.skill_scores = (skill_score_competitive, skill_score_climatological) self.report_error_metrics(errors) self.report_error_metrics({self.forecast_indicator: skill_score_climatological}) @@ -746,7 +746,7 @@ class PostProcessing(RunEnvironment): except (TypeError, AttributeError): return forecast if competitor is None else competitor - def calculate_error_metrics(self) -> Tuple[Dict, Dict, Dict]: + def calculate_error_metrics(self) -> Tuple[Dict, Dict, Dict, Dict]: """ Calculate error metrics and skill scores of NN forecast. @@ -759,6 +759,7 @@ class PostProcessing(RunEnvironment): path = self.data_store.get("forecast_path") all_stations = self.data_store.get("stations") skill_score_competitive = {} + skill_score_competitive_count = {} skill_score_climatological = {} errors = {} for station in all_stations: @@ -796,7 +797,7 @@ class PostProcessing(RunEnvironment): # skill score skill_score = statistics.SkillScores(combined, models=model_list, ahead_dim=self.ahead_dim) if external_data is not None: - skill_score_competitive[station] = skill_score.skill_scores() + skill_score_competitive[station], skill_score_competitive_count[station] = skill_score.skill_scores() internal_data = self._get_internal_data(station, path) if internal_data is not None: @@ -805,7 +806,21 @@ class PostProcessing(RunEnvironment): for model_type in errors.keys(): errors[model_type].update({"total": self.calculate_average_errors(errors[model_type])}) - return skill_score_competitive, skill_score_climatological, errors + skill_score_competitive.update({"total": self.calculate_average_skill_scores(skill_score_competitive, + skill_score_competitive_count)}) + return skill_score_competitive, skill_score_competitive_count, skill_score_climatological, errors + + @staticmethod + def calculate_average_skill_scores(scores, counts): + avg_skill_score = 0 + n_total = None + for vals in counts.values(): + n_total = vals if n_total is None else n_total + vals + for station, station_errors in scores.items(): + n = counts.get(station) + avg_skill_score = station_errors * n / n_total + avg_skill_score + return avg_skill_score + @staticmethod def calculate_average_errors(errors):