[SCM][TREND TEST] add linear margin model with respect to the scale and the shape.

experiment/meteo_france_SCM_study/visualization/studies_visualization/main_studies_visualizer.py

 from experiment.trend_analysis.abstract_score import MannKendall, WeigthedScore, MeanScore, MedianScore
-from experiment.trend_analysis.univariate_trend_test.abstract_gev_trend_test import GevLocationTrendTest
+from experiment.trend_analysis.univariate_trend_test.abstract_gev_trend_test import GevLocationTrendTest, \
+    GevScaleTrendTest, GevShapeTrendTest
 from experiment.trend_analysis.univariate_trend_test.abstract_trend_test import MannKendallTrendTest
 from experiment.meteo_france_SCM_study.safran.safran import ExtendedSafranTotalPrecip
 from experiment.meteo_france_SCM_study.visualization.studies_visualization.studies import Studies
@@ -45,7 +46,7 @@ def altitude_trends_significant():
     # altitudes = ALL_ALTITUDES[3:5]
     # altitudes = ALL_ALTITUDES[2:4]
     for study_class in SCM_STUDIES[:1]:
-        trend_test_classes = [MannKendallTrendTest, GevLocationTrendTest][1:]
+        trend_test_classes = [MannKendallTrendTest, GevLocationTrendTest, GevScaleTrendTest, GevShapeTrendTest][2:]
         visualizers = [StudyVisualizer(study, temporal_non_stationarity=True, verbose=False)
                        for study in study_iterator_global(study_classes=[study_class], only_first_one=only_first_one,
                                                           altitudes=altitudes)]

experiment/trend_analysis/univariate_trend_test/abstract_gev_trend_test.py

@@ -5,7 +5,8 @@ from sklearn.preprocessing import normalize
 from experiment.trend_analysis.univariate_trend_test.abstract_trend_test import AbstractTrendTest
 from extreme_estimator.estimator.margin_estimator.abstract_margin_estimator import LinearMarginEstimator
 from extreme_estimator.extreme_models.margin_model.temporal_linear_margin_model import StationaryStationModel, \
-    NonStationaryLocationStationModel
+    NonStationaryLocationStationModel, NonStationaryScaleStationModel, NonStationaryShapeStationModel
+from extreme_estimator.margin_fits.gev.gev_params import GevParams
 from spatio_temporal_dataset.coordinates.abstract_coordinates import AbstractCoordinates
 from spatio_temporal_dataset.coordinates.temporal_coordinates.abstract_temporal_coordinates import \
     AbstractTemporalCoordinates
@@ -55,4 +56,27 @@ class GevLocationTrendTest(AbstractGevTrendTest):
 
     @property
     def test_sign(self) -> int:
-        return np.sign(self.non_stationary_estimator.margin_function_fitted.mu1_temporal_trend)
+        return np.sign(self.non_stationary_estimator.margin_function_fitted.get_coef(GevParams.LOC,
+                                                                                     AbstractCoordinates.COORDINATE_T))
+
+
+class GevScaleTrendTest(AbstractGevTrendTest):
+
+    def __init__(self, years_after_change_point, maxima_after_change_point):
+        super().__init__(years_after_change_point, maxima_after_change_point, NonStationaryScaleStationModel)
+
+    @property
+    def test_sign(self) -> int:
+        return np.sign(self.non_stationary_estimator.margin_function_fitted.get_coef(GevParams.SCALE,
+                                                                                     AbstractCoordinates.COORDINATE_T))
+
+
+class GevShapeTrendTest(AbstractGevTrendTest):
+
+    def __init__(self, years_after_change_point, maxima_after_change_point):
+        super().__init__(years_after_change_point, maxima_after_change_point, NonStationaryShapeStationModel)
+
+    @property
+    def test_sign(self) -> int:
+        return np.sign(self.non_stationary_estimator.margin_function_fitted.get_coef(GevParams.SHAPE,
+                                                                                     AbstractCoordinates.COORDINATE_T))

extreme_estimator/extreme_models/margin_model/margin_function/linear_margin_function.py

@@ -82,6 +82,11 @@ class LinearMarginFunction(ParametricMarginFunction):
 
     # Properties for the location parameter
 
