diff --git a/projects/altitude_spatial_model/altitudes_fit/main_altitudes_studies.py b/projects/altitude_spatial_model/altitudes_fit/main_altitudes_studies.py index 60832aca40910268b742703d377e7a055c0a2d06..90277398940779d74f28c9e21fd385716631bc3a 100644 --- a/projects/altitude_spatial_model/altitudes_fit/main_altitudes_studies.py +++ b/projects/altitude_spatial_model/altitudes_fit/main_altitudes_studies.py @@ -11,7 +11,10 @@ from projects.altitude_spatial_model.altitudes_fit.one_fold_analysis.altitudes_s def plot_altitudinal_fit(studies): visualizer = AltitudesStudiesVisualizerForNonStationaryModels(studies=studies, model_classes=ALTITUDINAL_MODELS, - massif_names=['Belledonne']) + massif_names=['Belledonne'], + show=True) + visualizer.plot_mean() + visualizer.plot_relative_change() def plot_time_series(studies): diff --git a/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/altitudes_studies_visualizer_for_non_stationary_models.py b/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/altitudes_studies_visualizer_for_non_stationary_models.py index 9feb94a21b2512ca01c97769e1ec4b7e13044e70..0ea2b5e9bb366e6f74ffc26841e918b170f95d86 100644 --- a/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/altitudes_studies_visualizer_for_non_stationary_models.py +++ b/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/altitudes_studies_visualizer_for_non_stationary_models.py @@ -1,6 +1,12 @@ from typing import List +import matplotlib.pyplot as plt + +import numpy as np from extreme_data.meteo_france_data.scm_models_data.abstract_study import AbstractStudy +from extreme_data.meteo_france_data.scm_models_data.visualization.main_study_visualizer import \ + SCM_STUDY_CLASS_TO_ABBREVIATION +from extreme_data.meteo_france_data.scm_models_data.visualization.plot_utils import plot_against_altitude from extreme_data.meteo_france_data.scm_models_data.visualization.study_visualizer import StudyVisualizer from extreme_fit.model.margin_model.polynomial_margin_model.spatio_temporal_polynomial_model import \ AbstractSpatioTemporalPolynomialModel @@ -24,4 +30,31 @@ class AltitudesStudiesVisualizerForNonStationaryModels(StudyVisualizer): for massif_name in self.massif_names: dataset = studies.spatio_temporal_dataset(massif_name=massif_name) old_fold_fit = OneFoldFit(dataset, model_classes) - self.massif_name_to_one_fold_fit[massif_name] = old_fold_fit \ No newline at end of file + self.massif_name_to_one_fold_fit[massif_name] = old_fold_fit + + def plot_mean(self): + self.plot_general('mean') + + def plot_relative_change(self): + self.plot_general('relative_changes_in_the_mean') + + def plot_general(self, method_name): + ax = plt.gca() + min_altitude, *_, max_altitude = self.studies.altitudes + altitudes = np.linspace(min_altitude, max_altitude, num=50) + for massif_id, massif_name in enumerate(self.massif_names): + old_fold_fit = self.massif_name_to_one_fold_fit[massif_name] + values = old_fold_fit.__getattribute__(method_name)(altitudes) + plot_against_altitude(altitudes, ax, massif_id, massif_name, values) + # Plot settings + ax.legend(prop={'size': 7}, ncol=3) + moment = ' '.join(method_name.split('_')) + plot_name = '{} annual maxima of {}'.format(moment, SCM_STUDY_CLASS_TO_ABBREVIATION[self.studies.study_class]) + ax.set_ylabel('{} ({})'.format(plot_name, self.study.variable_unit), fontsize=15) + ax.set_xlabel('altitudes', fontsize=15) + # lim_down, lim_up = ax.get_ylim() + # lim_up += (lim_up - lim_down) / 3 + # ax.set_ylim([lim_down, lim_up]) + ax.tick_params(axis='both', which='major', labelsize=13) + self.studies.show_or_save_to_file(plot_name=plot_name, show=self.show) + ax.clear() diff --git a/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/one_fold_fit.py b/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/one_fold_fit.py index f1486e1037f4290cb077e2ac5f5871f59256c428..312caafb98d14da434a8beeeb3c2cbe5436276b3 100644 --- a/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/one_fold_fit.py +++ b/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/one_fold_fit.py @@ -1,3 +1,6 @@ +import numpy as np +from cached_property import cached_property + from extreme_fit.estimator.margin_estimator.utils import fitted_linear_margin_estimator_short from extreme_fit.model.margin_model.utils import MarginFitMethod @@ -18,3 +21,35 @@ class OneFoldFit(object): # Some display for estimator in self.model_class_to_estimator.values(): print(estimator.result_from_model_fit.aic) + + @cached_property + def best_estimator(self): + sorted_estimators = sorted([estimator for estimator in self.model_class_to_estimator.values()], + key=lambda e: e.result_from_model_fit.aic) + estimator_that_minimizes_aic = sorted_estimators[0] + return estimator_that_minimizes_aic + + @property + def best_function_from_fit(self): + return self.best_estimator.function_from_fit + + def mean(self, altitudes, year=2019): + return [self.get_mean(altitude, year) for altitude in altitudes] + + def get_mean(self, altitude, year): + coordinate = np.array([altitude, year]) + gev_params = self.best_function_from_fit.get_params(coordinate, is_transformed=False) + return gev_params.mean + + def relative_changes_in_the_mean(self, altitudes, year=2019, nb_years=50): + relative_changes = [] + for altitude in altitudes: + mean_after = self.get_mean(altitude, year) + mean_before = self.get_mean(altitude, year - nb_years) + relative_change = 100 * (mean_after - mean_before) / mean_before + relative_changes.append(relative_change) + return relative_changes + + + +