From e4b824178baa3d796cf676fa684b9c7a4eedf631 Mon Sep 17 00:00:00 2001 From: Le Roux Erwan <erwan.le-roux@irstea.fr> Date: Thu, 16 Jul 2020 15:47:48 +0200 Subject: [PATCH] [exceeding] modify plot settings for the paper. put "imap" instead of "map" to avoid computing 2 times the models --- extreme_data/eurocode_data/utils.py | 2 +- .../scm_models_data/abstract_study.py | 4 ++-- ...dy_visualizer_for_non_stationary_trends.py | 13 +++++++++++- extreme_trend/visualizers/utils.py | 5 +++-- .../main_snowfall_article.py | 2 +- .../section_data/main_eurocode_plot.py | 2 +- .../main_example_swe_total_plot.py | 2 +- .../main_result_trends_and_return_levels.py | 21 +++++++++++-------- .../plot_uncertainty_curves.py | 5 ++--- 9 files changed, 35 insertions(+), 21 deletions(-) diff --git a/extreme_data/eurocode_data/utils.py b/extreme_data/eurocode_data/utils.py index 3e5312ca..ed431a99 100644 --- a/extreme_data/eurocode_data/utils.py +++ b/extreme_data/eurocode_data/utils.py @@ -1,6 +1,6 @@ # Eurocode quantile str -EUROCODE_RETURN_LEVEL_STR = '50-year return level of GSL (kN $m^-2$)' +EUROCODE_RETURN_LEVEL_STR = '50-year return levels of GSL (kN $m^-2$)' # Eurocode quantile correspond to a 50 year return period EUROCODE_QUANTILE = 0.98 # Altitudes (between low and mid altitudes) < 2000m and should be > 200m diff --git a/extreme_data/meteo_france_data/scm_models_data/abstract_study.py b/extreme_data/meteo_france_data/scm_models_data/abstract_study.py index 43696208..c6212fef 100644 --- a/extreme_data/meteo_france_data/scm_models_data/abstract_study.py +++ b/extreme_data/meteo_france_data/scm_models_data/abstract_study.py @@ -589,8 +589,8 @@ class AbstractStudy(object): legend_elements = cls.get_legend_for_model_symbol(marker_style_to_label_name, markersize=7) ax.legend(handles=legend_elements[:], loc='upper right', prop={'size': 9}) # ax.legend(handles=legend_elements[4:], bbox_to_anchor=(0.01, 0.03), loc='lower left') - ax.annotate(' '.join(filled_marker_legend_list), - xy=(0.05, 0.015), xycoords='axes fraction', fontsize=7) + # ax.annotate(' '.join(filled_marker_legend_list), + # xy=(0.05, 0.015), xycoords='axes fraction', fontsize=7) if show: plt.show() diff --git a/extreme_trend/visualizers/study_visualizer_for_non_stationary_trends.py b/extreme_trend/visualizers/study_visualizer_for_non_stationary_trends.py index 36a4f22c..b07a5eca 100644 --- a/extreme_trend/visualizers/study_visualizer_for_non_stationary_trends.py +++ b/extreme_trend/visualizers/study_visualizer_for_non_stationary_trends.py @@ -58,12 +58,14 @@ class StudyVisualizerForNonStationaryTrends(StudyVisualizer): select_only_acceptable_shape_parameter=True, fit_gev_only_on_non_null_maxima=False, fit_only_time_series_with_ninety_percent_of_non_null_values=True, + select_only_model_that_pass_anderson_test=False, ): super().__init__(study, show, save_to_file, only_one_graph, only_first_row, vertical_kde_plot, year_for_kde_plot, plot_block_maxima_quantiles, temporal_non_stationarity, transformation_class, verbose, multiprocessing, complete_non_stationary_trend_analysis, normalization_under_one_observations) # Add some attributes + self.select_only_model_that_pass_anderson_test = select_only_model_that_pass_anderson_test self.fit_only_time_series_with_ninety_percent_of_non_null_values = fit_only_time_series_with_ninety_percent_of_non_null_values self.fit_gev_only_on_non_null_maxima = fit_gev_only_on_non_null_maxima self.select_only_acceptable_shape_parameter = select_only_acceptable_shape_parameter @@ -76,6 +78,8 @@ class StudyVisualizerForNonStationaryTrends(StudyVisualizer): self.uncertainty_massif_names = uncertainty_massif_names # Assign some default arguments if self.model_subsets_for_uncertainty is None: + # We regroup the results either with all the models, + # or only the stationary gumbel model which correspond to the building standards self.model_subsets_for_uncertainty = [ModelSubsetForUncertainty.stationary_gumbel, ModelSubsetForUncertainty.non_stationary_gumbel_and_gev][:] if self.uncertainty_methods is None: @@ -155,6 +159,7 @@ class StudyVisualizerForNonStationaryTrends(StudyVisualizer): @cached_property def massif_name_to_trend_test_tuple(self) -> Tuple[ Dict[str, AbstractGevTrendTest], Dict[str, AbstractGevTrendTest], Dict[str, AbstractGevTrendTest]]: + print('cached', self.altitude, id(self)) massif_name_to_trend_test_that_minimized_aic = {} massif_name_to_stationary_trend_test_that_minimized_aic = {} @@ -165,7 +170,13 @@ class StudyVisualizerForNonStationaryTrends(StudyVisualizer): # Extract the stationary or non-stationary model that minimized AIC trend_test_that_minimized_aic = sorted_trend_test[0] - massif_name_to_trend_test_that_minimized_aic[massif_name] = trend_test_that_minimized_aic + if self.select_only_model_that_pass_anderson_test and \ + (not trend_test_that_minimized_aic.goodness_of_fit_anderson_test): + print('not anderson') + else: + print('ok') + massif_name_to_trend_test_that_minimized_aic[massif_name] = trend_test_that_minimized_aic + # Extract the stationary model that minimized AIC stationary_trend_tests_that_minimized_aic = [t for t in sorted_trend_test if type(t) in [GumbelVersusGumbel, GevStationaryVersusGumbel]] diff --git a/extreme_trend/visualizers/utils.py b/extreme_trend/visualizers/utils.py index f056e22a..01c1b76b 100644 --- a/extreme_trend/visualizers/utils.py +++ b/extreme_trend/visualizers/utils.py @@ -24,7 +24,8 @@ def load_altitude_to_visualizer(altitudes, massif_names, model_subsets_for_uncer model_subsets_for_uncertainty=model_subsets_for_uncertainty, fit_method=fit_method, select_only_acceptable_shape_parameter=True, fit_gev_only_on_non_null_maxima=False, - fit_only_time_series_with_ninety_percent_of_non_null_values=True) + fit_only_time_series_with_ninety_percent_of_non_null_values=True, + select_only_model_that_pass_anderson_test=True, + ) altitude_to_visualizer[altitude] = study_visualizer - return altitude_to_visualizer 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 4fe9ff3a..7d5c0100 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 @@ -81,7 +81,7 @@ def intermediate_result(altitudes, massif_names=None, visualizers = list(altitude_to_visualizer.values()) if multiprocessing: with Pool(4) as p: - _ = p.map(compute_minimized_aic, visualizers) + _ = p.imap(compute_minimized_aic, visualizers) else: for visualizer in visualizers: _ = compute_minimized_aic(visualizer) diff --git a/projects/exceeding_snow_loads/section_data/main_eurocode_plot.py b/projects/exceeding_snow_loads/section_data/main_eurocode_plot.py index 0cdde265..850d0915 100644 --- a/projects/exceeding_snow_loads/section_data/main_eurocode_plot.py +++ b/projects/exceeding_snow_loads/section_data/main_eurocode_plot.py @@ -30,7 +30,7 @@ def main_eurocode_norms(ax=None, poster_plot=False): ax.set_yticks(list(range(0, 13, 2))) ax.set_xlim([0, 2000]) ax.set_xticks(list(range(0, 2001, 500))) - prop = {'size': 25} if poster_plot else {} + prop = {'size': 20} if poster_plot else {} ax.legend(prop=prop, loc='upper left') ax.grid() diff --git a/projects/exceeding_snow_loads/section_data/main_example_swe_total_plot.py b/projects/exceeding_snow_loads/section_data/main_example_swe_total_plot.py index 5719e15f..bb78ef8d 100644 --- a/projects/exceeding_snow_loads/section_data/main_example_swe_total_plot.py +++ b/projects/exceeding_snow_loads/section_data/main_example_swe_total_plot.py @@ -55,7 +55,7 @@ def max_graph_annual_maxima_poster(): last_plot = color == "magenta" label = '{} massif at {}m'.format(massif_name, altitude) tight_pad = {'h_pad': 0.2} - snow_abbreviation = 'max ' + snow_abbreviation + snow_abbreviation = 'annual maxima of ' + snow_abbreviation study_visualizer.visualize_max_graphs_poster(massif_name, altitude, snow_abbreviation, color, label, last_plot, ax, tight_pad=tight_pad, dpi=dpi_paper1_figure, diff --git a/projects/exceeding_snow_loads/section_results/main_result_trends_and_return_levels.py b/projects/exceeding_snow_loads/section_results/main_result_trends_and_return_levels.py index 61a2fe27..46f93e52 100644 --- a/projects/exceeding_snow_loads/section_results/main_result_trends_and_return_levels.py +++ b/projects/exceeding_snow_loads/section_results/main_result_trends_and_return_levels.py @@ -5,7 +5,8 @@ import matplotlib as mpl from projects.exceeding_snow_loads.checks.qqplot.plot_qqplot import \ plot_intensity_against_gumbel_quantile_for_3_examples, plot_full_diagnostic 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 +from projects.exceeding_snow_loads.section_results.plot_trend_curves import plot_trend_map, plot_trend_curves +from projects.exceeding_snow_loads.section_results.plot_uncertainty_histogram import plot_uncertainty_histogram mpl.rcParams['text.usetex'] = True mpl.rcParams['text.latex.preamble'] = [r'\usepackage{amsmath}'] @@ -33,6 +34,7 @@ def minor_result(altitude): def compute_minimized_aic(visualizer): + print(visualizer.altitude) _ = visualizer.massif_name_to_trend_test_that_minimized_aic return True @@ -58,22 +60,22 @@ def intermediate_result(altitudes, massif_names=None, for v in altitude_to_visualizer.values(): _ = v.study.year_to_variable_object # Compute minimized value efficiently - visualizers = list(altitude_to_visualizer.values()) + # visualizers = list() if multiprocessing: with Pool(NB_CORES) as p: - _ = p.map(compute_minimized_aic, visualizers) + _ = p.imap(compute_minimized_aic, altitude_to_visualizer.values()) else: - for visualizer in visualizers: + for visualizer in altitude_to_visualizer.values(): _ = compute_minimized_aic(visualizer) # Plots - # plot_trend_map(altitude_to_visualizer) - # plot_trend_curves(altitude_to_visualizer={a: v for a, v in altitude_to_visualizer.items() if a >= 900}) - # plot_uncertainty_massifs(altitude_to_visualizer) - # plot_uncertainty_histogram(altitude_to_visualizer) + plot_trend_map(altitude_to_visualizer) + plot_trend_curves(altitude_to_visualizer={a: v for a, v in altitude_to_visualizer.items() if a >= 900}) + plot_uncertainty_massifs(altitude_to_visualizer) + plot_uncertainty_histogram(altitude_to_visualizer) plot_selection_curves(altitude_to_visualizer) # uncertainty_interval_size(altitude_to_visualizer) - # plot_intensity_against_gumbel_quantile_for_3_examples(altitude_to_visualizer) + plot_intensity_against_gumbel_quantile_for_3_examples(altitude_to_visualizer) # plot_full_diagnostic(altitude_to_visualizer) @@ -89,6 +91,7 @@ def major_result(): # ModelSubsetForUncertainty.non_stationary_gumbel_and_gev] model_subsets_for_uncertainty = None # study_classes = [CrocusSnowLoad3Days, CrocusSnowLoad5Days, CrocusSnowLoad7Days][::-1] + altitudes = [300, 600] altitudes = paper_altitudes # altitudes = [900, 1800, 2700][:1] for study_class in study_classes: diff --git a/projects/exceeding_snow_loads/section_results/plot_uncertainty_curves.py b/projects/exceeding_snow_loads/section_results/plot_uncertainty_curves.py index 04403d6b..61ea9ccc 100644 --- a/projects/exceeding_snow_loads/section_results/plot_uncertainty_curves.py +++ b/projects/exceeding_snow_loads/section_results/plot_uncertainty_curves.py @@ -183,9 +183,8 @@ def plot_tdrl_bars(altitude_to_visualizer, ax, massif_name, valid_altitudes, leg # ax2.legend(handles=legend_elements, bbox_to_anchor=(0.93, 0.7), loc='upper right') # ax2.annotate("Filled symbol = significant trend ", xy=(0.85, 0.5), xycoords='axes fraction', fontsize=7) upper_legend_y = 0.55 - ax2.annotate('\n'.join(filled_marker_legend_list2), xy=(0.23, upper_legend_y - 0.2), xycoords='axes fraction', fontsize=fontsize) - ax2.annotate('Markers show selected model $\mathcal{M}_N$', xy=(0.02, upper_legend_y), xycoords='axes fraction', fontsize=fontsize) - print(legend_size) + # ax2.annotate('\n'.join(filled_marker_legend_list2), xy=(0.23, upper_legend_y - 0.2), xycoords='axes fraction', fontsize=fontsize) + # ax2.annotate('Markers show selected model $\mathcal{M}_N$', xy=(0.02, upper_legend_y), xycoords='axes fraction', fontsize=fontsize) ax2.legend(handles=legend_elements, prop={'size': legend_size}, loc='upper left', bbox_to_anchor=(0.00, upper_legend_y)) # for handle in lgnd.legendHandles: -- GitLab