+    def get_coef(self, gev_param_name, coef_name):
+        idx = 1 if coef_name in [AbstractCoordinates.COORDINATE_T] \
+            else AbstractCoordinates.COORDINATE_SPATIAL_NAMES.index(coef_name) + 2
+        return self.coef_dict[LinearCoef.coef_template_str(gev_param_name, coef_name).format(idx)]
+    
     @property
     def mu1_temporal_trend(self):
         return self.coef_dict[LinearCoef.coef_template_str(ExtremeParams.LOC, AbstractCoordinates.COORDINATE_T).format(1)]

extreme_estimator/extreme_models/margin_model/temporal_linear_margin_model.py

@@ -23,12 +23,21 @@ class TemporalLinearMarginModel(LinearMarginModel):
         # Gev Fit
         assert data.shape[1] == len(df_coordinates_temp.values)
         res = safe_run_r_estimator(function=r('gev.fit'), use_start=self.use_start_value,
-                                   xdat=ro.FloatVector(data[0]), y=df_coordinates_temp.values, mul=self.mul)
+                                   xdat=ro.FloatVector(data[0]), y=df_coordinates_temp.values, mul=self.mul,
+                                   sigl=self.sigl, shl=self.shl)
         return ResultFromIsmev(res, self.margin_function_start_fit.gev_param_name_to_dims)
 
     @property
     def mul(self):
-        return NotImplementedError
+        return get_null()
+
+    @property
+    def sigl(self):
+        return get_null()
+
+    @property
+    def shl(self):
+        return get_null()
 
 
 class StationaryStationModel(TemporalLinearMarginModel):
@@ -36,10 +45,6 @@ class StationaryStationModel(TemporalLinearMarginModel):
     def load_margin_functions(self, gev_param_name_to_dims=None):
         super().load_margin_functions({})
 
-    @property
-    def mul(self):
-        return get_null()
-
 
 class NonStationaryLocationStationModel(TemporalLinearMarginModel):
 
@@ -51,3 +56,21 @@ class NonStationaryLocationStationModel(TemporalLinearMarginModel):
         return 1
 
 
+class NonStationaryScaleStationModel(TemporalLinearMarginModel):
+
+    def load_margin_functions(self, gev_param_name_to_dims=None):
+        super().load_margin_functions({GevParams.SCALE: [self.coordinates.idx_temporal_coordinates]})
+
+    @property
+    def sigl(self):
+        return 1
+
+
+class NonStationaryShapeStationModel(TemporalLinearMarginModel):
+
+    def load_margin_functions(self, gev_param_name_to_dims=None):
+        super().load_margin_functions({GevParams.SHAPE: [self.coordinates.idx_temporal_coordinates]})
+
+    @property
+    def shl(self):
+        return 1

test/test_extreme_estimator/test_margin_fits/test_gev/test_gev_temporal.py

@@ -13,6 +13,7 @@ from spatio_temporal_dataset.coordinates.temporal_coordinates.abstract_temporal_
 from spatio_temporal_dataset.dataset.abstract_dataset import AbstractDataset
 from spatio_temporal_dataset.spatio_temporal_observations.abstract_spatio_temporal_observations import \
     AbstractSpatioTemporalObservations
+from test.test_utils import load_non_stationary_temporal_margin_models
 
 
 class TestGevTemporal(unittest.TestCase):
@@ -32,6 +33,8 @@ class TestGevTemporal(unittest.TestCase):
         df2 = pd.DataFrame(data=np.array(r['x_gev']), index=df.index)
         observations = AbstractSpatioTemporalObservations(df_maxima_gev=df2)
         self.dataset = AbstractDataset(observations=observations, coordinates=self.coordinates)
+        #
+        self.margin_models = load_non_stationary_temporal_margin_models(self.coordinates)
 
     def test_gev_temporal_margin_fit_stationary(self):
         # Create estimator
@@ -46,14 +49,14 @@ class TestGevTemporal(unittest.TestCase):
 
     def test_gev_temporal_margin_fit_nonstationary(self):
         # Create estimator
