From e0e5a0feff35b990aac5c961223259a368359e1b Mon Sep 17 00:00:00 2001 From: Le Roux Erwan <erwan.le-roux@irstea.fr> Date: Thu, 11 Jun 2020 17:52:01 +0200 Subject: [PATCH] [contrasting] add display of name for NON_STATIONARY_TREND_TEST_PAPER_2 --- .../visualization/plot_utils.py | 4 ++- .../visualization/study_visualizer.py | 10 ++++--- extreme_trend/abstract_gev_trend_test.py | 28 +++++++++++++++++-- .../gev_trend_test_two_parameters.py | 12 ++++---- .../main_snowfall_article.py | 10 ++++--- .../study_visualizer_for_mean_values.py | 23 +++++++++++++-- .../validation_plot.py | 18 ++++++------ projects/exceeding_snow_loads/utils.py | 11 +++++++- 8 files changed, 88 insertions(+), 28 deletions(-) diff --git a/extreme_data/meteo_france_data/scm_models_data/visualization/plot_utils.py b/extreme_data/meteo_france_data/scm_models_data/visualization/plot_utils.py index 1779e184..525f0cfc 100644 --- a/extreme_data/meteo_france_data/scm_models_data/visualization/plot_utils.py +++ b/extreme_data/meteo_france_data/scm_models_data/visualization/plot_utils.py @@ -23,7 +23,7 @@ def plot_against_altitude(altitudes, ax, massif_id, massif_name, values): def load_plot(cmap, graduation, label, massif_name_to_value, altitude, fit_method, add_x_label=True, - negative_and_positive_values=True): + negative_and_positive_values=True, massif_name_to_text=None): max_abs_change = max([abs(e) for e in massif_name_to_value.values()]) if negative_and_positive_values: ticks, labels = ticks_values_and_labels_for_percentages(graduation=graduation, max_abs_change=max_abs_change) @@ -54,6 +54,8 @@ def load_plot(cmap, graduation, label, massif_name_to_value, altitude, fit_metho ticks_values_and_labels=ticks_values_and_labels, label=label, fontsize_label=10, + massif_name_to_text=massif_name_to_text, + add_text=massif_name_to_text is not None, ) ax.get_xaxis().set_visible(True) ax.set_xticks([]) diff --git a/extreme_data/meteo_france_data/scm_models_data/visualization/study_visualizer.py b/extreme_data/meteo_france_data/scm_models_data/visualization/study_visualizer.py index 3288a89d..cbccba7b 100644 --- a/extreme_data/meteo_france_data/scm_models_data/visualization/study_visualizer.py +++ b/extreme_data/meteo_france_data/scm_models_data/visualization/study_visualizer.py @@ -709,20 +709,22 @@ class StudyVisualizer(VisualizationParameters): # PLot functions that should be common def plot_map(self, cmap, fit_method, graduation, label, massif_name_to_value, plot_name, add_x_label=True, - negative_and_positive_values=True): + negative_and_positive_values=True, massif_name_to_text=None): load_plot(cmap, graduation, label, massif_name_to_value, self.study.altitude, fitmethod_to_str(fit_method), - add_x_label=add_x_label, negative_and_positive_values=negative_and_positive_values) + add_x_label=add_x_label, negative_and_positive_values=negative_and_positive_values, + massif_name_to_text=massif_name_to_text) self.plot_name = plot_name # self.show_or_save_to_file(add_classic_title=False, tight_layout=True, no_title=True, dpi=500) self.show_or_save_to_file(add_classic_title=False, no_title=True, dpi=500) plt.close() def plot_abstract(self, massif_name_to_value, label, plot_name, fit_method='', graduation=10.0, cmap=plt.cm.bwr, - add_x_label=True, negative_and_positive_values=True): + add_x_label=True, negative_and_positive_values=True, massif_name_to_text=None): # Regroup the plot by altitudes plot_name1 = '{}/{}'.format(self.study.altitude, plot_name) # Regroup the plot by type of plot also plot_name2 = '{}/{}'.format(plot_name.split()[0], plot_name) for plot_name in [plot_name1, plot_name2]: - self.plot_map(cmap, fit_method, graduation, label, massif_name_to_value, plot_name, add_x_label, negative_and_positive_values) + self.plot_map(cmap, fit_method, graduation, label, massif_name_to_value, plot_name, add_x_label, negative_and_positive_values, + massif_name_to_text) diff --git a/extreme_trend/abstract_gev_trend_test.py b/extreme_trend/abstract_gev_trend_test.py index feea4fd4..46224e06 100644 --- a/extreme_trend/abstract_gev_trend_test.py +++ b/extreme_trend/abstract_gev_trend_test.py @@ -79,6 +79,31 @@ class AbstractGevTrendTest(object): test = cramer_von_mises_and_anderson_darling_tests_pvalues_for_gumbel_distribution(quantiles) return test[1] > self.SIGNIFICANCE_LEVEL + @property + def name(self): + if self.constrained_model_class is GumbelTemporalModel: + family = 'Gum' + else: + family = 'Gev' + + if self.unconstrained_estimator.margin_model.mul == 1: + family += '1' + else: + family += '0' + + if self.unconstrained_estimator.margin_model.sigl == 1: + family += '1' + else: + family += '0' + + if family.startswith('Gev'): + if self.unconstrained_estimator.margin_model.shl == 1: + family += '1' + else: + family += '0' + + return family + @property def degree_freedom_chi2(self) -> int: raise NotImplementedError @@ -225,7 +250,6 @@ class AbstractGevTrendTest(object): ax.set_ylabel("Non-zero annual maxima of GSL ({})".format(AbstractSnowLoadVariable.UNIT), fontsize=size) ax.legend(prop={'size': 17}) - def intensity_plot_wrt_standard_gumbel(self, massif_name, altitude, psnow): ax = plt.gca() sorted_maxima = sorted(self.maxima) @@ -446,7 +470,7 @@ class AbstractGevTrendTest(object): def unconstrained_average_mean_value(self, year_min, year_max) -> float: mean_values = [] - for year in range(year_min, year_max+1): + for year in range(year_min, year_max + 1): mean = self.get_unconstrained_gev_params(year).mean mean_values.append(mean) average_mean_value = np.mean(mean_values) diff --git a/extreme_trend/trend_test_two_parameters/gev_trend_test_two_parameters.py b/extreme_trend/trend_test_two_parameters/gev_trend_test_two_parameters.py index 20b6edef..7830fb4c 100644 --- a/extreme_trend/trend_test_two_parameters/gev_trend_test_two_parameters.py +++ b/extreme_trend/trend_test_two_parameters/gev_trend_test_two_parameters.py @@ -17,7 +17,8 @@ class GevTrendTestTwoParameters(AbstractGevTrendTest): def degree_freedom_chi2(self) -> int: return 2 -class GevTrendTestTwoParametersAgainstGev(AbstractGevTrendTest): + +class GevTrendTestTwoParametersAgainstGev(GevTrendTestTwoParameters): @classproperty def total_number_of_parameters_for_unconstrained_model(cls) -> int: @@ -25,7 +26,7 @@ class GevTrendTestTwoParametersAgainstGev(AbstractGevTrendTest): class GevLocationAndShapeTrendTest(GevTrendTestTwoParametersAgainstGev): - + def __init__(self, years, maxima, starting_year, constrained_model_class=StationaryTemporalModel, quantile_level=EUROCODE_QUANTILE, fit_method=MarginFitMethod.extremes_fevd_mle): super().__init__(years, maxima, starting_year, @@ -37,7 +38,6 @@ class GevLocationAndShapeTrendTest(GevTrendTestTwoParametersAgainstGev): class GevScaleAndShapeTrendTest(GevTrendTestTwoParametersAgainstGev): - def __init__(self, years, maxima, starting_year, constrained_model_class=StationaryTemporalModel, quantile_level=EUROCODE_QUANTILE, fit_method=MarginFitMethod.extremes_fevd_mle): super().__init__(years, maxima, starting_year, @@ -83,7 +83,8 @@ class GevLocationAndScaleTrendTest(GevTrendTestTwoParametersAgainstGev): class GevLocationAgainstGumbel(GevTrendTestTwoParameters, GevLocationTrendTest): - def __init__(self, years, maxima, starting_year, quantile_level=EUROCODE_QUANTILE, fit_method=MarginFitMethod.extremes_fevd_mle): + def __init__(self, years, maxima, starting_year, quantile_level=EUROCODE_QUANTILE, + fit_method=MarginFitMethod.extremes_fevd_mle): super().__init__(years, maxima, starting_year, quantile_level, GumbelTemporalModel, fit_method=fit_method) @classproperty @@ -101,7 +102,8 @@ class GevLocationAgainstGumbel(GevTrendTestTwoParameters, GevLocationTrendTest): class GevScaleAgainstGumbel(GevTrendTestTwoParameters, GevScaleTrendTest): - def __init__(self, years, maxima, starting_year, quantile_level=EUROCODE_QUANTILE, fit_method=MarginFitMethod.extremes_fevd_mle): + def __init__(self, years, maxima, starting_year, quantile_level=EUROCODE_QUANTILE, + fit_method=MarginFitMethod.extremes_fevd_mle): super().__init__(years, maxima, starting_year, quantile_level, GumbelTemporalModel, fit_method=fit_method) @classproperty diff --git a/projects/contrasting_trends_in_snow_loads/article2_snowfall_versus_time_and_altitude/main_snowfall_article.py b/projects/contrasting_trends_in_snow_loads/article2_snowfall_versus_time_and_altitude/main_snowfall_article.py index db109f24..4438841c 100644 --- a/projects/contrasting_trends_in_snow_loads/article2_snowfall_versus_time_and_altitude/main_snowfall_article.py +++ b/projects/contrasting_trends_in_snow_loads/article2_snowfall_versus_time_and_altitude/main_snowfall_article.py @@ -9,6 +9,7 @@ from projects.contrasting_trends_in_snow_loads.article2_snowfall_versus_time_and from projects.contrasting_trends_in_snow_loads.article2_snowfall_versus_time_and_altitude.study_visualizer_for_mean_values import \ StudyVisualizerForMeanValues from projects.contrasting_trends_in_snow_loads.article2_snowfall_versus_time_and_altitude.validation_plot import validation_plot +from projects.exceeding_snow_loads.section_results.plot_selection_curves import plot_selection_curves from projects.exceeding_snow_loads.section_results.plot_trend_curves import plot_trend_map mpl.rcParams['text.usetex'] = True @@ -75,6 +76,7 @@ def intermediate_result(altitudes, massif_names=None, validation_plot(altitude_to_visualizer, order_derivative=0) # validation_plot(altitude_to_visualizer, order_derivative=1) plot_snowfall_mean(altitude_to_visualizer) + plot_selection_curves(altitude_to_visualizer, paper1=False) # plot_snowfall_time_derivative_mean(altitude_to_visualizer) @@ -85,10 +87,10 @@ def major_result(): study_classes = [SafranSnowfall1Day] model_subsets_for_uncertainty = None altitudes = paper_altitudes - altitudes = [300, 600, 900, 1200, 1500, 1800, 2100, 2400, 2700, 3000, 3300, 3600, 3900] - # altitudes = [900, 1200] - draft_altitudes = [900, 1200, 1500, 1800, 2100, 2400, 2700, 3000] - altitudes = draft_altitudes + # altitudes = [300, 600, 900, 1200, 1500, 1800, 2100, 2400, 2700, 3000, 3300, 3600, 3900] + altitudes = [900, 1200, 1500, 1800] + # draft_altitudes = [900, 1200, 1500, 1800, 2100, 2400, 2700, 3000] + # altitudes = draft_altitudes for study_class in study_classes: intermediate_result(altitudes, massif_names, model_subsets_for_uncertainty, diff --git a/projects/contrasting_trends_in_snow_loads/article2_snowfall_versus_time_and_altitude/study_visualizer_for_mean_values.py b/projects/contrasting_trends_in_snow_loads/article2_snowfall_versus_time_and_altitude/study_visualizer_for_mean_values.py index ab0586f6..351a50c2 100644 --- a/projects/contrasting_trends_in_snow_loads/article2_snowfall_versus_time_and_altitude/study_visualizer_for_mean_values.py +++ b/projects/contrasting_trends_in_snow_loads/article2_snowfall_versus_time_and_altitude/study_visualizer_for_mean_values.py @@ -1,3 +1,4 @@ +from collections import Counter from typing import Dict import matplotlib.pyplot as plt @@ -37,9 +38,23 @@ class StudyVisualizerForMeanValues(StudyVisualizerForNonStationaryTrends): assert len(self.non_stationary_trend_test) == len(NON_STATIONARY_TREND_TEST_PAPER_2) def plot_abstract_fast(self, massif_name_to_value, label, graduation=10.0, cmap=plt.cm.coolwarm, add_x_label=True, - negative_and_positive_values=True): + negative_and_positive_values=True, add_text=False): + massif_name_to_text = self.massif_name_to_text if add_text else None super().plot_abstract(massif_name_to_value, label, label, self.fit_method, graduation, cmap, add_x_label, - negative_and_positive_values) + negative_and_positive_values, massif_name_to_text) + + @property + def massif_name_to_text(self): + d = {} + for m, t in self.massif_name_to_trend_test_that_minimized_aic.items(): + if t.is_significant: + name = '$\\textbf{' + name += t.name + name += '}$' + else: + name = '${}$'.format(t.name) + d[m] = name + return d # Override the main dict massif_name_to_trend_test_that_minimized_aic @@ -48,6 +63,10 @@ class StudyVisualizerForMeanValues(StudyVisualizerForNonStationaryTrends): return {m: t for m, t in super().massif_name_to_trend_test_that_minimized_aic.items() if t.goodness_of_fit_anderson_test} + @property + def super_massif_name_to_trend_test_that_minimized_aic(self): + return super().massif_name_to_trend_test_that_minimized_aic + # Study of the mean @cached_property diff --git a/projects/contrasting_trends_in_snow_loads/article2_snowfall_versus_time_and_altitude/validation_plot.py b/projects/contrasting_trends_in_snow_loads/article2_snowfall_versus_time_and_altitude/validation_plot.py index 99479248..053570e2 100644 --- a/projects/contrasting_trends_in_snow_loads/article2_snowfall_versus_time_and_altitude/validation_plot.py +++ b/projects/contrasting_trends_in_snow_loads/article2_snowfall_versus_time_and_altitude/validation_plot.py @@ -23,9 +23,9 @@ def validation_plot(altitude_to_visualizer: Dict[int, StudyVisualizerForMeanValu for altitude, visualizer in altitude_to_visualizer.items(): altitude_to_relative_differences[altitude] = plot_function(visualizer) study_visualizer.show_or_save_to_file(add_classic_title=False, dpi=500) - # Shoe plot with respect to the altitude. - plot_shoe_relative_differences_distribution(altitude_to_relative_differences, altitudes, study_visualizer, - order_derivative) + # # Shoe plot with respect to the altitude. + # plot_shoe_relative_differences_distribution(altitude_to_relative_differences, altitudes, study_visualizer, + # order_derivative) study_visualizer.show_or_save_to_file(add_classic_title=False, dpi=500) plt.close() @@ -51,13 +51,13 @@ def plot_shoe_relative_differences_distribution(altitude_to_relative_differences def plot_relative_difference_map_order_zero(visualizer: StudyVisualizerForMeanValues): study = visualizer.study label = ' mean annual maxima of {} ({})'.format(SCM_STUDY_CLASS_TO_ABBREVIATION[type(study)], study.variable_unit) - visualizer.plot_abstract_fast(massif_name_to_value=visualizer.massif_name_to_empirical_mean, - label='Empirical' + label, negative_and_positive_values=False) + # visualizer.plot_abstract_fast(massif_name_to_value=visualizer.massif_name_to_empirical_mean, + # label='Empirical' + label, negative_and_positive_values=False) visualizer.plot_abstract_fast(massif_name_to_value=visualizer.massif_name_to_model_mean, - label='Model' + label, negative_and_positive_values=False) - visualizer.plot_abstract_fast(massif_name_to_value=visualizer.massif_name_to_relative_difference_for_mean, - label='Relative difference of the model mean w.r.t. the empirical mean \n' - 'for the ' + label, graduation=1) + label='Model' + label, negative_and_positive_values=False, add_text=True) + # visualizer.plot_abstract_fast(massif_name_to_value=visualizer.massif_name_to_relative_difference_for_mean, + # label='Relative difference of the model mean w.r.t. the empirical mean \n' + # 'for the ' + label, graduation=1) return list(visualizer.massif_name_to_relative_difference_for_mean.values()) diff --git a/projects/exceeding_snow_loads/utils.py b/projects/exceeding_snow_loads/utils.py index 0faed07a..2594b132 100644 --- a/projects/exceeding_snow_loads/utils.py +++ b/projects/exceeding_snow_loads/utils.py @@ -31,9 +31,18 @@ NON_STATIONARY_TREND_TEST_PAPER_2 = [ # Gumbel models GumbelVersusGumbel, GumbelLocationTrendTest, GumbelScaleTrendTest, GumbelLocationAndScaleTrendTest, # GEV models with constant shape - GevVersusGev, GevLocationTrendTest, GevScaleTrendTest, GevLocationAndScaleTrendTest, + GevVersusGev, GevLocationTrendTest, GevScaleTrendTest, GevLocationAndScaleTrendTest, # GEV models with linear shape GevShapeTrendTest, GevLocationAndShapeTrendTest, GevScaleAndShapeTrendTest, GevLocationAndScaleAndShapeTrendTest ] +def get_trend_test_name(trend_test_class): + years = list(range(10)) + trend_test = trend_test_class(years, years, None) + return trend_test.name + + +if __name__ == '__main__': + for trend_test_class in NON_STATIONARY_TREND_TEST_PAPER_2: + print(get_trend_test_name(trend_test_class)) -- GitLab