From 08ca4b8843fc9ddfe55d8442017ddca42b2ad736 Mon Sep 17 00:00:00 2001
From: Le Roux Erwan <erwan.le-roux@irstea.fr>
Date: Thu, 25 Feb 2021 11:29:21 +0100
Subject: [PATCH] [refactor] remove anything related to total_aic for the
 second paper.

---
 .../altitudes_fit/main_altitudes_studies.py   | 31 +++----
 ...es_visualizer_for_non_stationary_models.py |  9 +-
 .../one_fold_analysis/one_fold_fit.py         | 77 ++--------------
 .../one_fold_analysis/plot_total_aic.py       | 90 -------------------
 .../altitudes_fit/plot_mean_and_std.py        |  8 --
 .../__init__.py                               |  0
 .../main_gap_altitudes_studies.py             |  0
 7 files changed, 25 insertions(+), 190 deletions(-)
 delete mode 100644 projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/plot_total_aic.py
 delete mode 100644 projects/altitude_spatial_model/altitudes_fit/plot_mean_and_std.py
 create mode 100644 projects/archive/gap_between_my_safran2019_and_safran_2019/__init__.py
 rename projects/{altitude_spatial_model/altitudes_fit => archive/gap_between_my_safran2019_and_safran_2019}/main_gap_altitudes_studies.py (100%)

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 e938e386..2209091f 100644
--- a/projects/altitude_spatial_model/altitudes_fit/main_altitudes_studies.py
+++ b/projects/altitude_spatial_model/altitudes_fit/main_altitudes_studies.py
@@ -2,35 +2,34 @@ import datetime
 import time
 from typing import List
 
-import matplotlib
 
+import matplotlib as mpl
+mpl.rcParams['text.usetex'] = True
+mpl.rcParams['text.latex.preamble'] = [r'\usepackage{amsmath}']
+
+import matplotlib
 matplotlib.use('Agg')
 
-from extreme_data.meteo_france_data.scm_models_data.safran.safran_max_snowf import SafranSnowfall2019, \
-    SafranSnowfall2020
 from projects.altitude_spatial_model.altitudes_fit.plots.plot_histogram_altitude_studies import \
     plot_shoe_plot_changes_against_altitude, plot_histogram_all_trends_against_altitudes, \
-    plot_shoe_plot_ratio_interval_size_against_altitude, plot_histogram_all_models_against_altitudes
+    plot_histogram_all_models_against_altitudes
 
 from extreme_fit.model.result_from_model_fit.result_from_extremes.abstract_extract_eurocode_return_level import \
     AbstractExtractEurocodeReturnLevel
 
-import matplotlib as mpl
+
 
 from extreme_fit.model.utils import set_seed_for_test
 from projects.altitude_spatial_model.altitudes_fit.plots.plot_coherence_curves import plot_coherence_curves
 
-mpl.rcParams['text.usetex'] = True
-mpl.rcParams['text.latex.preamble'] = [r'\usepackage{amsmath}']
+
 
 from projects.altitude_spatial_model.altitudes_fit.utils_altitude_studies_visualizer import load_visualizer_list
 
 from projects.altitude_spatial_model.altitudes_fit.one_fold_analysis.altitude_group import altitudes_for_groups
-from projects.altitude_spatial_model.altitudes_fit.one_fold_analysis.plot_total_aic import plot_individual_aic
 
 from extreme_data.meteo_france_data.scm_models_data.safran.safran import SafranSnowfall1Day, SafranSnowfall3Days, \
-    SafranSnowfall5Days, SafranSnowfall7Days, SafranDateFirstSnowfall, SafranPrecipitation1Day, \
-    SafranPrecipitation3Days, SafranSnowfallCenterOnDay1dayMeanRate, SafranSnowfallNotCenterOnDay1day
+    SafranSnowfall5Days, SafranSnowfall7Days
 from extreme_data.meteo_france_data.scm_models_data.utils import Season
 
 
@@ -96,13 +95,11 @@ def plot_visualizer(massif_names, visualizer):
     # visualizer.studies.plot_maxima_time_series(massif_names)
     # visualizer.studies.plot_maxima_time_series(['Vanoise'])
 
