From 3e4d965cfc2d0de7bd4c81953e6c97cba72f8a63 Mon Sep 17 00:00:00 2001 From: Le Roux Erwan <erwan.le-roux@irstea.fr> Date: Wed, 15 May 2019 15:12:51 +0200 Subject: [PATCH] [SCM][MAXIMA TREND] add score based on the difference between top 10 argmax and top 10 argmin --- .../main_study_visualizer.py | 34 +++++++---- .../study_visualization/study_visualizer.py | 58 ++++++++++++------- 2 files changed, 61 insertions(+), 31 deletions(-) diff --git a/experiment/meteo_france_SCM_study/visualization/study_visualization/main_study_visualizer.py b/experiment/meteo_france_SCM_study/visualization/study_visualization/main_study_visualizer.py index e5259df5..ace9884b 100644 --- a/experiment/meteo_france_SCM_study/visualization/study_visualization/main_study_visualizer.py +++ b/experiment/meteo_france_SCM_study/visualization/study_visualization/main_study_visualizer.py @@ -17,8 +17,15 @@ SCM_STUDIES = [SafranSnowfall, CrocusSwe, CrocusDepth] SCM_EXTENDED_STUDIES = [ExtendedSafranSnowfall, ExtendedCrocusSwe, ExtendedCrocusDepth] SCM_STUDY_TO_EXTENDED_STUDY = OrderedDict(zip(SCM_STUDIES, SCM_EXTENDED_STUDIES)) +ALL_ALTITUDES = [0, 300, 600, 900, 1200, 1500, 1800, 2100, 2400, 2700, 3000, 3300, 3600, 3900, 4200, 4500, 4800] +full_altitude_with_at_least_2_stations = [0, 300, 600, 900, 1200, 1500, 1800, 2100, 2400, 2700, 3000, 3300, 3600, 3900, + 4200] -def study_iterator_global(study_classes, only_first_one=False, both_altitude=False, verbose=True, altitudes=None) -> List[AbstractStudy]: +ALL_STUDIES = SCM_STUDIES + [SafranTemperature, SafranRainfall] + + +def study_iterator_global(study_classes, only_first_one=False, both_altitude=False, verbose=True, altitudes=None) -> \ +List[AbstractStudy]: for study_class in study_classes: for study in study_iterator(study_class, only_first_one, both_altitude, verbose, altitudes): yield study @@ -26,7 +33,8 @@ def study_iterator_global(study_classes, only_first_one=False, both_altitude=Fal break -def study_iterator(study_class, only_first_one=False, both_altitude=False, verbose=True, altitudes=None) -> List[AbstractStudy]: +def study_iterator(study_class, only_first_one=False, both_altitude=False, verbose=True, altitudes=None) -> List[ + AbstractStudy]: all_studies = [] is_safran_study = study_class in [SafranSnowfall, ExtendedSafranSnowfall] nb_days = [1] if is_safran_study else [1] @@ -93,15 +101,21 @@ def normal_visualization(temporal_non_stationarity=False): only_first_one = True # for study_class in SCM_STUDIES[:1]: for study_class in [CrocusDepth, SafranSnowfall, SafranRainfall, SafranTemperature][1:2]: - for study in study_iterator(study_class, only_first_one=only_first_one): + for study in study_iterator(study_class, only_first_one=only_first_one, altitudes=[300]): study_visualizer = StudyVisualizer(study, save_to_file=save_to_file, temporal_non_stationarity=temporal_non_stationarity) + study_visualizer.window_size_for_smoothing = 1 + study_visualizer.number_of_top_values = 10 + study_visualizer.visualize_all_mean_and_max_graphs() + # study_visualizer.visualize_independent_margin_fits(threshold=[None, 20, 40, 60][0]) # study_visualizer.visualize_annual_mean_values() - study_visualizer.visualize_all_mean_and_max_graphs() - # study_visualizer.visualize_linear_margin_fit(only_first_max_stable=None) +def all_normal_vizu(): + for study in study_iterator_global(study_classes=ALL_STUDIES, only_first_one=False, altitudes=ALL_ALTITUDES): + study_visualizer = StudyVisualizer(study, save_to_file=True, temporal_non_stationarity=True) + study_visualizer.visualize_all_mean_and_max_graphs() def complete_analysis(only_first_one=False): """An overview of everything that is possible with study OR extended study""" @@ -122,7 +136,8 @@ def trend_analysis(): save_to_file = False only_first_one = True short_altitudes = [300, 1200, 2100, 3000][:1] - full_altitude_with_at_least_2_stations = [0, 300, 600, 900, 1200, 1500, 1800, 2100, 2400, 2700, 3000, 3300, 3600, 3900, 4200][:] + full_altitude_with_at_least_2_stations = [0, 300, 600, 900, 1200, 1500, 1800, 2100, 2400, 2700, 3000, 3300, 3600, + 3900, 4200][:] durand_altitude = [1800] altitudes = durand_altitude normalization_class = [None, BetweenMinusOneAndOneNormalization, BetweenZeroAndOneNormalization][-1] @@ -141,17 +156,16 @@ def trend_analysis(): def main_run(): # normal_visualization(temporal_non_stationarity=True) - trend_analysis() - - + # trend_analysis() + all_normal_vizu() # annual_mean_vizu_compare_durand_study(safran=True, take_mean_value=True, altitude=2100) - # max_stable_process_vizu_compare_gaume_study(altitude=1800, nb_days=1) # extended_visualization() # complete_analysis() + if __name__ == '__main__': start = time.time() main_run() diff --git a/experiment/meteo_france_SCM_study/visualization/study_visualization/study_visualizer.py b/experiment/meteo_france_SCM_study/visualization/study_visualization/study_visualizer.py index 4f1eb446..af6b1a6e 100644 --- a/experiment/meteo_france_SCM_study/visualization/study_visualization/study_visualizer.py +++ b/experiment/meteo_france_SCM_study/visualization/study_visualization/study_visualizer.py @@ -27,6 +27,7 @@ from extreme_estimator.extreme_models.max_stable_model.abstract_max_stable_model from extreme_estimator.margin_fits.abstract_params import AbstractParams from extreme_estimator.margin_fits.gev.gev_params import GevParams from extreme_estimator.margin_fits.gev.gevmle_fit import GevMleFit +from extreme_estimator.margin_fits.gev.ismev_gev_fit import IsmevGevFit from extreme_estimator.margin_fits.gpd.gpd_params import GpdParams from extreme_estimator.margin_fits.gpd.gpdmle_fit import GpdMleFit from spatio_temporal_dataset.coordinates.spatial_coordinates.abstract_spatial_coordinates import \ @@ -41,7 +42,8 @@ from spatio_temporal_dataset.spatio_temporal_observations.abstract_spatio_tempor AbstractSpatioTemporalObservations from spatio_temporal_dataset.spatio_temporal_observations.annual_maxima_observations import AnnualMaxima from test.test_utils import load_test_max_stable_models -from utils import get_display_name_from_object_type, VERSION_TIME, float_to_str_with_only_some_significant_digits +from utils import get_display_name_from_object_type, VERSION_TIME, float_to_str_with_only_some_significant_digits, \ + cached_property BLOCK_MAXIMA_DISPLAY_NAME = 'block maxima ' @@ -83,7 +85,8 @@ class StudyVisualizer(object): self.year_for_kde_plot = year_for_kde_plot self.plot_block_maxima_quantiles = plot_block_maxima_quantiles - self.window_size_for_smoothing = 21 + self.window_size_for_smoothing = 1 # other value could be + self.number_of_top_values = 10 # 1 if we just want the maxima # PLOT ARGUMENTS self.show = False if self.save_to_file else show @@ -163,7 +166,7 @@ class StudyVisualizer(object): # PLOT FOR SEVERAL MASSIFS - def visualize_massif_graphs(self, visualize_function, specified_massif_names=None): + def visualize_massif_graphs(self, visualize_function, use_ordered_massif_names=False): if self.only_one_graph: fig, ax = plt.subplots(1, 1, figsize=self.figsize) visualize_function(ax, 0) @@ -177,11 +180,11 @@ class StudyVisualizer(object): ax = axes[massif_id] visualize_function(ax, massif_id) else: - if specified_massif_names is None: - massif_ids = list(range(len(self.study.study_massif_names))) - else: + if use_ordered_massif_names: massif_ids = [self.study.study_massif_names.index(massif_name) for massif_name in - specified_massif_names] + self.ordered_massif_names] + else: + massif_ids = list(range(len(self.study.study_massif_names))) for j, massif_id in enumerate(massif_ids): row_id, column_id = j // nb_columns, j % nb_columns ax = axes[row_id, column_id] @@ -324,18 +327,29 @@ class StudyVisualizer(object): all_massif_data = np.sort(all_massif_data) return all_massif_data - # - - def visualize_all_mean_and_max_graphs(self): - # Compute the order of massif names + @cached_property + def massif_name_to_score(self): + # Ordered massif by scores massif_name_to_score = {} for massif_id, massif_name in enumerate(self.study.study_massif_names): - score = self.smooth_maxima_x_y(massif_id)[1].argmax() - massif_name_to_score[massif_name] = score - ordered_massif_names = sorted(self.study.study_massif_names[:], key=lambda s: massif_name_to_score[s]) - self.visualize_massif_graphs(self.visualize_mean_and_max_graph, specified_massif_names=ordered_massif_names) - self.plot_name = ' smoothing values temporally with sliding window of size {}'.format( - self.window_size_for_smoothing) + years, smooth_maxima = self.smooth_maxima_x_y(massif_id) + start_year = years[0] + sorted_indices = [i for i, e in sorted(enumerate(smooth_maxima), key=lambda s: s[1])] + mean_max_year = np.mean(sorted_indices[-self.number_of_top_values:]) + start_year + mean_min_years = np.mean(sorted_indices[:self.number_of_top_values]) + start_year + score = mean_max_year - mean_min_years + massif_name_to_score[massif_name] = (score, mean_max_year, mean_min_years) + return massif_name_to_score + + @cached_property + def ordered_massif_names(self): + return sorted(self.study.study_massif_names[:], key=lambda s: self.massif_name_to_score[s][0]) + + def visualize_all_mean_and_max_graphs(self): + self.visualize_massif_graphs(self.visualize_mean_and_max_graph, use_ordered_massif_names=True) + self.plot_name = ' smoothing values temporally with sliding window of size {}' \ + 'and {} top values taking into account'.format( + self.window_size_for_smoothing, self.number_of_top_values) self.show_or_save_to_file() def visualize_mean_and_max_graph(self, ax, massif_id): @@ -355,9 +369,11 @@ class StudyVisualizer(object): x, y = average_smoothing_with_sliding_window(x, y, window_size_for_smoothing=self.window_size_for_smoothing) ax.plot(x, y, color=color_mean) ax.set_ylabel('mean'.format(self.window_size_for_smoothing), color=color_mean) - ax.set_xlabel('year') - ax.set_title(self.study.study_massif_names[massif_id]) - ax.xaxis.set_ticks(x[2::10]) + massif_name = self.study.study_massif_names[massif_id] + title = massif_name + title += ' {}={}-{}'.format(*[round(e, 1) for e in list(self.massif_name_to_score[massif_name])]) + ax.set_title(title) + ax.xaxis.set_ticks(x[2::20]) def smooth_maxima_x_y(self, massif_id): if massif_id not in self.massif_id_to_smooth_maxima: @@ -536,7 +552,7 @@ class StudyVisualizer(object): @property def df_gev_mle(self) -> pd.DataFrame: # Fit a margin_fits on each massif - massif_to_gev_mle = {massif_name: GevMleFit(self.df_maxima_gev.loc[massif_name]).gev_params.summary_serie + massif_to_gev_mle = {massif_name: IsmevGevFit(self.df_maxima_gev.loc[massif_name]).gev_params.summary_serie for massif_name in self.study.study_massif_names} return pd.DataFrame(massif_to_gev_mle, columns=self.study.study_massif_names) -- GitLab