[contrasting] add time derivative for the linear models. add one_axis_param_function.py

extreme_data/meteo_france_data/scm_models_data/visualization/plot_utils.py

@@ -17,7 +17,8 @@ def plot_against_altitude(altitudes, ax, massif_id, massif_name, values):
     colors = list(mcolors.TABLEAU_COLORS)
     colors[-3:-1] = []  # remove gray and olive
     color = colors[massif_id % 8]
-    ax.plot(altitudes, values, color=color, linewidth=2, label=massif_name, linestyle=linestyle)
+    massif_name_str = ' '.join(massif_name.split('_'))
+    ax.plot(altitudes, values, color=color, linewidth=2, label=massif_name_str, linestyle=linestyle)
 
 
 def load_plot(cmap, graduation, label, massif_name_to_value, altitude, fit_method, add_x_label=True):

extreme_fit/function/param_function/one_axis_param_function.py

+import numpy as np
+
+
+class AbstractOneAxisParamFunction(object):
+
+    def __init__(self, dim: int, coordinates: np.ndarray):
+        self.dim = dim
+        self.t_min = coordinates[:, dim].min()
+        self.t_max = coordinates[:, dim].max()
+
+    def get_param_value(self, coordinate: np.ndarray, assert_range=True) -> float:
+        t = coordinate[self.dim]
+        if assert_range:
+            assert self.t_min <= t <= self.t_max, '{} is out of bounds ({}, {})'.format(t, self.t_min, self.t_max)
+        return self._get_param_value(t)
+
+    def _get_param_value(self, t) -> float:
+        raise NotImplementedError
+
+
+class LinearOneAxisParamFunction(AbstractOneAxisParamFunction):
+
+    def __init__(self, dim: int, coordinates: np.ndarray, coef: float = 0.01):
+        super().__init__(dim, coordinates)
+        self.coef = coef
+
+    def _get_param_value(self, t) -> float:
+        return t * self.coef
+
+
+class QuadraticOneAxisParamFunction(AbstractOneAxisParamFunction):
+
+    def __init__(self, dim: int, coordinates: np.ndarray, coef1: float = 0.01, coef2: float = 0.01):
+        super().__init__(dim, coordinates)
+        self.coef1 = coef1
+        self.coef2 = coef2
+
+    def _get_param_value(self, t) -> float:
+        return np.power(t, 2) * self.coef2 + t * self.coef1

extreme_fit/function/param_function/param_function.py

 from typing import List
 import numpy as np
 from extreme_fit.function.param_function.linear_coef import LinearCoef
+from extreme_fit.function.param_function.one_axis_param_function import LinearOneAxisParamFunction
 from extreme_fit.function.param_function.spline_coef import SplineCoef
 
 
@@ -22,20 +23,8 @@ class ConstantParamFunction(AbstractParamFunction):
     def get_param_value(self, coordinate: np.ndarray) -> float:
         return self.constant
 
-
-class LinearOneAxisParamFunction(AbstractParamFunction):
-
-    def __init__(self, dim: int, coordinates: np.ndarray, coef: float = 0.01):
-        self.dim = dim
-        self.t_min = coordinates[:, dim].min()
-        self.t_max = coordinates[:, dim].max()
-        self.coef = coef
-
-    def get_param_value(self, coordinate: np.ndarray) -> float:
-        t = coordinate[self.dim]
-        if self.OUT_OF_BOUNDS_ASSERT:
-            assert self.t_min <= t <= self.t_max, '{} is out of bounds ({}, {})'.format(t, self.t_min, self.t_max)
-        return t * self.coef
+    def get_first_derivative_param_value(self, coordinate: np.ndarray, dim: int) -> float:
+        return 0
 
 
 class LinearParamFunction(AbstractParamFunction):
@@ -48,14 +37,20 @@ class LinearParamFunction(AbstractParamFunction):
             param_function = LinearOneAxisParamFunction(dim=dim, coordinates=coordinates,
                                                         coef=self.linear_coef.get_coef(idx=dim))
             self.linear_one_axis_param_functions.append(param_function)
+        self.dim_to_linear_one_axis_param_function = dict(zip(dims, self.linear_one_axis_param_functions))
 
     def get_param_value(self, coordinate: np.ndarray) -> float:
         # Add the intercept and the value with respect to each axis
         gev_param_value = self.linear_coef.intercept
         for linear_one_axis_param_function in self.linear_one_axis_param_functions:
-            gev_param_value += linear_one_axis_param_function.get_param_value(coordinate)
+            gev_param_value += linear_one_axis_param_function.get_param_value(coordinate, self.OUT_OF_BOUNDS_ASSERT)
         return gev_param_value
 
+    def get_first_derivative_param_value(self, coordinate: np.ndarray, dim: int) -> float:
+        return self.dim_to_linear_one_axis_param_function[dim].coef
+
+
+
 
 class SplineParamFunction(AbstractParamFunction):
 

projects/contrasting_trends_in_snow_loads/snowfall_versus_time_and_altitude/main_snowfall_article.py