-    # Plot the results for the model that minimizes the individual aic
-    plot_individual_aic(visualizer)
-
-    # Plot the results for the model that minimizes the total aic
-    # plot_total_aic(model_classes, visualizer)
-    pass
-
+    visualizer.plot_shape_map()
+    visualizer.plot_moments()
+    visualizer.plot_qqplots()
+    # for std in [True, False]:
+    #     visualizer.studies.plot_mean_maxima_against_altitude(std=std)
 
 if __name__ == '__main__':
     main()
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 9c3b7f9f..ffe76956 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
@@ -191,7 +191,7 @@ class AltitudesStudiesVisualizerForNonStationaryModels(StudyVisualizer):
             if self.temporal_covariate_for_fit is AnomalyTemperatureTemporalCovariate else ' in {}'
         str_for_last_year = str_for_last_year.format(self.first_one_fold_fit.covariate_after, **d_temperature)
         moment = moment.replace('moment', '{}{}'.format(OneFoldFit.get_moment_str(order=order), str_for_last_year))
-        plot_name = '{}{} '.format(OneFoldFit.folder_for_plots, moment)
+        plot_name = '{} '.format(moment)
 
         if 'change' in method_name:
             plot_name = plot_name.replace(str_for_last_year, '')
@@ -264,7 +264,7 @@ class AltitudesStudiesVisualizerForNonStationaryModels(StudyVisualizer):
         moment = ' '.join(method_name.split('_'))
         moment = moment.replace('moment',
                                 '{} in {}'.format(OneFoldFit.get_moment_str(order=order), OneFoldFit.last_year))
