[STRENGTH EVOLUTION] fix strength computation for the scale test.

experiment/paper1_steps/hard extreme evolution - annual maxima/main1_good_stationary_gev_fit.py

@@ -6,7 +6,7 @@ from experiment.meteo_france_data.scm_models_data.visualization.study_visualizat
 
 
 def maxima_analysis():
-    save_to_file = True
+    save_to_file = False
     only_first_one = False
     durand_altitude = [900, 1500, 2100, 2700]
     altitudes = durand_altitude

experiment/paper1_steps/hard extreme evolution - annual maxima/main3_non_stationary_strength_evolution.py

@@ -2,6 +2,8 @@ import time
 
 from experiment.meteo_france_data.scm_models_data.visualization.hypercube_visualization.altitude_year_hypercube_visualizer import \
     Altitude_Hypercube_Year_Visualizer
+from experiment.trend_analysis.univariate_test.abstract_gev_change_point_test import GevLocationChangePointTest, \
+    GevScaleChangePointTest
 
 """
 Visualize the 0.99 quantile initial value and its evolution
@@ -12,25 +14,24 @@ from experiment.paper1_steps.utils import get_full_altitude_visualizer, FULL_ALT
 def main_fast_spatial_risk_evolution():
     for altitude in FULL_ALTITUDES[-1:]:
         vizualiser = get_full_altitude_visualizer(Altitude_Hypercube_Year_Visualizer, altitude=altitude,
-                                                  exact_starting_year=1958, reduce_strength_array=True)
+                                                  exact_starting_year=1958, reduce_strength_array=True,
+                                                  trend_test_class=GevScaleChangePointTest)
         vizualiser.save_to_file = False
         vizualiser.visualize_massif_trend_test_one_altitude()
 
 
-def main_full_spatial_repartition():
+def main_full_spatial_risk_evolution():
     for altitude in FULL_ALTITUDES[:]:
-        # Compute for the most likely starting year
-        # vizualiser = get_full_altitude_visualizer(Altitude_Hypercube_Year_Visualizer, altitude=altitude)
-        # vizualiser.visualize_massif_trend_test_one_altitude()
-        # Compute the trend for a linear trend
-        vizualiser = get_full_altitude_visualizer(Altitude_Hypercube_Year_Visualizer, altitude=altitude,
-                                                  exact_starting_year=1958)
-        vizualiser.visualize_massif_trend_test_one_altitude()
+        for trend_test_class in [GevLocationChangePointTest, GevScaleChangePointTest][1:]:
+            vizualiser = get_full_altitude_visualizer(Altitude_Hypercube_Year_Visualizer, altitude=altitude,
+                                                      exact_starting_year=1958, reduce_strength_array=True,
+                                                      trend_test_class=trend_test_class)
+            vizualiser.visualize_massif_trend_test_one_altitude()
 
 
 def main_run():
-    # main_full_spatial_repartition()
-    main_fast_spatial_risk_evolution()
+    main_full_spatial_risk_evolution()
+    # main_fast_spatial_risk_evolution()
 
 
 if __name__ == '__main__':

experiment/paper1_steps/utils.py

@@ -14,13 +14,13 @@ FULL_ALTITUDES = [900, 1200, 1500, 1800, 2100, 2400, 2700, 3000]
 
 
 def get_full_altitude_visualizer(altitude_hypercube_class, exact_starting_year=None, altitude=900,
-                                           reduce_strength_array=False):
-    altitudes, first_starting_year, last_starting_year, nb_data_reduced_for_speed, only_first_one, save_to_file, trend_test_class = get_full_parameters(
+                                 reduce_strength_array=False,
+                                 trend_test_class = GevLocationChangePointTest):
+    altitudes, first_starting_year, last_starting_year, nb_data_reduced_for_speed, only_first_one, save_to_file, _ = get_full_parameters(
         altitude=altitude)
     if exact_starting_year is not None:
         first_starting_year, last_starting_year = None, None
     study_classes = [CrocusRecentSwe]
-    trend_test_class = GevLocationChangePointTest
     visualizer = load_altitude_visualizer(altitude_hypercube_class, altitudes, last_starting_year,
                                           nb_data_reduced_for_speed, only_first_one, save_to_file, study_classes,
                                           trend_test_class, first_starting_year=first_starting_year,

experiment/trend_analysis/univariate_test/abstract_gev_change_point_test.py

@@ -169,6 +169,11 @@ class GevScaleChangePointTest(AbstractGevChangePointTest):
         super().__init__(years, maxima, starting_year,
                          NonStationaryScaleStationModel, GevParams.SCALE)
 
+    def _slope_strength(self):
+        return self.non_stationary_constant_gev_params.quantile_strength_evolution_ratio(
+            p=self.quantile_for_strength,
+            sigma1=self.non_stationary_linear_coef)
+
 
 class GevShapeChangePointTest(AbstractGevChangePointTest):
 

extreme_estimator/margin_fits/gev/gev_params.py

@@ -58,7 +58,8 @@ class GevParams(ExtremeParams):
         initial_quantile = self.quantile(p)
         quantity_increased = mu1
         if sigma1 != 0:
-            quantity_increased += (sigma1 / self.shape) * (1 - (- np.float_power(np.log(p), -self.shape)))
+            power = np.float_power(- np.log(p), -self.shape)
+            quantity_increased -= (sigma1 / self.shape) * (1 - power)
         return quantity_increased / initial_quantile
 
     # Compute some indicators (such as the mean and the variance)