From 823d5c8bcdbd148d35af8f8fd36c858bc5fe5488 Mon Sep 17 00:00:00 2001
From: Le Roux Erwan <erwan.le-roux@irstea.fr>
Date: Mon, 25 Nov 2019 11:26:50 +0100
Subject: [PATCH] [PAPER 1] set starting year to None to avoid bug due to
missing data for the year 1958 for instance
---
.../eurocode_visualizer.py | 5 +++--
.../main_result_trends_and_return_levels.py | 17 ++++++++++++-----
...tudy_visualizer_for_non_stationary_trends.py | 9 ++++++++-
.../eurocode_return_level_uncertainties.py | 2 +-
4 files changed, 24 insertions(+), 9 deletions(-)
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 6b887c9c..78bd644d 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
@@ -37,8 +37,8 @@ def plot_uncertainty_massifs(altitude_to_visualizer: Dict[int, StudyVisualizerFo
:return:
"""
visualizer = list(altitude_to_visualizer.values())[0]
- # Subdivide massif names in group of 5
- m = 4
+ # Subdivide massif names in group of 3
+ m = 3
uncertainty_massif_names = visualizer.uncertainty_massif_names
n = (len(uncertainty_massif_names) // m) + 1
for i in list(range(n))[:]:
@@ -140,6 +140,7 @@ def plot_single_uncertainty_massif_and_non_stationary_context(ax, massif_name, n
title = '{} ({} Eurocodes area) with a {} model'.format(massif_name_str, eurocode_region_str,
non_stationary_context)
ax.set_title(title)
+ ax.set_xticks(altitudes)
ax.set_ylabel('50-year return level of SL (kN $m^-2$)')
ax.set_xlabel('Altitude (m)')
ax.grid()
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 172b0396..8554af5c 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
@@ -67,9 +67,16 @@ def major_result():
if __name__ == '__main__':
# minor_result(altitude=1800)
- intermediate_result(altitudes=[1500, 1800], massif_names=['Chartreuse'],
- uncertainty_methods=[ConfidenceIntervalMethodFromExtremes.ci_mle],
- non_stationary_uncertainty=[False])
+ # intermediate_result(altitudes=[1500, 1800], massif_names=['Chartreuse'],
+ # uncertainty_methods=[ConfidenceIntervalMethodFromExtremes.ci_mle],
+ # non_stationary_uncertainty=[False])
+ # intermediate_result(altitudes=[1500, 1800], massif_names=None,
+ # uncertainty_methods=[ConfidenceIntervalMethodFromExtremes.ci_mle],
+ # non_stationary_uncertainty=[False])
+ # intermediate_result(altitudes=[300, 600, 900, 1200, 1500, 1800], massif_names=None,
+ # uncertainty_methods=[ConfidenceIntervalMethodFromExtremes.ci_mle],
+ # non_stationary_uncertainty=[False])
intermediate_result(altitudes=[300, 600, 900, 1200, 1500, 1800], massif_names=None,
- uncertainty_methods=[ConfidenceIntervalMethodFromExtremes.ci_mle],
- non_stationary_uncertainty=[False])
+ uncertainty_methods=[ConfidenceIntervalMethodFromExtremes.ci_mle,
+ ConfidenceIntervalMethodFromExtremes.ci_bayes],
+ non_stationary_uncertainty=[False, True])
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/result_trends_and_return_levels/study_visualizer_for_non_stationary_trends.py
index e4f21054..b80f46ac 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/result_trends_and_return_levels/study_visualizer_for_non_stationary_trends.py
@@ -6,6 +6,7 @@ from typing import Dict, Tuple
import numpy as np
from cached_property import cached_property
+from experiment.eurocode_data.utils import EUROCODE_QUANTILE
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
@@ -70,6 +71,12 @@ class StudyVisualizerForNonStationaryTrends(StudyVisualizer):
return {m: np.count_nonzero(maxima) / len(maxima) for m, (_, maxima) in
self.massif_name_to_years_and_maxima.items()}
+ @cached_property
+ def massif_name_to_eurocode_quantile_in_practice(self):
+ """Due to missing data, the the eurocode quantile which 0.98 if we have all the data
+ correspond in practice to the quantile psnow x 0.98 of the data where there is snow"""
+ return {m: p * EUROCODE_QUANTILE for m, p in self.massif_name_to_psnow.items()}
+
@cached_property
def massif_name_to_non_null_years_and_maxima(self):
d = {}
@@ -80,7 +87,7 @@ class StudyVisualizerForNonStationaryTrends(StudyVisualizer):
@cached_property
def massif_name_to_minimized_aic_non_stationary_trend_test(self) -> Dict[str, AbstractGevTrendTest]:
- starting_year = 1958
+ starting_year = None
massif_name_to_trend_test_that_minimized_aic = {}
for massif_name, (x, y) in self.massif_name_to_non_null_years_and_maxima.items():
non_stationary_trend_test = [t(x, y, starting_year) for t in self.non_stationary_trend_test]
diff --git a/extreme_fit/model/result_from_model_fit/result_from_extremes/eurocode_return_level_uncertainties.py b/extreme_fit/model/result_from_model_fit/result_from_extremes/eurocode_return_level_uncertainties.py
index 99599dca..71d7e183 100644
--- a/extreme_fit/model/result_from_model_fit/result_from_extremes/eurocode_return_level_uncertainties.py
+++ b/extreme_fit/model/result_from_model_fit/result_from_extremes/eurocode_return_level_uncertainties.py
@@ -49,7 +49,7 @@ class EurocodeConfidenceIntervalFromExtremes(object):
else:
fit_method = TemporalMarginFitMethod.extremes_fevd_mle
# Fitted estimator
- fitted_estimator = fitted_linear_margin_estimator(model_class, coordinates, dataset, starting_year=1958,
+ fitted_estimator = fitted_linear_margin_estimator(model_class, coordinates, dataset, starting_year=None,
fit_method=fit_method)
# Load object from result from extremes
return cls.from_estimator_extremes(fitted_estimator, ci_method, temporal_covariate)
--
GitLab