[contrasting] change name of the model. add quadratic models. improve plots

extreme_data/meteo_france_data/scm_models_data/visualization/plot_utils.py

@@ -80,4 +80,4 @@ def load_plot(cmap, graduation, label, massif_name_to_value, altitude, fit_metho
     ax.set_xticks([])
     if add_x_label:
         ax.set_xlabel('Altitude = {}m'.format(altitude), fontsize=15)
-        ax.set_title('Fit method is {}'.format(fit_method))
+        # ax.set_title('Fit method is {}'.format(fit_method))

extreme_fit/model/margin_model/polynomial_margin_model/gev_altitudinal_models.py

@@ -34,19 +34,37 @@ class AbstractAltitudinalModel(AbstractSpatioTemporalPolynomialModel):
 
     @property
     def name_str(self):
-        name = ''
-        for coordinate_name, idx in zip(['s', 't'], [self.coordinates.idx_x_coordinates, self.coordinates.idx_temporal_coordinates]):
-            name += coordinate_name
-            for param_name in GevParams.PARAM_NAMES:
-                name += self.dim_to_str_number(param_name, idx)
-        if isinstance(self, AbstractAddCrossTerm):
-            name += 'sxt'
-            for c in [AbstractAddCrossTermForLocation, AbstractAddCrossTermForScale, AbstractAddCrossTermForShape]:
-                if isinstance(self, c):
-                    name += '1'
-                else:
-                    name += '0'
-        return name
+        s = ''
+        if self.dim_to_str_number(GevParams.SHAPE, self.coordinates.idx_x_coordinates) == '1':
+            s += '\\zeta_s'
+        if self.dim_to_str_number(GevParams.LOC, self.coordinates.idx_temporal_coordinates) in ['1', '2']:
+            s += '\\mu_t'
+            if self.dim_to_str_number(GevParams.LOC, self.coordinates.idx_temporal_coordinates) == '2':
+                s += '^2'
+        if self.dim_to_str_number(GevParams.SCALE, self.coordinates.idx_temporal_coordinates) in ['1', '2']:
+            s += '\\sigma_t'
+            if self.dim_to_str_number(GevParams.LOC, self.coordinates.idx_temporal_coordinates) == '2':
+                s += '^2'
+        if len(s) == 0:
+            s = '0'
+        return '$ \\boldsymbol{ \\mathcal{M}_{%s} }$' % s
+
+
+    # @property
+    # def name_str(self):
+    #     name = ''
+    #     for coordinate_name, idx in zip(['s', 't'], [self.coordinates.idx_x_coordinates, self.coordinates.idx_temporal_coordinates]):
+    #         name += coordinate_name
+    #         for param_name in GevParams.PARAM_NAMES:
+    #             name += self.dim_to_str_number(param_name, idx)
+    #     if isinstance(self, AbstractAddCrossTerm):
+    #         name += 'sxt'
+    #         for c in [AbstractAddCrossTermForLocation, AbstractAddCrossTermForScale, AbstractAddCrossTermForShape]:
+    #             if isinstance(self, c):
+    #                 name += '1'
+    #             else:
+    #                 name += '0'
+    #     return name
 
 
 class StationaryAltitudinal(AbstractAltitudinalModel):
@@ -134,6 +152,7 @@ class NonStationaryAltitudinalLocationQuadraticScaleLinear(AbstractAltitudinalMo
 class AbstractAddCrossTerm(AbstractAltitudinalModel):
     pass
 
+
 class AbstractAddCrossTermForLocation(AbstractAddCrossTerm):
 
     # @property
@@ -180,7 +199,8 @@ class AbstractAddCrossTermForShape(AbstractAddCrossTerm):
         if GevParams.SHAPE in d:
             insert_index = 1 if self.coordinates.idx_x_coordinates == d[GevParams.SHAPE][0][0] else 0
             d[GevParams.SHAPE].insert(insert_index,
-                                      ((self.coordinates.idx_x_coordinates, self.coordinates.idx_temporal_coordinates), 1))
+                                      ((self.coordinates.idx_x_coordinates, self.coordinates.idx_temporal_coordinates),
+                                       1))
         else:
             d[GevParams.SHAPE] = [((self.coordinates.idx_x_coordinates, self.coordinates.idx_temporal_coordinates), 1)]
         return d

extreme_fit/model/margin_model/polynomial_margin_model/models_based_on_pariwise_analysis/gev_quadratic_location_and_scale_with_constant_shape.py

+from extreme_fit.distribution.gev.gev_params import GevParams
+from extreme_fit.model.margin_model.polynomial_margin_model.gev_altitudinal_models import AbstractAltitudinalModel
+
+
+class AltitudinalShapeConstantTimeLocationQuadratic(AbstractAltitudinalModel):
+
+    @property
+    def param_name_to_list_dim_and_degree(self):
+        return {
+            GevParams.LOC: [(self.coordinates.idx_x_coordinates, 1), (self.coordinates.idx_temporal_coordinates, 2)],
+            GevParams.SCALE: [(self.coordinates.idx_x_coordinates, 1)]
+        }
+
+
+class AltitudinalShapeConstantTimeLocationQuadraticScaleLinear(AbstractAltitudinalModel):
+
+    @property
+    def param_name_to_list_dim_and_degree(self):
+        return {
+            GevParams.LOC: [(self.coordinates.idx_x_coordinates, 1), (self.coordinates.idx_temporal_coordinates, 2)],
+            GevParams.SCALE: [(self.coordinates.idx_x_coordinates, 1), (self.coordinates.idx_temporal_coordinates, 1)]
+        }
+
+
+class AltitudinalShapeConstantTimeLocationQuadraticScaleQuadratic(AbstractAltitudinalModel):
+
+    @property
+    def param_name_to_list_dim_and_degree(self):
+        return {
+            GevParams.LOC: [(self.coordinates.idx_x_coordinates, 1), (self.coordinates.idx_temporal_coordinates, 2)],
+            GevParams.SCALE: [(self.coordinates.idx_x_coordinates, 1), (self.coordinates.idx_temporal_coordinates, 2)]
+        }
+
+
+class AltitudinalShapeConstantTimeLocationLinearScaleQuadratic(AbstractAltitudinalModel):
+
+    @property
+    def param_name_to_list_dim_and_degree(self):
+        return {
+            GevParams.LOC: [(self.coordinates.idx_x_coordinates, 1), (self.coordinates.idx_temporal_coordinates, 1)],
+            GevParams.SCALE: [(self.coordinates.idx_x_coordinates, 1), (self.coordinates.idx_temporal_coordinates, 2)]
+        }
+
+
+class AltitudinalShapeConstantTimeScaleQuadratic(AbstractAltitudinalModel):
+
+    @property
+    def param_name_to_list_dim_and_degree(self):
+        return {
+            GevParams.LOC: [(self.coordinates.idx_x_coordinates, 1)],
+            GevParams.SCALE: [(self.coordinates.idx_x_coordinates, 1), (self.coordinates.idx_temporal_coordinates, 2)]
+        }

extreme_fit/model/margin_model/polynomial_margin_model/models_based_on_pariwise_analysis/gev_quadratic_location_and_scale_with_linear_shape.py

+from extreme_fit.distribution.gev.gev_params import GevParams
+from extreme_fit.model.margin_model.polynomial_margin_model.gev_altitudinal_models import AbstractAltitudinalModel
+
+
+class AltitudinalShapeLinearTimeLocationQuadratic(AbstractAltitudinalModel):
+
+    @property
+    def param_name_to_list_dim_and_degree(self):
+        return {
+            GevParams.LOC: [(self.coordinates.idx_x_coordinates, 1), (self.coordinates.idx_temporal_coordinates, 2)],
+            GevParams.SCALE: [(self.coordinates.idx_x_coordinates, 1)],
+            GevParams.SHAPE: [(self.coordinates.idx_x_coordinates, 1)]
+        }
+
+
+class AltitudinalShapeLinearTimeLocationQuadraticScaleLinear(AbstractAltitudinalModel):
+
+    @property
+    def param_name_to_list_dim_and_degree(self):
+        return {
+            GevParams.LOC: [(self.coordinates.idx_x_coordinates, 1), (self.coordinates.idx_temporal_coordinates, 2)],
+            GevParams.SCALE: [(self.coordinates.idx_x_coordinates, 1), (self.coordinates.idx_temporal_coordinates, 1)],
+            GevParams.SHAPE: [(self.coordinates.idx_x_coordinates, 1)]
+        }
+
+
+class AltitudinalShapeLinearTimeLocationQuadraticScaleQuadratic(AbstractAltitudinalModel):
+
+    @property
+    def param_name_to_list_dim_and_degree(self):
+        return {
+            GevParams.LOC: [(self.coordinates.idx_x_coordinates, 1), (self.coordinates.idx_temporal_coordinates, 2)],
+            GevParams.SCALE: [(self.coordinates.idx_x_coordinates, 1), (self.coordinates.idx_temporal_coordinates, 2)],
+            GevParams.SHAPE: [(self.coordinates.idx_x_coordinates, 1)]
+        }
+
+
+class AltitudinalShapeLinearTimeLocationLinearScaleQuadratic(AbstractAltitudinalModel):
+
+    @property
+    def param_name_to_list_dim_and_degree(self):
+        return {
+            GevParams.LOC: [(self.coordinates.idx_x_coordinates, 1), (self.coordinates.idx_temporal_coordinates, 1)],
+            GevParams.SCALE: [(self.coordinates.idx_x_coordinates, 1), (self.coordinates.idx_temporal_coordinates, 2)],
+            GevParams.SHAPE: [(self.coordinates.idx_x_coordinates, 1)]
+        }
+
+
+class AltitudinalShapeLinearTimeScaleQuadratic(AbstractAltitudinalModel):
+
+    @property
+    def param_name_to_list_dim_and_degree(self):
+        return {
+            GevParams.LOC: [(self.coordinates.idx_x_coordinates, 1)],
+            GevParams.SCALE: [(self.coordinates.idx_x_coordinates, 1), (self.coordinates.idx_temporal_coordinates, 2)],
+            GevParams.SHAPE: [(self.coordinates.idx_x_coordinates, 1)]
+        }

extreme_fit/model/margin_model/polynomial_margin_model/utils.py

@@ -29,6 +29,14 @@ from extreme_fit.model.margin_model.polynomial_margin_model.gumbel_altitudinal_m
     NonStationaryGumbelAltitudinalLocationLinearScaleLinearCrossTermForLocation, \
     NonStationaryGumbelAltitudinalLocationQuadraticScaleLinearCrossTermForLocation, \
     NonStationaryGumbelAltitudinalScaleLinearCrossTermForLocation, NonStationaryGumbelCrossTermForLocation
+from extreme_fit.model.margin_model.polynomial_margin_model.models_based_on_pariwise_analysis.gev_quadratic_location_and_scale_with_constant_shape import \
+    AltitudinalShapeConstantTimeLocationQuadratic, AltitudinalShapeConstantTimeLocationQuadraticScaleLinear, \
+    AltitudinalShapeConstantTimeLocationQuadraticScaleQuadratic, \
+    AltitudinalShapeConstantTimeLocationLinearScaleQuadratic, AltitudinalShapeConstantTimeScaleQuadratic
+from extreme_fit.model.margin_model.polynomial_margin_model.models_based_on_pariwise_analysis.gev_quadratic_location_and_scale_with_linear_shape import \
+    AltitudinalShapeLinearTimeLocationQuadratic, AltitudinalShapeLinearTimeLocationQuadraticScaleLinear, \
+    AltitudinalShapeLinearTimeLocationQuadraticScaleQuadratic, AltitudinalShapeLinearTimeLocationLinearScaleQuadratic, \
+    AltitudinalShapeLinearTimeScaleQuadratic
 from extreme_fit.model.margin_model.polynomial_margin_model.models_based_on_pariwise_analysis.gev_with_constant_shape_wrt_altitude import \
     AltitudinalShapeConstantTimeLocShapeLinear, \
     AltitudinalShapeConstantTimeLocScaleLinear, AltitudinalShapeConstantTimeScaleShapeLinear, \
@@ -112,6 +120,20 @@ ALTITUDINAL_GEV_MODELS_BASED_ON_POINTWISE_ANALYSIS = [
     # AltitudinalShapeLinearTimeScaleShapeLinear,
     # AltitudinalShapeLinearTimeLocScaleShapeLinear,
 
+    # With a quadratic term
+
+    # AltitudinalShapeConstantTimeLocationQuadratic,
+    # AltitudinalShapeConstantTimeLocationQuadraticScaleLinear,
+    # AltitudinalShapeConstantTimeLocationQuadraticScaleQuadratic,
+    # AltitudinalShapeConstantTimeLocationLinearScaleQuadratic,
+    # AltitudinalShapeConstantTimeScaleQuadratic,
+    #
+    # AltitudinalShapeLinearTimeLocationQuadratic,
+    # AltitudinalShapeLinearTimeLocationQuadraticScaleLinear,
+    # AltitudinalShapeLinearTimeLocationQuadraticScaleQuadratic,
+    # AltitudinalShapeLinearTimeLocationLinearScaleQuadratic,
+    # AltitudinalShapeLinearTimeScaleQuadratic,
+
 ]
 
 ALTITUDINAL_GEV_MODELS_LOCATION_ONLY_SCALE_ALTITUDES = [

projects/altitude_spatial_model/altitudes_fit/main_altitudes_studies.py

+
+import time
 from typing import List
 
 import matplotlib as mpl
-import matplotlib.pyplot as plt
-import numpy as np
-
-from extreme_data.meteo_france_data.adamont_data.abstract_simulation_study import SimulationStudy
-from extreme_data.meteo_france_data.adamont_data.snowfall_simulation import SafranSnowfallSimulationRCP85
-from extreme_fit.model.margin_model.polynomial_margin_model.utils import ALTITUDINAL_GEV_MODELS, \
-    ALTITUDINAL_GEV_MODELS_LOCATION_QUADRATIC_MINIMUM, ALTITUDINAL_GEV_MODELS_LOCATION_ONLY_SCALE_ALTITUDES, \
-    ALTITUDINAL_GEV_MODELS_LOCATION, ALTITUDINAL_GEV_MODELS_BASED_ON_POINTWISE_ANALYSIS
-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_total_aic, \
-    plot_individual_aic
-from spatio_temporal_dataset.coordinates.temporal_coordinates.temperature_covariate import MeanAlpsTemperatureCovariate
 
 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
+from extreme_fit.model.margin_model.polynomial_margin_model.utils import \
+    ALTITUDINAL_GEV_MODELS_BASED_ON_POINTWISE_ANALYSIS
+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, SafranPrecipitation1Day, SafranPrecipitation3Days, \
@@ -47,12 +38,13 @@ def main():
     study_classes = [SafranSnowfall1Day, SafranSnowfall3Days, SafranPrecipitation1Day
                         , SafranPrecipitation3Days][:1]
     altitudes = [1800, 2100, 2400]
-    study_classes = [SafranSnowfall1Day, SafranSnowfall3Days][:1]
+    study_classes = [SafranSnowfall1Day, SafranSnowfall3Days, SafranSnowfall5Days, SafranSnowfall7Days][:1]
+    # study_classes = [SafranPrecipitation1Day][:1]
 
     # Common parameters
     # altitudes = [600, 900, 1200, 1500, 1800, 2100, 2400, 2700, 3000, 3300, 3600]
     massif_names = None
-    massif_names = ['Mercantour', 'Vercors', 'Ubaye']
+    # massif_names = ['Mercantour', 'Vercors', 'Ubaye']
     seasons = [Season.annual, Season.winter, Season.spring, Season.automn][:1]
 
     main_loop(altitudes_for_groups, massif_names, seasons, study_classes)
@@ -72,6 +64,8 @@ def main_loop(altitudes_list, massif_names, seasons, study_classes):
             visualizer_list = load_visualizer_list(season, study_class, altitudes_list, massif_names)
             for visualizer in visualizer_list:
                 plots(massif_names, season, visualizer)
+            del visualizer_list
+            time.sleep(2)
 
 
 def load_visualizer_list(season, study_class, altitudes_list, massif_names):
@@ -79,7 +73,6 @@ def load_visualizer_list(season, study_class, altitudes_list, massif_names):
     visualizer_list = []
     # Load all studies
     for altitudes in altitudes_list:
-        print('here', altitudes)
         studies = AltitudesStudies(study_class, altitudes, season=season)
         visualizer = AltitudesStudiesVisualizerForNonStationaryModels(studies=studies,
                                                                       model_classes=model_classes,

projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/altitudes_studies_visualizer_for_non_stationary_models.py

@@ -61,8 +61,8 @@ class AltitudesStudiesVisualizerForNonStationaryModels(StudyVisualizer):
         self._method_name_and_order_to_max_abs = {}
         self._max_abs_for_shape = None
 
-    moment_names = ['moment', 'changes_in_the_moment', 'relative_changes_in_the_moment'][2:]
-    orders = [1, 2, None][2:]
+    moment_names = ['moment', 'changes_for_moment', 'relative_changes_for_moment'][:]
+    orders = [1, 2, None][:]
 
     @property
     def massif_name_to_one_fold_fit(self) -> Dict[str, OneFoldFit]:
@@ -106,10 +106,13 @@ class AltitudesStudiesVisualizerForNonStationaryModels(StudyVisualizer):
         # Plot settings
         moment = ' '.join(method_name.split('_'))
         moment = moment.replace('moment', '{} in 2019'.format(OneFoldFit.get_moment_str(order=order)))
-        plot_name = '{}{} annual maxima of {}'.format(OneFoldFit.folder_for_plots, moment,
+        plot_name = '{}{} for annual maxima of {}'.format(OneFoldFit.folder_for_plots, moment,
                                                       SCM_STUDY_CLASS_TO_ABBREVIATION[
                                                           self.studies.study_class])
-        ylabel = '{} ({})'.format(plot_name, self.study.variable_unit)
+        if 'change' in method_name:
+            plot_name += '\nbetween {} and {}'.format(2019-50, 2019)
+        parenthesis = self.study.variable_unit if 'relative' not in method_name else '\%'
+        ylabel = '{} ({})'.format(plot_name, parenthesis)
 
         # Plot the map
         a_change_is_displayed = self.moment_names.index(method_name) > 0
@@ -117,9 +120,10 @@ class AltitudesStudiesVisualizerForNonStationaryModels(StudyVisualizer):
                       label=ylabel, massif_name_to_value=massif_name_to_value,
                       plot_name=plot_name, add_x_label=True,
                       negative_and_positive_values=a_change_is_displayed,
-                      massif_name_to_text=None, altitude=self.altitude_group.reference_altitude,
+                      altitude=self.altitude_group.reference_altitude,
                       add_colorbar=self.add_colorbar,
-                      max_abs_change=self.method_name_and_order_to_max_abs(method_name, order)
+                      max_abs_change=self.method_name_and_order_to_max_abs(method_name, order),
+                      massif_name_to_text=self.massif_name_to_name,
                       )
 
     @property

projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/one_fold_fit.py

@@ -24,11 +24,13 @@ from spatio_temporal_dataset.dataset.abstract_dataset import AbstractDataset
 class OneFoldFit(object):
     SIGNIFICANCE_LEVEL = 0.05
     best_estimator_minimizes_total_aic = False
+    return_period = 50
 
     def __init__(self, massif_name: str, dataset: AbstractDataset, models_classes,
                  fit_method=MarginFitMethod.extremes_fevd_mle, temporal_covariate_for_fit=None,
                  altitude_class=DefaultAltitudeGroup,
-                 only_models_that_pass_anderson_test=True):
+                 only_models_that_pass_anderson_test=True,
+                 ):
         self.only_models_that_pass_anderson_test = only_models_that_pass_anderson_test
         self.altitude_group = altitude_class()
         self.massif_name = massif_name
@@ -57,11 +59,11 @@ class OneFoldFit(object):
     @classmethod
     def get_moment_str(cls, order):
         if order == 1:
-            return 'Mean'
+            return 'mean'
         elif order == 2:
-            return 'Std'
+            return 'std'
         elif order is None:
-            return 'Return level'
+            return '{}-year return levels'.format(cls.return_period)
 
     def get_moment(self, altitude, year, order=1):
         gev_params = self.get_gev_params(altitude, year)
@@ -70,7 +72,7 @@ class OneFoldFit(object):
         elif order == 2:
             return gev_params.std
         elif order is None:
-            return gev_params.return_level(return_period=50)
+            return gev_params.return_level(return_period=self.return_period)
         else:
             raise NotImplementedError
 
@@ -82,7 +84,7 @@ class OneFoldFit(object):
     def moment(self, altitudes, year=2019, order=1):
         return [self.get_moment(altitude, year, order) for altitude in altitudes]
 
-    def changes_in_the_moment(self, altitudes, year=2019, nb_years=50, order=1):
+    def changes_for_moment(self, altitudes, year=2019, nb_years=50, order=1):
         changes = []
         for altitude in altitudes:
             mean_after = self.get_moment(altitude, year, order)
@@ -91,7 +93,7 @@ class OneFoldFit(object):
             changes.append(change)
         return changes
 
-    def relative_changes_in_the_moment(self, altitudes, year=2019, nb_years=50, order=1):
+    def relative_changes_for_moment(self, altitudes, year=2019, nb_years=50, order=1):
         relative_changes = []
         for altitude in altitudes:
             mean_after = self.get_moment(altitude, year, order)