-        plot_name = '{}Model {} annual maxima of {}'.format(OneFoldFit.folder_for_plots, moment,
+        plot_name = 'Model {} 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)
@@ -321,8 +321,7 @@ class AltitudesStudiesVisualizerForNonStationaryModels(StudyVisualizer):
         negative_and_positive_values = (max(values) > 0) and (min(values) < 0)
         raise NotImplementedError
         self.plot_map(massif_name_to_value=massif_name_to_best_coef,
-                      label='{}Coef/{} plot for {} {}'.format(OneFoldFit.folder_for_plots,
-                                                              coef_name,
+                      label='Coef/{} plot for {} {}'.format(coef_name,
                                                               SCM_STUDY_CLASS_TO_ABBREVIATION[
                                                                   type(self.study)],
                                                               self.study.variable_unit),
@@ -387,7 +386,7 @@ class AltitudesStudiesVisualizerForNonStationaryModels(StudyVisualizer):
             ax.set_ylabel('{} of {} in {} ({})'.format(ylabel, SCM_STUDY_CLASS_TO_ABBREVIATION[type(self.study)],
                                                        massif_name.replace('_', ' '), self.study.variable_unit))
             peak_year_folder = 'Peak year ' + ylabel
-            plot_name = '{}{}/Peak year for {}'.format(OneFoldFit.folder_for_plots, peak_year_folder,
+            plot_name = '{}/Peak year for {}'.format(peak_year_folder,
                                                        massif_name.replace('_', ''))
             self.studies.show_or_save_to_file(plot_name=plot_name, show=self.show, no_title=True, tight_layout=True)
             plt.close()
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 a0451c44..525e969e 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
@@ -43,7 +43,6 @@ from spatio_temporal_dataset.spatio_temporal_observations.annual_maxima_observat
 
 class OneFoldFit(object):
     SIGNIFICANCE_LEVEL = 0.05
-    best_estimator_minimizes_total_aic = False
     return_period = 100
     quantile_level = 1 - (1 / return_period)
     nb_years = 60
@@ -75,21 +74,12 @@ class OneFoldFit(object):
         # Compute sorted estimators indirectly
         _ = self.has_at_least_one_valid_model
 
-        # Best estimator definition
-        self.best_estimator_class_for_total_aic = None
-        # Cached object
-        self._folder_to_goodness = {}
-
     def fitted_linear_margin_estimator(self, model_class, dataset):
         return fitted_linear_margin_estimator_short(model_class=model_class,
                                                     dataset=dataset,
                                                     fit_method=self.fit_method,
                                                     temporal_covariate_for_fit=self.temporal_covariate_for_fit)
 
-    @classproperty
-    def folder_for_plots(cls):
-        return 'Total aic/' if cls.best_estimator_minimizes_total_aic else ''
-
     @classmethod
     def get_moment_str(cls, order):
         if order == 1:
@@ -159,7 +149,7 @@ class OneFoldFit(object):
         s = ' between +${}^o\mathrm{C}$ and +${}^o\mathrm{C}$' if self.temporal_covariate_for_fit is AnomalyTemperatureTemporalCovariate \
             else ' between {} and {}'
         s = s.format(self.covariate_before, self.covariate_after,
-                                 **d_temperature)
+                     **d_temperature)
         return s
 
     def relative_changes_of_moment(self, altitudes, order=1):
@@ -184,7 +174,7 @@ class OneFoldFit(object):
             for e in estimators:
                 coordinate_values_for_the_fit = e.coordinates_for_nllh(Split.all)
                 coordinate_values_for_the_result = [np.array([self.altitude_group.reference_altitude, c])
-                                                             for c in self._covariate_before_and_after]
+                                                    for c in self._covariate_before_and_after]
                 coordinate_values_to_check = list(coordinate_values_for_the_fit) + coordinate_values_for_the_result
                 has_undefined_parameters = False
                 for coordinate in coordinate_values_to_check:
@@ -218,8 +208,6 @@ class OneFoldFit(object):
         if self.only_models_that_pass_goodness_of_fit_test:
             return [e for e in self.sorted_estimators
                     if self.goodness_of_fit_test(e)
-                    # and self.sensitivity_of_fit_test_top_maxima(e)
-                    # and self.sensitivity_of_fit_test_last_years(e)
                     ]
         else:
             if not self.remove_physically_implausible_models:
@@ -236,16 +224,12 @@ class OneFoldFit(object):
 
     @property
     def best_estimator(self):
-        if self.best_estimator_minimizes_total_aic and self.best_estimator_class_for_total_aic is not None:
-            return self.model_class_to_estimator[self.best_estimator_class_for_total_aic]
+        if self.has_at_least_one_valid_model:
+            best_estimator = self.sorted_estimators_with_stationary[0]
+            return best_estimator
         else:
-            # With stationary
-            if self.has_at_least_one_valid_model:
-                best_estimator = self.sorted_estimators_with_stationary[0]
-                return best_estimator
-            else:
-                raise ValueError('This object should not have been called because '
-                                 'has_at_least_one_valid_model={}'.format(self.has_at_least_one_valid_model))
+            raise ValueError('This object should not have been called because '
+                             'has_at_least_one_valid_model={}'.format(self.has_at_least_one_valid_model))
 
     @property
     def best_margin_model(self):
@@ -322,52 +306,6 @@ class OneFoldFit(object):
             # Bootstrap based significance
             return self.cached_results_from_bootstrap[0]
 
-    # @property
-    # def goodness_of_fit_anderson_test(self):
-    #     if self.folder_for_plots in self._folder_to_goodness:
-    #         return self._folder_to_goodness[self.folder_for_plots]
-    #     else:
-    #         estimator = self.best_estimator
-    #         goodness_of_fit_anderson_test = self.goodness_of_fit_test(estimator)
-    #         if not goodness_of_fit_anderson_test:
-    #             print('{} with {} does not pass the anderson test for model {}'.format(self.massif_name,
-    #                                                                                    self.folder_for_plots,
-    #                                                                                    type(self.best_margin_model)))
-    #         self._folder_to_goodness[self.folder_for_plots] = goodness_of_fit_anderson_test
-    #         return goodness_of_fit_anderson_test
-
-    def sensitivity_of_fit_test_top_maxima(self, estimator: LinearMarginEstimator):
-        # Build the dataset without the maxima for each altitude
-        new_dataset = AbstractDataset.remove_top_maxima(self.dataset.observations,
-                                                        self.dataset.coordinates)
-        # Fit the new estimator
-        new_estimator = fitted_linear_margin_estimator_short(model_class=type(estimator.margin_model),
-                                                             dataset=new_dataset,
-                                                             fit_method=self.fit_method,
-                                                             temporal_covariate_for_fit=self.temporal_covariate_for_fit)
-        # Compare sign of change
-        has_not_opposite_sign = self.sign_of_change(estimator.function_from_fit) * self.sign_of_change(
-            new_estimator.function_from_fit) >= 0
-        return has_not_opposite_sign
-
-    def sensitivity_of_fit_test_last_years(self, estimator: LinearMarginEstimator):
-        # Build the dataset without the maxima for each altitude
-        new_dataset = AbstractDataset.remove_last_maxima(self.dataset.observations,
-                                                         self.dataset.coordinates,
-                                                         nb_years=10)
-        # Fit the new estimator
-        model_class = type(estimator.margin_model)
-        new_estimator = fitted_linear_margin_estimator_short(model_class=model_class,
-                                                             dataset=new_dataset,
-                                                             fit_method=self.fit_method,
-                                                             temporal_covariate_for_fit=self.temporal_covariate_for_fit)
-        # Compare sign of change
-        has_not_opposite_sign = self.sign_of_change(estimator.function_from_fit) * self.sign_of_change(
-            new_estimator.function_from_fit) >= 0
-        # if not has_not_opposite_sign:
-        # print('Last years', self.massif_name, model_class, self.sign_of_change(estimator), self.sign_of_change(new_estimator))
-        return has_not_opposite_sign
-
     def sign_of_change(self, function_from_fit):
         return_levels = []
         for temporal_covariate in self._covariate_before_and_after:
@@ -445,7 +383,6 @@ class OneFoldFit(object):
 
         return is_significant, altitude_and_year_to_return_level_mean_estimate, altitude_and_year_to_return_level_confidence_interval
 
-
     @property
     def bootstrap_fitted_functions_from_fit(self):
         print('nb of bootstrap for confidence interval=', AbstractExtractEurocodeReturnLevel.NB_BOOTSTRAP)
diff --git a/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/plot_total_aic.py b/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/plot_total_aic.py
deleted file mode 100644
index 2a50df32..00000000
--- a/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/plot_total_aic.py
+++ /dev/null
@@ -1,90 +0,0 @@
-import matplotlib as mpl
-import matplotlib.pyplot as plt
-import numpy as np
-
-from projects.altitude_spatial_model.altitudes_fit.one_fold_analysis.altitudes_studies_visualizer_for_non_stationary_models import \
-    AltitudesStudiesVisualizerForNonStationaryModels
-
-mpl.rcParams['text.usetex'] = True
-mpl.rcParams['text.latex.preamble'] = [r'\usepackage{amsmath}']
-
-from extreme_data.meteo_france_data.scm_models_data.visualization.utils import create_adjusted_axes
-from projects.altitude_spatial_model.altitudes_fit.one_fold_analysis.one_fold_fit import OneFoldFit
-from projects.exceeding_snow_loads.utils import dpi_paper1_figure
-
-
-def plots(visualizer: AltitudesStudiesVisualizerForNonStationaryModels):
-    visualizer.plot_shape_map()
-    visualizer.plot_moments()
-    visualizer.plot_qqplots()
-    # for std in [True, False]:
-    #     visualizer.studies.plot_mean_maxima_against_altitude(std=std)
-
-
-def plot_individual_aic(visualizer):
-    OneFoldFit.best_estimator_minimizes_total_aic = False
-    plots(visualizer)
-
-
-def plot_total_aic(model_classes, visualizer):
-    # Compute the mean AIC for each model_class
-    OneFoldFit.best_estimator_minimizes_total_aic = True
-    model_class_to_total_aic = {model_class: 0 for model_class in model_classes}
-    model_class_to_name_str = {}
-    # model_class_to_aic_scores = {model_class: [] for model_class in model_classes}
-    # model_class_to_n = {model_class: [] for model_class in model_classes}
-    for one_fold_fit in visualizer.massif_name_to_one_fold_fit.values():
-        for model_class, estimator in one_fold_fit.model_class_to_estimator_with_finite_aic.items():
-            model_class_to_total_aic[model_class] += estimator.aic()
-            model_class_to_name_str[model_class] = estimator.margin_model.name_str
-            # model_class_to_aic_scores[model_class].append(estimator.aic())
-            # model_class_to_n[model_class].append(estimator.n())
-    # model_class_to_mean_aic_score = {model_class: np.array(aic_scores).sum() / np.array(model_class_to_n[model_class]).sum()
-    #                                  for model_class, aic_scores in model_class_to_aic_scores.items()}
-    # print(model_class_to_mean_aic_score)
-    sorted_model_class = sorted(model_classes, key=lambda m: model_class_to_total_aic[m])
-    sorted_scores = [model_class_to_total_aic[model_class] for model_class in sorted_model_class]
-    sorted_labels = [model_class_to_name_str[model_class] for model_class in sorted_model_class]
-    print(sorted_model_class)
-    print(sorted_scores)
-    print(sorted_labels)
-    best_model_class_for_total_aic = sorted_model_class[0]
-    for one_fold_fit in visualizer.massif_name_to_one_fold_fit.values():
-        one_fold_fit.best_estimator_class_for_total_aic = best_model_class_for_total_aic
-    # Plot the ranking of the model based on their total aic
-    plot_total_aic_repartition(visualizer, sorted_labels, sorted_scores)
-    # Plot the results for the model that minimizes the mean aic
-    plots(visualizer)
-
-
-def plot_total_aic_repartition(visualizer, labels, scores):
-    """
-    Plot a single trend curves
-    :return:
-    """
-    scores = np.array(scores)
-    ax = create_adjusted_axes(1, 1)
-
-    # parameters
-    width = 3
-    size = 30
-    legend_fontsize = 30
-    labelsize = 10
-    linewidth = 3
-    x = [2 * width * (i + 1) for i in range(len(labels))]
-    ax.bar(x, scores, width=width, color='grey', edgecolor='grey', label='Total Aic',
-           linewidth=linewidth)
-    ax.legend(loc='upper right', prop={'size': size})
-    ax.set_ylabel(' Total AIC score \n '
-                  'i.e. sum of AIC score for all massifs ', fontsize=legend_fontsize)
-    ax.set_xlabel('Models', fontsize=legend_fontsize)
-    ax.tick_params(axis='both', which='major', labelsize=labelsize)
-    ax.set_xticks(x)
-    ax.set_ylim([scores.min(), scores.max()])
-    ax.yaxis.grid()
-    ax.set_xticklabels(labels)
-
-    # Save plot
-    visualizer.plot_name = 'Total AIC ranking'
-    visualizer.show_or_save_to_file(no_title=True, dpi=dpi_paper1_figure)
-    plt.close()
diff --git a/projects/altitude_spatial_model/altitudes_fit/plot_mean_and_std.py b/projects/altitude_spatial_model/altitudes_fit/plot_mean_and_std.py
deleted file mode 100644
index 20bc2fc7..00000000
--- a/projects/altitude_spatial_model/altitudes_fit/plot_mean_and_std.py
+++ /dev/null
@@ -1,8 +0,0 @@
-from extreme_data.meteo_france_data.scm_models_data.safran.safran import SafranSnowfall1Day
-from projects.altitude_spatial_model.altitudes_fit.altitudes_studies import AltitudesStudies
-
-if __name__ == '__main__':
-    studies = AltitudesStudies(study_class=SafranSnowfall1Day,
-                     altitudes=[600, 900, 1200, 1500, 1800, 2100, 2400, 2700, 3000, 3300, 3600])
-    for std in [True, False]:
-        studies.plot_mean_maxima_against_altitude(std=std)
\ No newline at end of file
diff --git a/projects/archive/gap_between_my_safran2019_and_safran_2019/__init__.py b/projects/archive/gap_between_my_safran2019_and_safran_2019/__init__.py
new file mode 100644
index 00000000..e69de29b
diff --git a/projects/altitude_spatial_model/altitudes_fit/main_gap_altitudes_studies.py b/projects/archive/gap_between_my_safran2019_and_safran_2019/main_gap_altitudes_studies.py
similarity index 100%
rename from projects/altitude_spatial_model/altitudes_fit/main_gap_altitudes_studies.py
rename to projects/archive/gap_between_my_safran2019_and_safran_2019/main_gap_altitudes_studies.py
-- 
GitLab