[PAPER 1] set starting year to None to avoid bug due to missing data for the year 1958 for instance

experiment/paper_past_snow_loads/result_trends_and_return_levels/eurocode_visualizer.py

@@ -37,8 +37,8 @@ def plot_uncertainty_massifs(altitude_to_visualizer: Dict[int, StudyVisualizerFo
     :return:
     """
     visualizer = list(altitude_to_visualizer.values())[0]
-    # Subdivide massif names in group of 5
-    m = 4
+    # Subdivide massif names in group of 3
+    m = 3
     uncertainty_massif_names = visualizer.uncertainty_massif_names
     n = (len(uncertainty_massif_names) // m) + 1
     for i in list(range(n))[:]:
@@ -140,6 +140,7 @@ def plot_single_uncertainty_massif_and_non_stationary_context(ax, massif_name, n
     title = '{} ({} Eurocodes area) with a {} model'.format(massif_name_str, eurocode_region_str,
                                                             non_stationary_context)
     ax.set_title(title)
+    ax.set_xticks(altitudes)
     ax.set_ylabel('50-year return level of SL (kN $m^-2$)')
     ax.set_xlabel('Altitude (m)')
     ax.grid()

experiment/paper_past_snow_loads/result_trends_and_return_levels/main_result_trends_and_return_levels.py

@@ -67,9 +67,16 @@ def major_result():
 
 if __name__ == '__main__':
     # minor_result(altitude=1800)
-    intermediate_result(altitudes=[1500, 1800], massif_names=['Chartreuse'],
-                        uncertainty_methods=[ConfidenceIntervalMethodFromExtremes.ci_mle],
-                        non_stationary_uncertainty=[False])
+    # intermediate_result(altitudes=[1500, 1800], massif_names=['Chartreuse'],
+    #                     uncertainty_methods=[ConfidenceIntervalMethodFromExtremes.ci_mle],
+    #                     non_stationary_uncertainty=[False])
+    # intermediate_result(altitudes=[1500, 1800], massif_names=None,
+    #                     uncertainty_methods=[ConfidenceIntervalMethodFromExtremes.ci_mle],
+    #                     non_stationary_uncertainty=[False])
+    # intermediate_result(altitudes=[300, 600, 900, 1200, 1500, 1800], massif_names=None,
+    #                     uncertainty_methods=[ConfidenceIntervalMethodFromExtremes.ci_mle],
+    #                     non_stationary_uncertainty=[False])
     intermediate_result(altitudes=[300, 600, 900, 1200, 1500, 1800], massif_names=None,
-                        uncertainty_methods=[ConfidenceIntervalMethodFromExtremes.ci_mle],
-                        non_stationary_uncertainty=[False])
+                        uncertainty_methods=[ConfidenceIntervalMethodFromExtremes.ci_mle,
+                                             ConfidenceIntervalMethodFromExtremes.ci_bayes],
+                        non_stationary_uncertainty=[False, True])

experiment/paper_past_snow_loads/result_trends_and_return_levels/study_visualizer_for_non_stationary_trends.py

@@ -6,6 +6,7 @@ from typing import Dict, Tuple
 import numpy as np
 from cached_property import cached_property
 
+from experiment.eurocode_data.utils import EUROCODE_QUANTILE
 from experiment.meteo_france_data.scm_models_data.abstract_study import AbstractStudy
 from experiment.meteo_france_data.scm_models_data.visualization.study_visualization.study_visualizer import \
     StudyVisualizer
@@ -70,6 +71,12 @@ class StudyVisualizerForNonStationaryTrends(StudyVisualizer):
         return {m: np.count_nonzero(maxima) / len(maxima) for m, (_, maxima) in
                 self.massif_name_to_years_and_maxima.items()}
 
+    @cached_property
+    def massif_name_to_eurocode_quantile_in_practice(self):
+        """Due to missing data, the the eurocode quantile which 0.98 if we have all the data
+        correspond in practice to the quantile psnow x 0.98 of the data where there is snow"""
+        return {m: p * EUROCODE_QUANTILE for m, p in self.massif_name_to_psnow.items()}
+
     @cached_property
     def massif_name_to_non_null_years_and_maxima(self):
         d = {}
@@ -80,7 +87,7 @@ class StudyVisualizerForNonStationaryTrends(StudyVisualizer):
 
     @cached_property
     def massif_name_to_minimized_aic_non_stationary_trend_test(self) -> Dict[str, AbstractGevTrendTest]:
-        starting_year = 1958
+        starting_year = None
         massif_name_to_trend_test_that_minimized_aic = {}
         for massif_name, (x, y) in self.massif_name_to_non_null_years_and_maxima.items():
             non_stationary_trend_test = [t(x, y, starting_year) for t in self.non_stationary_trend_test]

extreme_fit/model/result_from_model_fit/result_from_extremes/eurocode_return_level_uncertainties.py

@@ -49,7 +49,7 @@ class EurocodeConfidenceIntervalFromExtremes(object):
         else:
             fit_method = TemporalMarginFitMethod.extremes_fevd_mle
         # Fitted estimator
-        fitted_estimator = fitted_linear_margin_estimator(model_class, coordinates, dataset, starting_year=1958,
+        fitted_estimator = fitted_linear_margin_estimator(model_class, coordinates, dataset, starting_year=None,
                                                           fit_method=fit_method)
         # Load object from result from extremes
         return cls.from_estimator_extremes(fitted_estimator, ci_method, temporal_covariate)