From bdcaf002c989e6e1346988fd255642ca6e28e19d Mon Sep 17 00:00:00 2001
From: Le Roux Erwan <erwan.le-roux@irstea.fr>
Date: Thu, 30 Apr 2020 11:57:04 +0200
Subject: [PATCH] [contrasting] add a validation for the change in the ratio
---
extreme_trend/abstract_gev_trend_test.py | 5 +-
.../main_snowfall_article.py | 7 +--
.../study_visualizer_for_mean_values.py | 50 +++++++++++++++----
.../validation_plot.py | 32 +++++++++---
4 files changed, 73 insertions(+), 21 deletions(-)
diff --git a/extreme_trend/abstract_gev_trend_test.py b/extreme_trend/abstract_gev_trend_test.py
index eac300d0..ce17e3c3 100644
--- a/extreme_trend/abstract_gev_trend_test.py
+++ b/extreme_trend/abstract_gev_trend_test.py
@@ -385,10 +385,9 @@ class AbstractGevTrendTest(object):
# Mean value
- @property
- def unconstrained_average_mean_value(self) -> float:
+ def unconstrained_average_mean_value(self, year_min, year_max) -> float:
mean_values = []
- for year in self.years:
+ for year in range(year_min, year_max+1):
mean = self.get_unconstrained_gev_params(year).mean
mean_values.append(mean)
average_mean_value = np.mean(mean_values)
diff --git a/projects/contrasting_trends_in_snow_loads/snowfall_versus_time_and_altitude/main_snowfall_article.py b/projects/contrasting_trends_in_snow_loads/snowfall_versus_time_and_altitude/main_snowfall_article.py
index 8c409087..f7137c76 100644
--- a/projects/contrasting_trends_in_snow_loads/snowfall_versus_time_and_altitude/main_snowfall_article.py
+++ b/projects/contrasting_trends_in_snow_loads/snowfall_versus_time_and_altitude/main_snowfall_article.py
@@ -72,9 +72,10 @@ def intermediate_result(altitudes, massif_names=None,
_ = compute_minimized_aic(visualizer)
# Plots
- validation_plot(altitude_to_visualizer)
- plot_snowfall_mean(altitude_to_visualizer)
- plot_snowfall_time_derivative_mean(altitude_to_visualizer)
+ validation_plot(altitude_to_visualizer, order_derivative=0)
+ validation_plot(altitude_to_visualizer, order_derivative=1)
+ # plot_snowfall_mean(altitude_to_visualizer)
+ # plot_snowfall_time_derivative_mean(altitude_to_visualizer)
def major_result():
diff --git a/projects/contrasting_trends_in_snow_loads/snowfall_versus_time_and_altitude/study_visualizer_for_mean_values.py b/projects/contrasting_trends_in_snow_loads/snowfall_versus_time_and_altitude/study_visualizer_for_mean_values.py
index 9fd53430..ca2cfce6 100644
--- a/projects/contrasting_trends_in_snow_loads/snowfall_versus_time_and_altitude/study_visualizer_for_mean_values.py
+++ b/projects/contrasting_trends_in_snow_loads/snowfall_versus_time_and_altitude/study_visualizer_for_mean_values.py
@@ -27,6 +27,13 @@ class StudyVisualizerForMeanValues(StudyVisualizerForNonStationaryTrends):
non_stationary_trend_test_to_marker, fit_method, select_only_acceptable_shape_parameter,
fit_gev_only_on_non_null_maxima, fit_only_time_series_with_ninety_percent_of_non_null_values)
+
+ def plot_abstract_fast(self, massif_name_to_value, label, graduation=10.0, cmap=plt.cm.coolwarm, add_x_label=True,
+ negative_and_positive_values=True):
+ super().plot_abstract(massif_name_to_value, label, label, self.fit_method, graduation, cmap, add_x_label, negative_and_positive_values)
+
+ # Study of the mean
+
@cached_property
def massif_name_to_empirical_mean(self):
massif_name_to_empirical_value = {}
@@ -36,23 +43,48 @@ class StudyVisualizerForMeanValues(StudyVisualizerForNonStationaryTrends):
@cached_property
def massif_name_to_model_mean(self):
- massif_name_to_parameter_value = {}
+ massif_name_to_model_value = {}
for massif_name, trend_test in self.massif_name_to_trend_test_that_minimized_aic.items():
- parameter_value = trend_test.unconstrained_average_mean_value
- massif_name_to_parameter_value[massif_name] = parameter_value
- return massif_name_to_parameter_value
+ parameter_value = trend_test.unconstrained_average_mean_value(self.study.year_min, self.study.year_max)
+ massif_name_to_model_value[massif_name] = parameter_value
+ return massif_name_to_model_value
@cached_property
def massif_name_to_relative_difference_for_mean(self):
massif_name_to_relative_difference = {}
for massif_name in self.massif_name_to_trend_test_that_minimized_aic.keys():
e = self.massif_name_to_empirical_mean[massif_name]
- p = self.massif_name_to_model_mean[massif_name]
- relative_diference = 100 * (p-e) / e
+ m = self.massif_name_to_model_mean[massif_name]
+ relative_diference = 100 * (m-e) / e
massif_name_to_relative_difference[massif_name] = relative_diference
return massif_name_to_relative_difference
- def plot_abstract_fast(self, massif_name_to_value, label, graduation=10.0, cmap=plt.cm.coolwarm, add_x_label=True,
- negative_and_positive_values=True):
- super().plot_abstract(massif_name_to_value, label, label, self.fit_method, graduation, cmap, add_x_label, negative_and_positive_values)
+ # Study of the change in the mean
+ @cached_property
+ def massif_name_to_change_ratio_in_empirical_mean(self):
+ massif_name_to_empirical_value = {}
+ for massif_name, maxima in self.study.massif_name_to_annual_maxima.items():
+ index = self.study.ordered_years.index(1989)
+ maxima_before, maxima_after = maxima[:index+1], maxima[index+1:]
+ massif_name_to_empirical_value[massif_name] = np.mean(maxima_after) / np.mean(maxima_before)
+ return massif_name_to_empirical_value
+
+ @cached_property
+ def massif_name_to_change_ratio_in_model_mean(self):
+ massif_name_to_parameter_value = {}
+ for massif_name, trend_test in self.massif_name_to_trend_test_that_minimized_aic.items():
+ model_mean_before = trend_test.unconstrained_average_mean_value(year_min=self.study.year_min, year_max=1989)
+ model_mean_after = trend_test.unconstrained_average_mean_value(year_min=1990, year_max=self.study.year_max)
+ massif_name_to_parameter_value[massif_name] = model_mean_after / model_mean_before
+ return massif_name_to_parameter_value
+
+ @cached_property
+ def massif_name_to_relative_difference_for_change_ratio_in_mean(self):
+ massif_name_to_relative_difference = {}
+ for massif_name in self.massif_name_to_trend_test_that_minimized_aic.keys():
+ e = self.massif_name_to_change_ratio_in_empirical_mean[massif_name]
+ m = self.massif_name_to_change_ratio_in_model_mean[massif_name]
+ relative_diference = 100 * (m-e) / e
+ massif_name_to_relative_difference[massif_name] = relative_diference
+ return massif_name_to_relative_difference
\ No newline at end of file
diff --git a/projects/contrasting_trends_in_snow_loads/snowfall_versus_time_and_altitude/validation_plot.py b/projects/contrasting_trends_in_snow_loads/snowfall_versus_time_and_altitude/validation_plot.py
index 2da8bda5..4847c939 100644
--- a/projects/contrasting_trends_in_snow_loads/snowfall_versus_time_and_altitude/validation_plot.py
+++ b/projects/contrasting_trends_in_snow_loads/snowfall_versus_time_and_altitude/validation_plot.py
@@ -10,34 +10,42 @@ from projects.ogorman.gorman_figures.figure1.study_visualizer_for_double_station
StudyVisualizerForReturnLevelChange
-def validation_plot(altitude_to_visualizer: Dict[int, StudyVisualizerForMeanValues]):
+def validation_plot(altitude_to_visualizer: Dict[int, StudyVisualizerForMeanValues], order_derivative=0):
# Plot the mean empirical, the mean parametric and the relative difference between the two
altitudes = list(altitude_to_visualizer.keys())
study_visualizer = list(altitude_to_visualizer.values())[0]
altitude_to_relative_differences = {}
+ if order_derivative == 0:
+ plot_function = plot_relative_difference_map_order_zero
+ else:
+ plot_function = plot_relative_difference_map_order_one
# Plot map for the repartition of the difference
for altitude, visualizer in altitude_to_visualizer.items():
- altitude_to_relative_differences[altitude] = plot_relative_difference_map(visualizer)
+ altitude_to_relative_differences[altitude] = plot_function(visualizer)
study_visualizer.show_or_save_to_file(add_classic_title=False, dpi=500)
# Shoe plot with respect to the altitude.
- plot_shoe_relative_differences_distribution(altitude_to_relative_differences, altitudes)
+ plot_shoe_relative_differences_distribution(altitude_to_relative_differences, altitudes, study_visualizer, order_derivative)
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):
+def plot_shoe_relative_differences_distribution(altitude_to_relative_differences, altitudes, visualizer, order_derivative):
+ study = visualizer.study
ax = plt.gca()
width = 150
ax.boxplot([altitude_to_relative_differences[a] for a in altitudes], positions=altitudes, widths=width)
ax.set_xlim([min(altitudes) - width, max(altitudes) + width])
- ylabel = 'Relative difference of the model mean w.r.t. the empirical mean (\%)'
+ moment = '' if order_derivative == 0 else 'time derivative of '
+ ylabel = 'Relative difference of the {} model mean \n' \
+ 'w.r.t. the {}empirical mean of {} (\%)'.format(moment, moment, SCM_STUDY_CLASS_TO_ABBREVIATION[type(study)])
+ visualizer.plot_trends = ylabel
ax.set_ylabel(ylabel)
ax.set_xlabel('Altitude (m)')
ax.legend()
ax.grid()
-def plot_relative_difference_map(visualizer: StudyVisualizerForMeanValues):
+def plot_relative_difference_map_order_zero(visualizer: StudyVisualizerForMeanValues):
study = visualizer.study
label = ' mean annual maxima of {} ({})'.format(SCM_STUDY_CLASS_TO_ABBREVIATION[type(study)], study.variable_unit)
visualizer.plot_abstract_fast(massif_name_to_value=visualizer.massif_name_to_empirical_mean,
@@ -48,3 +56,15 @@ def plot_relative_difference_map(visualizer: StudyVisualizerForMeanValues):
label='Relative difference of the model mean w.r.t. the empirical mean \n'
'for the ' + label, graduation=1)
return list(visualizer.massif_name_to_relative_difference_for_mean.values())
+
+def plot_relative_difference_map_order_one(visualizer: StudyVisualizerForMeanValues):
+ study = visualizer.study
+ label = ' time derivative of mean annual maxima of {} ({})'.format(SCM_STUDY_CLASS_TO_ABBREVIATION[type(study)], study.variable_unit)
+ visualizer.plot_abstract_fast(massif_name_to_value=visualizer.massif_name_to_change_ratio_in_empirical_mean,
+ label='Empirical' + label, negative_and_positive_values=False, graduation=0.5)
+ visualizer.plot_abstract_fast(massif_name_to_value=visualizer.massif_name_to_change_ratio_in_model_mean,
+ label='Model' + label, negative_and_positive_values=False, graduation=0.5)
+ visualizer.plot_abstract_fast(massif_name_to_value=visualizer.massif_name_to_relative_difference_for_change_ratio_in_mean,
+ label='Relative difference of the model mean w.r.t. the empirical mean \n'
+ 'for the ' + label, graduation=1)
+ return list(visualizer.massif_name_to_relative_difference_for_mean.values())
\ No newline at end of file
--
GitLab