-        margin_model = NonStationaryLocationStationModel(self.coordinates)
-        estimator = LinearMarginEstimator(self.dataset, margin_model)
-        estimator.fit()
-        self.assertNotEqual(estimator.margin_function_fitted.mu1_temporal_trend, 0.0)
-        # Checks that parameters returned are indeed different
-        mle_params_estimated_year1 = estimator.margin_function_fitted.get_gev_params(np.array([1])).to_dict()
-        mle_params_estimated_year3 = estimator.margin_function_fitted.get_gev_params(np.array([3])).to_dict()
-        self.assertNotEqual(mle_params_estimated_year1, mle_params_estimated_year3)
+        for margin_model in self.margin_models:
+            # margin_model = NonStationaryLocationStationModel(self.coordinates)
+            estimator = LinearMarginEstimator(self.dataset, margin_model)
+            estimator.fit()
+            # Checks that parameters returned are indeed different
+            mle_params_estimated_year1 = estimator.margin_function_fitted.get_gev_params(np.array([1])).to_dict()
+            mle_params_estimated_year3 = estimator.margin_function_fitted.get_gev_params(np.array([3])).to_dict()
+            self.assertNotEqual(mle_params_estimated_year1, mle_params_estimated_year3)
 
     def test_gev_temporal_margin_fit_nonstationary_with_start_point(self):
         # Create estimator
@@ -80,7 +83,6 @@ class TestGevTemporal(unittest.TestCase):
         estimator2 = self.fit_non_stationary_estimator(starting_point=28)
         mu1_estimator1 = estimator1.margin_function_fitted.mu1_temporal_trend
         mu1_estimator2 = estimator2.margin_function_fitted.mu1_temporal_trend
-        print(mu1_estimator1, mu1_estimator2)
         self.assertNotEqual(mu1_estimator1, mu1_estimator2)
 
 

test/test_utils.py

@@ -9,6 +9,8 @@ from extreme_estimator.estimator.full_estimator.abstract_full_estimator import S
 from extreme_estimator.estimator.max_stable_estimator.abstract_max_stable_estimator import MaxStableEstimator
 from extreme_estimator.extreme_models.margin_model.linear_margin_model import LinearAllParametersAllDimsMarginModel, \
     ConstantMarginModel
+from extreme_estimator.extreme_models.margin_model.temporal_linear_margin_model import \
+    NonStationaryLocationStationModel, NonStationaryScaleStationModel, NonStationaryShapeStationModel
 from extreme_estimator.extreme_models.max_stable_model.abstract_max_stable_model import \
     AbstractMaxStableModelWithCovarianceFunction, CovarianceFunction
 from extreme_estimator.extreme_models.max_stable_model.max_stable_models import Smith, BrownResnick, Schlather, \
@@ -37,10 +39,15 @@ TEST_3D_SPATIAL_COORDINATES = [AlpsStation3DCoordinatesWithAnisotropy]
 TEST_TEMPORAL_COORDINATES = [ConsecutiveTemporalCoordinates]
 TEST_SPATIO_TEMPORAL_COORDINATES = [UniformSpatioTemporalCoordinates, LinSpaceSpatial2DSpatioTemporalCoordinates]
 TEST_MARGIN_TYPES = [ConstantMarginModel, LinearAllParametersAllDimsMarginModel][:]
+TEST_NON_STATIONARY_TEMPORAL_MARGIN_TYPES = [NonStationaryLocationStationModel, NonStationaryScaleStationModel, NonStationaryShapeStationModel]
 TEST_MAX_STABLE_ESTIMATOR = [MaxStableEstimator]
 TEST_FULL_ESTIMATORS = [SmoothMarginalsThenUnitaryMsp, FullEstimatorInASingleStepWithSmoothMargin][:]
 
 
+def load_non_stationary_temporal_margin_models(coordinates):
+    return [margin_class(coordinates=coordinates) for margin_class in TEST_NON_STATIONARY_TEMPORAL_MARGIN_TYPES]
+
+
 def load_test_full_estimators(dataset, margin_model, max_stable_model):
     return [full_estimator(dataset=dataset, margin_model=margin_model, max_stable_model=max_stable_model) for
             full_estimator in TEST_FULL_ESTIMATORS]
@@ -89,7 +96,8 @@ def load_test_temporal_coordinates(nb_steps, train_split_ratio=None):
             TEST_TEMPORAL_COORDINATES]
 
 
-def load_test_spatiotemporal_coordinates(nb_points, nb_steps, train_split_ratio=None, transformation_class: type = None):
+def load_test_spatiotemporal_coordinates(nb_points, nb_steps, train_split_ratio=None,
+                                         transformation_class: type = None):
     return [coordinate_class.from_nb_points_and_nb_steps(nb_points=nb_points, nb_steps=nb_steps,
                                                          train_split_ratio=train_split_ratio,
                                                          transformation_class=transformation_class)