From b9ba4e91f1a908f1fe3a22ef26950f16cfecc7b9 Mon Sep 17 00:00:00 2001
From: Le Roux Erwan <erwan.le-roux@irstea.fr>
Date: Thu, 28 Nov 2019 17:16:33 +0100
Subject: [PATCH] [PAPER 1] refactor folder names in the paper snow load
folder. add mani gelman convergence test. add the possibility to change the
number of iterations for the bayesian fevd estimation
---
.../main_study_visualizer.py | 2 +-
.../__init__.py | 0
.../gelman_convergence_test.py | 5 +++
.../main_bayesian_mcmc.py | 10 +++--
.../main_gelman_convergence_test.py | 39 +++++++++++++++++++
.../__init__.py | 0
.../main_mle_diagnosis.py | 0
.../__init__.py | 0
.../crocus_study_comparison_with_eurocode.py | 0
.../main_comparison_with_eurocode.py | 2 +-
.../paper_past_snow_loads/paper_utils.py | 19 +++++++++
.../eurocode_visualizer.py | 7 +---
.../main_result_trends_and_return_levels.py | 28 ++++---------
...dy_visualizer_for_non_stationary_trends.py | 28 ++++++++-----
.../trend_analysis/univariate_test/utils.py | 4 +-
.../abstract_temporal_linear_margin_model.py | 6 ++-
16 files changed, 106 insertions(+), 44 deletions(-)
rename experiment/paper_past_snow_loads/{result_data_comparison_with_eurocode => check_mcmc_convergence_for_return_levels}/__init__.py (100%)
create mode 100644 experiment/paper_past_snow_loads/check_mcmc_convergence_for_return_levels/gelman_convergence_test.py
rename experiment/paper_past_snow_loads/{result_mcmc_check_for_return_levels => check_mcmc_convergence_for_return_levels}/main_bayesian_mcmc.py (92%)
create mode 100644 experiment/paper_past_snow_loads/check_mcmc_convergence_for_return_levels/main_gelman_convergence_test.py
rename experiment/paper_past_snow_loads/{result_mcmc_check_for_return_levels => check_mle_convergence_for_trends}/__init__.py (100%)
rename experiment/paper_past_snow_loads/{result_mle_check_for_trends => check_mle_convergence_for_trends}/main_mle_diagnosis.py (100%)
rename experiment/paper_past_snow_loads/{result_mle_check_for_trends => discussion_data_comparison_with_eurocode}/__init__.py (100%)
rename experiment/paper_past_snow_loads/{result_data_comparison_with_eurocode => discussion_data_comparison_with_eurocode}/crocus_study_comparison_with_eurocode.py (100%)
rename experiment/paper_past_snow_loads/{result_data_comparison_with_eurocode => discussion_data_comparison_with_eurocode}/main_comparison_with_eurocode.py (96%)
create mode 100644 experiment/paper_past_snow_loads/paper_utils.py
rename experiment/paper_past_snow_loads/{result_trends_and_return_levels => }/study_visualizer_for_non_stationary_trends.py (91%)
diff --git a/experiment/meteo_france_data/scm_models_data/visualization/study_visualization/main_study_visualizer.py b/experiment/meteo_france_data/scm_models_data/visualization/study_visualization/main_study_visualizer.py
index 147dc88f..67ac17dd 100644
--- a/experiment/meteo_france_data/scm_models_data/visualization/study_visualization/main_study_visualizer.py
+++ b/experiment/meteo_france_data/scm_models_data/visualization/study_visualization/main_study_visualizer.py
@@ -3,7 +3,7 @@ from typing import List
from experiment.meteo_france_data.scm_models_data.visualization.study_visualization.study_visualizer import \
StudyVisualizer
-from experiment.paper_past_snow_loads.result_data_comparison_with_eurocode.crocus_study_comparison_with_eurocode import \
+from experiment.paper_past_snow_loads.discussion_data_comparison_with_eurocode.crocus_study_comparison_with_eurocode import \
CrocusDifferenceSnowLoad, \
CrocusSnowDensityAtMaxofSwe, CrocusDifferenceSnowLoadRescaledAndEurocodeToSeeSynchronization, \
CrocusSnowDepthDifference, CrocusSnowDepthAtMaxofSwe
diff --git a/experiment/paper_past_snow_loads/result_data_comparison_with_eurocode/__init__.py b/experiment/paper_past_snow_loads/check_mcmc_convergence_for_return_levels/__init__.py
similarity index 100%
rename from experiment/paper_past_snow_loads/result_data_comparison_with_eurocode/__init__.py
rename to experiment/paper_past_snow_loads/check_mcmc_convergence_for_return_levels/__init__.py
diff --git a/experiment/paper_past_snow_loads/check_mcmc_convergence_for_return_levels/gelman_convergence_test.py b/experiment/paper_past_snow_loads/check_mcmc_convergence_for_return_levels/gelman_convergence_test.py
new file mode 100644
index 00000000..7df11604
--- /dev/null
+++ b/experiment/paper_past_snow_loads/check_mcmc_convergence_for_return_levels/gelman_convergence_test.py
@@ -0,0 +1,5 @@
+
+def compute_gelman_convergence_value(non_null_years_and_maxima, mcmc_iterations, model_class, nb_chains):
+ return 1.0
+
+# rbeta(10000, 6, 9) - 0.5
\ No newline at end of file
diff --git a/experiment/paper_past_snow_loads/result_mcmc_check_for_return_levels/main_bayesian_mcmc.py b/experiment/paper_past_snow_loads/check_mcmc_convergence_for_return_levels/main_bayesian_mcmc.py
similarity index 92%
rename from experiment/paper_past_snow_loads/result_mcmc_check_for_return_levels/main_bayesian_mcmc.py
rename to experiment/paper_past_snow_loads/check_mcmc_convergence_for_return_levels/main_bayesian_mcmc.py
index f8057130..1a2a7d47 100644
--- a/experiment/paper_past_snow_loads/result_mcmc_check_for_return_levels/main_bayesian_mcmc.py
+++ b/experiment/paper_past_snow_loads/check_mcmc_convergence_for_return_levels/main_bayesian_mcmc.py
@@ -1,7 +1,7 @@
import seaborn as sns
import matplotlib.pyplot as plt
-from experiment.meteo_france_data.scm_models_data.crocus.crocus import CrocusSwe3Days
+from experiment.meteo_france_data.scm_models_data.crocus.crocus import CrocusSwe3Days, CrocusSnowLoadTotal
from experiment.meteo_france_data.scm_models_data.visualization.utils import create_adjusted_axes
from experiment.trend_analysis.univariate_test.utils import load_temporal_coordinates_and_dataset, \
fitted_linear_margin_estimator
@@ -31,7 +31,8 @@ def main_drawing_bayesian():
colors = ['r', 'g', 'b', 'tab:orange']
gev_param_name_to_color = dict(zip(GevParams.PARAM_NAMES, colors))
gev_param_name_to_ax = dict(
- zip(GevParams.PARAM_NAMES, [ax_trajectories, ax_trajectories, ax_trajectories_inverted]))
+ zip(GevParams.PARAM_NAMES, [ax_trajectories, ax_trajectories, ax_trajectories]))
+ # zip(GevParams.PARAM_NAMES, [ax_trajectories, ax_trajectories, ax_trajectories_inverted]))
gev_param_name_to_label = {n: GevParams.greek_letter_from_gev_param_name(n) for n in GevParams.PARAM_NAMES}
for j, gev_param_name in enumerate(GevParams.PARAM_NAMES[:]):
label = gev_param_name_to_label[gev_param_name]
@@ -77,13 +78,14 @@ def main_drawing_bayesian():
def get_return_level_bayesian_example():
- maxima, years = CrocusSwe3Days(altitude=1800).annual_maxima_and_years('Vercors')
+ maxima, years = CrocusSnowLoadTotal(altitude=1800).annual_maxima_and_years('Vercors')
# plt.plot(years, maxima)
# plt.show()
model_class = StationaryTemporalModel
coordinates, dataset = load_temporal_coordinates_and_dataset(maxima, years)
fitted_estimator = fitted_linear_margin_estimator(model_class, coordinates, dataset, starting_year=1958,
- fit_method=TemporalMarginFitMethod.extremes_fevd_bayesian)
+ fit_method=TemporalMarginFitMethod.extremes_fevd_bayesian,
+ nb_iterations_for_bayesian_fit=100000)
return_level_bayesian = ExtractEurocodeReturnLevelFromMyBayesianExtremes(estimator=fitted_estimator,
ci_method=ConfidenceIntervalMethodFromExtremes.my_bayes,
temporal_covariate=2017)
diff --git a/experiment/paper_past_snow_loads/check_mcmc_convergence_for_return_levels/main_gelman_convergence_test.py b/experiment/paper_past_snow_loads/check_mcmc_convergence_for_return_levels/main_gelman_convergence_test.py
new file mode 100644
index 00000000..d435acf8
--- /dev/null
+++ b/experiment/paper_past_snow_loads/check_mcmc_convergence_for_return_levels/main_gelman_convergence_test.py
@@ -0,0 +1,39 @@
+import pandas as pd
+from experiment.paper_past_snow_loads.paper_utils import paper_altitudes, paper_study_classes, \
+ load_altitude_to_visualizer
+from extreme_fit.model.margin_model.linear_margin_model.temporal_linear_margin_models import StationaryTemporalModel
+
+
+def gelman_convergence_test(mcmc_iterations, model_class, altitudes, study_class, nb_chains=3, massif_names=None):
+ """ Test if a given number of MCMC iterations is enough for the convergence of each parameter of the model_class
+ for every time series with non-null values present in the study_classes
+
+ Ideally, return a DataFrame with altitude as columns, and massif name as index that contain the R score
+ Then we could compute the max R that should ideally be below 1.2
+ """
+
+ altitude_to_visualizer = load_altitude_to_visualizer(altitudes, massif_names=massif_names,
+ non_stationary_uncertainty=None,
+ study_class=study_class, uncertainty_methods=None)
+ altitude_to_d = {}
+ for altitude, vizu in altitude_to_visualizer.items():
+ altitude_to_d[altitude] = vizu.massif_name_to_gelman_convergence_value(mcmc_iterations, model_class, nb_chains)
+ massif_names = list(altitude_to_d[altitudes[0]].keys())
+ df = pd.DataFrame(altitude_to_d, index=massif_names, columns=altitudes)
+ return df
+
+
+
+"""
+test gelman for the 4 types of models
+and the for the 3 variables considered: GSL, GSL from eurocode, GLS in 3 days
+"""
+
+if __name__ == '__main__':
+ mcmc_iterations = 1000
+ df = gelman_convergence_test(mcmc_iterations=mcmc_iterations, altitudes=paper_altitudes[:1],
+ study_class=paper_study_classes[0], model_class=StationaryTemporalModel,
+ massif_names=['Chartreuse'])
+ print(mcmc_iterations)
+ print(df.head())
+ print('Overall maxima:', df.max().max())
diff --git a/experiment/paper_past_snow_loads/result_mcmc_check_for_return_levels/__init__.py b/experiment/paper_past_snow_loads/check_mle_convergence_for_trends/__init__.py
similarity index 100%
rename from experiment/paper_past_snow_loads/result_mcmc_check_for_return_levels/__init__.py
rename to experiment/paper_past_snow_loads/check_mle_convergence_for_trends/__init__.py
diff --git a/experiment/paper_past_snow_loads/result_mle_check_for_trends/main_mle_diagnosis.py b/experiment/paper_past_snow_loads/check_mle_convergence_for_trends/main_mle_diagnosis.py
similarity index 100%
rename from experiment/paper_past_snow_loads/result_mle_check_for_trends/main_mle_diagnosis.py
rename to experiment/paper_past_snow_loads/check_mle_convergence_for_trends/main_mle_diagnosis.py
diff --git a/experiment/paper_past_snow_loads/result_mle_check_for_trends/__init__.py b/experiment/paper_past_snow_loads/discussion_data_comparison_with_eurocode/__init__.py
similarity index 100%
rename from experiment/paper_past_snow_loads/result_mle_check_for_trends/__init__.py
rename to experiment/paper_past_snow_loads/discussion_data_comparison_with_eurocode/__init__.py
diff --git a/experiment/paper_past_snow_loads/result_data_comparison_with_eurocode/crocus_study_comparison_with_eurocode.py b/experiment/paper_past_snow_loads/discussion_data_comparison_with_eurocode/crocus_study_comparison_with_eurocode.py
similarity index 100%
rename from experiment/paper_past_snow_loads/result_data_comparison_with_eurocode/crocus_study_comparison_with_eurocode.py
rename to experiment/paper_past_snow_loads/discussion_data_comparison_with_eurocode/crocus_study_comparison_with_eurocode.py
diff --git a/experiment/paper_past_snow_loads/result_data_comparison_with_eurocode/main_comparison_with_eurocode.py b/experiment/paper_past_snow_loads/discussion_data_comparison_with_eurocode/main_comparison_with_eurocode.py
similarity index 96%
rename from experiment/paper_past_snow_loads/result_data_comparison_with_eurocode/main_comparison_with_eurocode.py
rename to experiment/paper_past_snow_loads/discussion_data_comparison_with_eurocode/main_comparison_with_eurocode.py
index f266e948..3a54757d 100644
--- a/experiment/paper_past_snow_loads/result_data_comparison_with_eurocode/main_comparison_with_eurocode.py
+++ b/experiment/paper_past_snow_loads/discussion_data_comparison_with_eurocode/main_comparison_with_eurocode.py
@@ -6,7 +6,7 @@ from experiment.meteo_france_data.scm_models_data.visualization.study_visualizat
StudyVisualizer
import matplotlib.pyplot as plt
-from experiment.paper_past_snow_loads.result_data_comparison_with_eurocode.crocus_study_comparison_with_eurocode import \
+from experiment.paper_past_snow_loads.discussion_data_comparison_with_eurocode.crocus_study_comparison_with_eurocode import \
CrocusDifferenceSnowLoad, \
CrocusSnowDensityAtMaxofSwe, CrocusDifferenceSnowLoadRescaledAndEurocodeToSeeSynchronization, \
CrocusSnowDepthAtMaxofSwe, CrocusSnowDepthDifference
diff --git a/experiment/paper_past_snow_loads/paper_utils.py b/experiment/paper_past_snow_loads/paper_utils.py
new file mode 100644
index 00000000..7af16183
--- /dev/null
+++ b/experiment/paper_past_snow_loads/paper_utils.py
@@ -0,0 +1,19 @@
+from collections import OrderedDict
+
+from experiment.meteo_france_data.scm_models_data.crocus.crocus import CrocusSnowLoadTotal, CrocusSnowLoadEurocode, \
+ CrocusSnowLoad3Days
+from experiment.paper_past_snow_loads.study_visualizer_for_non_stationary_trends import \
+ StudyVisualizerForNonStationaryTrends
+
+paper_altitudes = [300, 600, 900, 1200, 1500, 1800, 2100, 2400, 2700]
+paper_study_classes = [CrocusSnowLoadTotal, CrocusSnowLoadEurocode, CrocusSnowLoad3Days]
+
+
+def load_altitude_to_visualizer(altitudes, massif_names, non_stationary_uncertainty, study_class, uncertainty_methods):
+ altitude_to_visualizer = OrderedDict()
+ for altitude in altitudes:
+ altitude_to_visualizer[altitude] = StudyVisualizerForNonStationaryTrends(
+ study=study_class(altitude=altitude), multiprocessing=True, save_to_file=True,
+ uncertainty_massif_names=massif_names, uncertainty_methods=uncertainty_methods,
+ non_stationary_contexts=non_stationary_uncertainty)
+ return altitude_to_visualizer
diff --git a/experiment/paper_past_snow_loads/result_trends_and_return_levels/eurocode_visualizer.py b/experiment/paper_past_snow_loads/result_trends_and_return_levels/eurocode_visualizer.py
index 73921d0a..5d6bbedf 100644
--- a/experiment/paper_past_snow_loads/result_trends_and_return_levels/eurocode_visualizer.py
+++ b/experiment/paper_past_snow_loads/result_trends_and_return_levels/eurocode_visualizer.py
@@ -1,15 +1,12 @@
-from typing import Dict, List, Tuple
+from typing import Dict
-import matplotlib.pyplot as plt
import numpy as np
from experiment.eurocode_data.utils import EUROCODE_RETURN_LEVEL_STR, EUROCODE_ALTITUDES
-from experiment.paper_past_snow_loads.result_trends_and_return_levels.study_visualizer_for_non_stationary_trends import \
+from experiment.paper_past_snow_loads.study_visualizer_for_non_stationary_trends import \
StudyVisualizerForNonStationaryTrends
from extreme_fit.model.result_from_model_fit.result_from_extremes.abstract_extract_eurocode_return_level import \
AbstractExtractEurocodeReturnLevel
-from extreme_fit.model.result_from_model_fit.result_from_extremes.eurocode_return_level_uncertainties import \
- EurocodeConfidenceIntervalFromExtremes
from experiment.eurocode_data.massif_name_to_departement import massif_name_to_eurocode_region
from experiment.meteo_france_data.scm_models_data.visualization.utils import create_adjusted_axes
from root_utils import get_display_name_from_object_type
diff --git a/experiment/paper_past_snow_loads/result_trends_and_return_levels/main_result_trends_and_return_levels.py b/experiment/paper_past_snow_loads/result_trends_and_return_levels/main_result_trends_and_return_levels.py
index b0bfdf88..4476a9ff 100644
--- a/experiment/paper_past_snow_loads/result_trends_and_return_levels/main_result_trends_and_return_levels.py
+++ b/experiment/paper_past_snow_loads/result_trends_and_return_levels/main_result_trends_and_return_levels.py
@@ -1,21 +1,14 @@
-from collections import OrderedDict
-import os.path as op
-
import matplotlib as mpl
import matplotlib.pyplot as plt
-from experiment.meteo_france_data.plot.create_shifted_cmap import get_shifted_map
-from experiment.meteo_france_data.scm_models_data.visualization.study_visualization.study_visualizer import \
- StudyVisualizer
+from experiment.paper_past_snow_loads.paper_utils import paper_study_classes, paper_altitudes
from experiment.paper_past_snow_loads.result_trends_and_return_levels.eurocode_visualizer import \
plot_uncertainty_massifs
-from experiment.meteo_france_data.scm_models_data.crocus.crocus import CrocusSnowLoadTotal, CrocusSnowLoadEurocode, \
- CrocusSnowLoad3Days
-from experiment.paper_past_snow_loads.result_trends_and_return_levels.study_visualizer_for_non_stationary_trends import \
+from experiment.meteo_france_data.scm_models_data.crocus.crocus import CrocusSnowLoadTotal
+from experiment.paper_past_snow_loads.study_visualizer_for_non_stationary_trends import \
StudyVisualizerForNonStationaryTrends
from extreme_fit.model.result_from_model_fit.result_from_extremes.confidence_interval_method import \
ConfidenceIntervalMethodFromExtremes
-from root_utils import VERSION_TIME
mpl.rcParams['text.usetex'] = True
mpl.rcParams['text.latex.preamble'] = [r'\usepackage{amsmath}']
@@ -42,12 +35,8 @@ def intermediate_result(altitudes, massif_names=None,
:return:
"""
# Load altitude to visualizer
- altitude_to_visualizer = OrderedDict()
- for altitude in altitudes:
- altitude_to_visualizer[altitude] = StudyVisualizerForNonStationaryTrends(
- study=study_class(altitude=altitude), multiprocessing=True, save_to_file=True,
- uncertainty_massif_names=massif_names, uncertainty_methods=uncertainty_methods,
- non_stationary_contexts=non_stationary_uncertainty)
+ altitude_to_visualizer = load_altitude_to_visualizer(altitudes, massif_names, non_stationary_uncertainty,
+ study_class, uncertainty_methods)
# Plot trends
max_abs_tdrl = max([visualizer.max_abs_tdrl for visualizer in altitude_to_visualizer.values()])
for visualizer in altitude_to_visualizer.values():
@@ -57,15 +46,14 @@ def intermediate_result(altitudes, massif_names=None,
return altitude_to_visualizer
+
def major_result():
- altitudes = [300, 600, 900, 1200, 1500, 1800, 2100, 2400, 2700]
uncertainty_methods = [ConfidenceIntervalMethodFromExtremes.my_bayes,
ConfidenceIntervalMethodFromExtremes.ci_mle][:]
massif_names = None
non_stationary_uncertainty = [False, True][:]
- study_classes = [CrocusSnowLoadTotal, CrocusSnowLoadEurocode, CrocusSnowLoad3Days]
- for study_class in study_classes[2:]:
- intermediate_result(altitudes, massif_names, non_stationary_uncertainty, uncertainty_methods, study_class)
+ for study_class in paper_study_classes[2:]:
+ intermediate_result(paper_altitudes, massif_names, non_stationary_uncertainty, uncertainty_methods, study_class)
if __name__ == '__main__':
diff --git a/experiment/paper_past_snow_loads/result_trends_and_return_levels/study_visualizer_for_non_stationary_trends.py b/experiment/paper_past_snow_loads/study_visualizer_for_non_stationary_trends.py
similarity index 91%
rename from experiment/paper_past_snow_loads/result_trends_and_return_levels/study_visualizer_for_non_stationary_trends.py
rename to experiment/paper_past_snow_loads/study_visualizer_for_non_stationary_trends.py
index c8ee7c7b..57d8456a 100644
--- a/experiment/paper_past_snow_loads/result_trends_and_return_levels/study_visualizer_for_non_stationary_trends.py
+++ b/experiment/paper_past_snow_loads/study_visualizer_for_non_stationary_trends.py
@@ -11,6 +11,8 @@ from experiment.meteo_france_data.plot.create_shifted_cmap import get_shifted_ma
from experiment.meteo_france_data.scm_models_data.abstract_study import AbstractStudy
from experiment.meteo_france_data.scm_models_data.visualization.study_visualization.study_visualizer import \
StudyVisualizer
+from experiment.paper_past_snow_loads.check_mcmc_convergence_for_return_levels.gelman_convergence_test import \
+ compute_gelman_convergence_value
from experiment.trend_analysis.abstract_score import MeanScore
from experiment.trend_analysis.univariate_test.abstract_gev_trend_test import AbstractGevTrendTest
from experiment.trend_analysis.univariate_test.gev_trend_test_one_parameter import GevScaleTrendTest, \
@@ -147,7 +149,7 @@ class StudyVisualizerForNonStationaryTrends(StudyVisualizer):
def massif_name_to_tdrl_value(self):
return {m: t.time_derivative_of_return_level for m, t in
self.massif_name_to_minimized_aic_non_stationary_trend_test.items()}
-
+
@cached_property
def cmap(self):
return get_shifted_map(-self._max_abs_tdrl, self._max_abs_tdrl)
@@ -169,8 +171,8 @@ class StudyVisualizerForNonStationaryTrends(StudyVisualizer):
# Part 2 - Uncertainty return level plot
- def massif_name_to_model_class(self, massif_name, non_stationary_model):
- if not non_stationary_model:
+ def massif_name_and_non_stationary_context_to_model_class(self, massif_name, non_stationary_context):
+ if not non_stationary_context:
return StationaryTemporalModel
else:
return self.massif_name_to_minimized_aic_non_stationary_trend_test[massif_name].unconstrained_model_class
@@ -179,16 +181,12 @@ class StudyVisualizerForNonStationaryTrends(StudyVisualizer):
def nb_contexts(self):
return len(self.non_stationary_contexts)
- @property
- def nb_uncertainty_method(self):
- return len(self.uncertainty_methods)
-
- def all_massif_name_to_eurocode_uncertainty_for_minimized_aic_model_class(self, ci_method, non_stationary_model) \
+ def all_massif_name_to_eurocode_uncertainty_for_minimized_aic_model_class(self, ci_method, non_stationary_context) \
-> Dict[str, Tuple[int, EurocodeConfidenceIntervalFromExtremes]]:
# Compute for the uncertainty massif names
arguments = [
[self.massif_name_to_non_null_years_and_maxima[m],
- self.massif_name_to_model_class(m, non_stationary_model),
+ self.massif_name_and_non_stationary_context_to_model_class(m, non_stationary_context),
ci_method, self.effective_temporal_covariate,
self.massif_name_to_eurocode_quantile_level_in_practice[m]
] for m in self.uncertainty_massif_names]
@@ -219,3 +217,15 @@ class StudyVisualizerForNonStationaryTrends(StudyVisualizer):
def model_name_to_uncertainty_method_to_ratio_above_eurocode(self):
assert self.uncertainty_massif_names == self.study.study_massif_names
+
+ # Some checks with Gelman convergence diagnosis
+
+ def massif_name_to_gelman_convergence_value(self, mcmc_iterations, model_class, nb_chains):
+ arguments = [(self.massif_name_to_non_null_years_and_maxima[m], mcmc_iterations, model_class, nb_chains)
+ for m in self.uncertainty_massif_names]
+ if self.multiprocessing:
+ with Pool(NB_CORES) as p:
+ res = p.starmap(compute_gelman_convergence_value, arguments)
+ else:
+ res = [compute_gelman_convergence_value(*argument) for argument in arguments]
+ return dict(zip(self.uncertainty_massif_names, res))
diff --git a/experiment/trend_analysis/univariate_test/utils.py b/experiment/trend_analysis/univariate_test/utils.py
index 75e4f7ff..cd37c902 100644
--- a/experiment/trend_analysis/univariate_test/utils.py
+++ b/experiment/trend_analysis/univariate_test/utils.py
@@ -20,8 +20,8 @@ def load_temporal_coordinates_and_dataset(maxima, years):
return coordinates, dataset
-def fitted_linear_margin_estimator(model_class, coordinates, dataset, starting_year, fit_method):
- model = model_class(coordinates, starting_point=starting_year, fit_method=fit_method)
+def fitted_linear_margin_estimator(model_class, coordinates, dataset, starting_year, fit_method, **model_kwargs):
+ model = model_class(coordinates, starting_point=starting_year, fit_method=fit_method, **model_kwargs)
estimator = LinearMarginEstimator(dataset, model)
estimator.fit()
return estimator
diff --git a/extreme_fit/model/margin_model/linear_margin_model/abstract_temporal_linear_margin_model.py b/extreme_fit/model/margin_model/linear_margin_model/abstract_temporal_linear_margin_model.py
index 851b9f54..161033e9 100644
--- a/extreme_fit/model/margin_model/linear_margin_model/abstract_temporal_linear_margin_model.py
+++ b/extreme_fit/model/margin_model/linear_margin_model/abstract_temporal_linear_margin_model.py
@@ -24,8 +24,10 @@ class AbstractTemporalLinearMarginModel(LinearMarginModel):
"""Linearity only with respect to the temporal coordinates"""
def __init__(self, coordinates: AbstractCoordinates, use_start_value=False, params_start_fit=None,
- params_sample=None, starting_point=None, fit_method=TemporalMarginFitMethod.is_mev_gev_fit):
+ params_sample=None, starting_point=None, fit_method=TemporalMarginFitMethod.is_mev_gev_fit,
+ nb_iterations_for_bayesian_fit=5000):
super().__init__(coordinates, use_start_value, params_start_fit, params_sample, starting_point)
+ self.nb_iterations_for_bayesian_fit = nb_iterations_for_bayesian_fit
assert isinstance(fit_method, TemporalMarginFitMethod)
self.fit_method = fit_method
@@ -73,7 +75,7 @@ class AbstractTemporalLinearMarginModel(LinearMarginModel):
method='Bayesian',
priorFun="fevdPriorCustom",
priorParams=r.list(q=r.c(6), p=r.c(9)),
- iter=5000,
+ iter=self.nb_iterations_for_bayesian_fit,
**r_type_argument_kwargs
)
return ResultFromBayesianExtremes(res, self.margin_function_start_fit.gev_param_name_to_dims)
--
GitLab