diff --git a/extreme_fit/estimator/margin_estimator/utils.py b/extreme_fit/estimator/margin_estimator/utils.py
index 6867ff6b2d62b705abeb871eaf73dea43b6ae330..99f454f339c80075beea62f904fabcb51addb02b 100644
--- a/extreme_fit/estimator/margin_estimator/utils.py
+++ b/extreme_fit/estimator/margin_estimator/utils.py
@@ -14,12 +14,14 @@ from spatio_temporal_dataset.dataset.abstract_dataset import AbstractDataset
from spatio_temporal_dataset.spatio_temporal_observations.annual_maxima_observations import AnnualMaxima
-def fitted_linear_margin_estimator_short(model_class, dataset, fit_method, **model_kwargs) -> LinearMarginEstimator:
- return fitted_linear_margin_estimator(model_class, dataset.coordinates, dataset, None, fit_method, **model_kwargs)
+def fitted_linear_margin_estimator_short(model_class, dataset, fit_method, drop_duplicates=None, **model_kwargs) -> LinearMarginEstimator:
+ return fitted_linear_margin_estimator(model_class, dataset.coordinates, dataset, None, fit_method, drop_duplicates, **model_kwargs)
-def fitted_linear_margin_estimator(model_class, coordinates, dataset, starting_year, fit_method, **model_kwargs):
+def fitted_linear_margin_estimator(model_class, coordinates, dataset, starting_year, fit_method, drop_duplicates=None, **model_kwargs):
model = model_class(coordinates, starting_point=starting_year, fit_method=fit_method, **model_kwargs)
+ if drop_duplicates is not None:
+ model.drop_duplicates = drop_duplicates
estimator = LinearMarginEstimator(dataset, model)
estimator.fit()
return estimator
diff --git a/extreme_trend/ensemble_fit/independent_ensemble_fit/one_fold_fit_merge.py b/extreme_trend/ensemble_fit/independent_ensemble_fit/one_fold_fit_merge.py
index 48e012b6813f356150e4cf25e78f6f9559b03e90..b81a2e3bfaf4449efa1203e288d69229806fa283 100644
--- a/extreme_trend/ensemble_fit/independent_ensemble_fit/one_fold_fit_merge.py
+++ b/extreme_trend/ensemble_fit/independent_ensemble_fit/one_fold_fit_merge.py
@@ -22,14 +22,16 @@ class OneFoldFitMerge(OneFoldFit):
def get_moment(self, altitude, temporal_covariate, order=1):
return self.merge_function([o.get_moment(altitude, temporal_covariate, order) for o in self.one_fold_fit_list])
- def changes_of_moment(self, altitudes, order=1):
- all_changes = [o.changes_of_moment(altitudes, order) for o in self.one_fold_fit_list]
+ def changes_of_moment(self, altitudes, order=1, covariate_before=None, covariate_after=None):
+ all_changes = [o.changes_of_moment(altitudes, order, covariate_before, covariate_after)
+ for o in self.one_fold_fit_list]
merged_changes = list(self.merge_function(np.array(all_changes), axis=0))
assert len(all_changes[0]) == len(merged_changes)
return merged_changes
- def relative_changes_of_moment(self, altitudes, order=1):
- all_relative_changes = [o.relative_changes_of_moment(altitudes, order) for o in self.one_fold_fit_list]
+ def relative_changes_of_moment(self, altitudes, order=1, covariate_before=None, covariate_after=None):
+ all_relative_changes = [o.relative_changes_of_moment(altitudes, order, covariate_before, covariate_after)
+ for o in self.one_fold_fit_list]
merged_relative_changes = list(self.merge_function(np.array(all_relative_changes), axis=0))
assert len(all_relative_changes[0]) == len(merged_relative_changes)
return merged_relative_changes
diff --git a/extreme_trend/ensemble_fit/together_ensemble_fit/visualizer_non_stationary_ensemble.py b/extreme_trend/ensemble_fit/together_ensemble_fit/visualizer_non_stationary_ensemble.py
index 4d8b1999c9a5569b1cdb509f741288b79ac2fc77..fbe621fbae9506bd6d97ae69acce0da2d25ff8db 100644
--- a/extreme_trend/ensemble_fit/together_ensemble_fit/visualizer_non_stationary_ensemble.py
+++ b/extreme_trend/ensemble_fit/together_ensemble_fit/visualizer_non_stationary_ensemble.py
@@ -39,8 +39,14 @@ class VisualizerNonStationaryEnsemble(AltitudesStudiesVisualizerForNonStationary
dataset = studies.spatio_temporal_dataset(massif_name=massif_name, massif_altitudes=massif_altitudes)
df_coordinates_list.append(dataset.coordinates.df_coordinates(add_climate_informations=True))
df_maxima_gev_list.append(dataset.observations.df_maxima_gev)
- observations = AbstractSpatioTemporalObservations(df_maxima_gev=pd.concat(df_maxima_gev_list, axis=0))
- coordinates = AbstractCoordinates(df=pd.concat(df_coordinates_list, axis=0),
+
+ index = pd.RangeIndex(0, sum([len(df) for df in df_maxima_gev_list]))
+ df_maxima_gev = pd.concat(df_maxima_gev_list, axis=0)
+ df_maxima_gev.index = index
+ observations = AbstractSpatioTemporalObservations(df_maxima_gev=df_maxima_gev)
+ df = pd.concat(df_coordinates_list, axis=0)
+ df.index = index
+ coordinates = AbstractCoordinates(df=df,
slicer_class=type(dataset.slicer))
dataset = AbstractDataset(observations=observations, coordinates=coordinates)
return dataset
diff --git a/extreme_trend/ensemble_fit/visualizer_for_projection_ensemble.py b/extreme_trend/ensemble_fit/visualizer_for_projection_ensemble.py
index 04934a36ba4bfad3a0ce3ef144986b1212715b33..0a91aee99b01e747e24088891743dfea5712afcb 100644
--- a/extreme_trend/ensemble_fit/visualizer_for_projection_ensemble.py
+++ b/extreme_trend/ensemble_fit/visualizer_for_projection_ensemble.py
@@ -13,10 +13,13 @@ from extreme_data.meteo_france_data.scm_models_data.altitudes_studies import Alt
from extreme_trend.ensemble_fit.abstract_ensemble_fit import AbstractEnsembleFit
from extreme_trend.ensemble_fit.independent_ensemble_fit.independent_ensemble_fit import IndependentEnsembleFit
from extreme_trend.ensemble_fit.together_ensemble_fit.together_ensemble_fit import TogetherEnsembleFit
-from extreme_trend.one_fold_fit.altitude_group import get_altitude_class_from_altitudes
+from extreme_trend.one_fold_fit.altitude_group import get_altitude_class_from_altitudes, \
+ get_altitude_group_from_altitudes
from extreme_trend.one_fold_fit.plots.plot_histogram_altitude_studies import \
plot_histogram_all_trends_against_altitudes, plot_shoe_plot_changes_against_altitude
from extreme_trend.one_fold_fit.utils_altitude_studies_visualizer import compute_and_assign_max_abs
+from projects.projected_extreme_snowfall.results.plot_relative_change_in_return_level import \
+ plot_relative_dynamic_in_return_level
class VisualizerForProjectionEnsemble(object):
@@ -48,9 +51,9 @@ class VisualizerForProjectionEnsemble(object):
assert len(years) == 2, years
# Load all studies
- altitude_class_to_gcm_couple_to_studies = OrderedDict()
+ altitude_group_to_gcm_couple_to_studies = OrderedDict()
for altitudes in altitudes_list:
- altitude_class = get_altitude_class_from_altitudes(altitudes)
+ altitude_group = get_altitude_group_from_altitudes(altitudes)
gcm_rcm_couple_to_studies = {}
for gcm_rcm_couple in gcm_rcm_couples:
if gcm_to_year_min_and_year_max is None:
@@ -69,11 +72,11 @@ class VisualizerForProjectionEnsemble(object):
gcm_rcm_couple_to_studies[gcm_rcm_couple] = studies
if len(gcm_rcm_couple_to_studies) == 0:
print('No valid studies for the following couples:', self.gcm_rcm_couples)
- altitude_class_to_gcm_couple_to_studies[altitude_class] = gcm_rcm_couple_to_studies
+ altitude_group_to_gcm_couple_to_studies[altitude_group] = gcm_rcm_couple_to_studies
# Load ensemble fit
- self.altitude_class_to_ensemble_class_to_ensemble_fit = OrderedDict()
- for altitude_class, gcm_rcm_couple_to_studies in altitude_class_to_gcm_couple_to_studies.items():
+ self.altitude_group_to_ensemble_class_to_ensemble_fit = OrderedDict()
+ for altitude_group, gcm_rcm_couple_to_studies in altitude_group_to_gcm_couple_to_studies.items():
ensemble_class_to_ensemble_fit = {}
for ensemble_fit_class in ensemble_fit_classes:
ensemble_fit = ensemble_fit_class(massif_names, gcm_rcm_couple_to_studies, model_classes,
@@ -82,7 +85,7 @@ class VisualizerForProjectionEnsemble(object):
confidence_interval_based_on_delta_method,
remove_physically_implausible_models)
ensemble_class_to_ensemble_fit[ensemble_fit_class] = ensemble_fit
- self.altitude_class_to_ensemble_class_to_ensemble_fit[altitude_class] = ensemble_class_to_ensemble_fit
+ self.altitude_group_to_ensemble_class_to_ensemble_fit[altitude_group] = ensemble_class_to_ensemble_fit
@property
def has_elevation_non_stationarity(self):
@@ -99,7 +102,7 @@ class VisualizerForProjectionEnsemble(object):
plot_shoe_plot_changes_against_altitude(self.massif_names, visualizer_list, relative=relative,
with_significance=with_significance)
else:
- print('here plot for one altitude visualizer')
+ plot_relative_dynamic_in_return_level(self.massif_names, visualizer_list)
def plot(self):
# Set limit for the plot
@@ -142,7 +145,7 @@ class VisualizerForProjectionEnsemble(object):
"""Return the ordered ensemble fit for a given ensemble class (in the order of the altitudes)"""
return [ensemble_class_to_ensemble_fit[ensemble_class]
for ensemble_class_to_ensemble_fit
- in self.altitude_class_to_ensemble_class_to_ensemble_fit.values()]
+ in self.altitude_group_to_ensemble_class_to_ensemble_fit.values()]
def plot_preliminary_first_part(self):
if self.massif_names is None:
diff --git a/extreme_trend/ensemble_fit/visualizer_for_sensitivity.py b/extreme_trend/ensemble_fit/visualizer_for_sensitivity.py
index 37150b16ed99b34e9ad2a74b9358110c904cc673..a33668bb854abeb3f1d0c76d4f94e8f004dc8c80 100644
--- a/extreme_trend/ensemble_fit/visualizer_for_sensitivity.py
+++ b/extreme_trend/ensemble_fit/visualizer_for_sensitivity.py
@@ -235,12 +235,12 @@ class VisualizerForSensivity(object):
merge_visualizer_str):
if merge_visualizer_str in [AbstractEnsembleFit.Median_merge, AbstractEnsembleFit.Mean_merge]:
independent_ensemble_fit = \
- visualizer_projection.altitude_class_to_ensemble_class_to_ensemble_fit[altitude_class][
+ visualizer_projection.altitude_group_to_ensemble_class_to_ensemble_fit[altitude_class][
IndependentEnsembleFit]
merge_visualizer = independent_ensemble_fit.merge_function_name_to_visualizer[merge_visualizer_str]
else:
together_ensemble_fit = \
- visualizer_projection.altitude_class_to_ensemble_class_to_ensemble_fit[altitude_class][
+ visualizer_projection.altitude_group_to_ensemble_class_to_ensemble_fit[altitude_class][
TogetherEnsembleFit]
merge_visualizer = together_ensemble_fit.visualizer
merge_visualizer.studies.study.gcm_rcm_couple = (merge_visualizer_str, "merge")
diff --git a/extreme_trend/one_fold_fit/one_fold_fit.py b/extreme_trend/one_fold_fit/one_fold_fit.py
index bb70613465a5da3c8e31878f5d89e2eb49dc274f..19396a19dc5e25671e3758096f842d1f10260452 100644
--- a/extreme_trend/one_fold_fit/one_fold_fit.py
+++ b/extreme_trend/one_fold_fit/one_fold_fit.py
@@ -77,7 +77,8 @@ class OneFoldFit(object):
return fitted_linear_margin_estimator_short(model_class=model_class,
dataset=dataset,
fit_method=self.fit_method,
- temporal_covariate_for_fit=self.temporal_covariate_for_fit)
+ temporal_covariate_for_fit=self.temporal_covariate_for_fit,
+ drop_duplicates=False)
@classmethod
def get_moment_str(cls, order):
@@ -115,11 +116,15 @@ class OneFoldFit(object):
def relative_change_in_return_level_for_reference_altitude(self) -> float:
return self.relative_changes_of_moment(altitudes=[self.altitude_plot], order=None)[0]
- def changes_of_moment(self, altitudes, order=1):
+ def changes_of_moment(self, altitudes, order=1, covariate_before=None, covariate_after=None):
+ if covariate_before is None:
+ covariate_before = self.covariate_before
+ if covariate_after is None:
+ covariate_after = self.covariate_after
changes = []
for altitude in altitudes:
- mean_after = self.get_moment(altitude, self.covariate_after, order)
- mean_before = self.get_moment(altitude, self.covariate_before, order)
+ mean_after = self.get_moment(altitude, covariate_after, order)
+ mean_before = self.get_moment(altitude, covariate_before, order)
change = mean_after - mean_before
changes.append(change)
return changes
@@ -151,11 +156,11 @@ class OneFoldFit(object):
**d_temperature)
return s
- def relative_changes_of_moment(self, altitudes, order=1):
+ def relative_changes_of_moment(self, altitudes, order=1, covariate_before=None, covariate_after=None):
relative_changes = []
for altitude in altitudes:
- mean_after = self.get_moment(altitude, self.covariate_after, order)
- mean_before = self.get_moment(altitude, self.covariate_before, order)
+ mean_after = self.get_moment(altitude, covariate_after, order)
+ mean_before = self.get_moment(altitude, covariate_before, order)
relative_change = 100 * (mean_after - mean_before) / mean_before
relative_changes.append(relative_change)
return relative_changes
diff --git a/projects/projected_extreme_snowfall/results/main_projections_ensemble.py b/projects/projected_extreme_snowfall/results/main_projections_ensemble.py
index 2069ed808f519b32c04d49d8d7d963ef9b3958f0..79cf3832b21f6a4c87ca1ba6819780731c0e47de 100644
--- a/projects/projected_extreme_snowfall/results/main_projections_ensemble.py
+++ b/projects/projected_extreme_snowfall/results/main_projections_ensemble.py
@@ -43,7 +43,7 @@ from extreme_data.meteo_france_data.scm_models_data.utils import Season
def main():
start = time.time()
study_class = AdamontSnowfall
- ensemble_fit_classes = [IndependentEnsembleFit, TogetherEnsembleFit][:1]
+ ensemble_fit_classes = [IndependentEnsembleFit, TogetherEnsembleFit][1:]
temporal_covariate_for_fit = [TimeTemporalCovariate,
AnomalyTemperatureWithSplineTemporalCovariate][1]
set_seed_for_test()
diff --git a/projects/projected_extreme_snowfall/results/plot_relative_change_in_return_level.py b/projects/projected_extreme_snowfall/results/plot_relative_change_in_return_level.py
new file mode 100644
index 0000000000000000000000000000000000000000..34a9db8d44d4efab239b94a86382534ac3318823
--- /dev/null
+++ b/projects/projected_extreme_snowfall/results/plot_relative_change_in_return_level.py
@@ -0,0 +1,42 @@
+from typing import List
+import matplotlib.pyplot as plt
+
+import numpy as np
+
+from extreme_trend.one_fold_fit.altitudes_studies_visualizer_for_non_stationary_models import \
+ AltitudesStudiesVisualizerForNonStationaryModels
+from spatio_temporal_dataset.coordinates.temporal_coordinates.temperature_covariate import \
+ AnomalyTemperatureWithSplineTemporalCovariate
+
+
+def plot_relative_dynamic_in_return_level(massif_names, visualizer_list: List[
+ AltitudesStudiesVisualizerForNonStationaryModels]):
+ ax = plt.gca()
+ print(len(visualizer_list))
+ visualizer = visualizer_list[0]
+ assert len(massif_names) == 1
+ assert visualizer.temporal_covariate_for_fit is AnomalyTemperatureWithSplineTemporalCovariate
+ massif_name = massif_names[0]
+ for v in visualizer_list:
+ plot_curve(ax, massif_name, v)
+
+ ax.set_xlabel('Anomaly of global temperature w.r.t. pre-industrial levels (K)')
+ ax.set_ylabel('Relative change in return levels (\%)')
+
+ ax.legend()
+ visualizer.plot_name = 'dynamic of return level'
+ visualizer.show_or_save_to_file(add_classic_title=False, no_title=True)
+
+ plt.close()
+
+
+def plot_curve(ax, massif_name, visualizer: AltitudesStudiesVisualizerForNonStationaryModels):
+ temperatures_list = np.linspace(1, 4, num=4)
+ one_fold_fit = visualizer.massif_name_to_one_fold_fit[massif_name]
+ return_levels = [one_fold_fit.relative_changes_of_moment([None], order=None,
+ covariate_before=1,
+ covariate_after=t)[0] for t in temperatures_list]
+ label = '{} m'.format(visualizer.altitude_group.reference_altitude)
+ ax.plot(temperatures_list, return_levels, label=label)
+
+