@@ -4,7 +4,8 @@ from multiprocessing.pool import Pool
 import matplotlib as mpl
 
 from extreme_data.meteo_france_data.scm_models_data.safran.safran import SafranSnowfall1Day
-from projects.contrasting_trends_in_snow_loads.snowfall_versus_time_and_altitude.snowfall_plot import plot_snowfall_mean
+from projects.contrasting_trends_in_snow_loads.snowfall_versus_time_and_altitude.snowfall_plot import \
+    plot_snowfall_mean, plot_snowfall_time_derivative_mean
 from projects.contrasting_trends_in_snow_loads.snowfall_versus_time_and_altitude.study_visualizer_for_mean_values import \
     StudyVisualizerForMeanValues
 from projects.contrasting_trends_in_snow_loads.snowfall_versus_time_and_altitude.validation_plot import validation_plot
@@ -71,8 +72,9 @@ def intermediate_result(altitudes, massif_names=None,
             _ = compute_minimized_aic(visualizer)
 
     # Plots
-    validation_plot(altitude_to_visualizer)
-    plot_snowfall_mean(altitude_to_visualizer)
+    # validation_plot(altitude_to_visualizer)
+    # plot_snowfall_mean(altitude_to_visualizer)
+    plot_snowfall_time_derivative_mean(altitude_to_visualizer)
 
 
 def major_result():

projects/contrasting_trends_in_snow_loads/snowfall_versus_time_and_altitude/snowfall_plot.py

@@ -20,6 +20,30 @@ def fit_linear_regression(x, y):
     return a, b, r2_score
 
 
+def plot_snowfall_time_derivative_mean(altitude_to_visualizer: Dict[int, StudyVisualizerForMeanValues]):
+    visualizer = list(altitude_to_visualizer.values())[0]
+    # Plot the curve for the evolution of the mean
+    massif_name_to_a, massif_name_to_b, massif_name_to_r2_score = plot_mean(altitude_to_visualizer, derivative=True)
+    # Augmentation every km
+    massif_name_to_augmentation_every_km = {m: a * 1000 for m, a in massif_name_to_a.items()}
+    visualizer.plot_abstract_fast(massif_name_to_augmentation_every_km,
+                                  label='Augmentation of time derivative of mean annual maxima for every km of elevation ()',
+                                  add_x_label=False)
+    # Value at 2000 m
+    massif_name_to_mean_at_2000 = {m: a * 2000 + massif_name_to_b[m] for m, a in massif_name_to_a.items()}
+    visualizer.plot_abstract_fast(massif_name_to_mean_at_2000, label='Time derivative of mean annual maxima of at 2000 m ()',
+                                  add_x_label=False)
+    # Altitude for the change of dynamic
+    massif_name_to_altitude_change_dynamic = {m: - massif_name_to_b[m] / a for m, a in massif_name_to_a.items()}
+    # Keep only those that are in a reasonable range
+    massif_name_to_altitude_change_dynamic = {m: d for m, d in massif_name_to_altitude_change_dynamic.items()
+                                              if 0 < d < 3000}
+    visualizer.plot_abstract_fast(massif_name_to_altitude_change_dynamic, label='Altitude for the change of dynamic ()',
+                                  add_x_label=False, graduation=500)
+    # R2 score
+    visualizer.plot_abstract_fast(massif_name_to_r2_score, label='r2', graduation=0.1, add_x_label=False)
+
+
 def plot_snowfall_mean(altitude_to_visualizer: Dict[int, StudyVisualizerForMeanValues]):
     visualizer = list(altitude_to_visualizer.values())[0]
     # Plot the curve for the evolution of the mean
@@ -60,6 +84,8 @@ def plot_mean(altitude_to_visualizer: Dict[int, StudyVisualizerForMeanValues], d
                 values = [t.mean_value(year=year) for t in trend_tests]
             massif_name_to_linear_regression_result[massif_name] = fit_linear_regression(altitudes_massif, values)
             plot_values_against_altitudes(ax, altitudes_massif, massif_id, massif_name, moment, study, values, visualizer)
+    ax.legend(prop={'size': 7}, ncol=3)
+    visualizer.show_or_save_to_file(dpi=500, add_classic_title=False)
     plt.close()
 
     return [{m: t[i][0] if i == 0 else t[i]
@@ -76,4 +102,4 @@ def plot_values_against_altitudes(ax, altitudes, massif_id, massif_name, moment,
     # ax.set_ylim([lim_down, lim_up])
     ax.tick_params(axis='both', which='major', labelsize=13)
     visualizer.plot_name = plot_name
-    visualizer.show_or_save_to_file(dpi=500, add_classic_title=False)
+

projects/contrasting_trends_in_snow_loads/snowfall_versus_time_and_altitude/validation_plot.py

@@ -20,6 +20,7 @@ def validation_plot(altitude_to_visualizer: Dict[int, StudyVisualizerForMeanValu
     # Shoe plot with respect to the altitude.
     plot_shoe_relative_differences_distribution(altitude_to_relative_differences, altitudes)
     study_visualizer.show_or_save_to_file(add_classic_title=False, dpi=500)
+    plt.close()
 
 
 def plot_shoe_relative_differences_distribution(altitude_to_relative_differences, altitudes):