From f8fdcd0c11e39cda82832a12e17ed7d0891b604f Mon Sep 17 00:00:00 2001
From: Le Roux Erwan <erwan.le-roux@irstea.fr>
Date: Thu, 21 Feb 2019 15:46:35 +0100
Subject: [PATCH] [EXTREME PARAMETERS] create an attribute undefined_parameters
in ExtremeParameters object, to avoid crash when a shape is < 0 for point
while we want to display
---
experiment/meteo_france_SCM_study/main_visualize.py | 11 ++++++-----
.../safran/safran_visualizer.py | 4 ++++
.../margin_function/abstract_margin_function.py | 1 +
extreme_estimator/margin_fits/extreme_params.py | 7 ++++++-
extreme_estimator/margin_fits/gev/gev_params.py | 11 +++++++++--
.../margin_fits/plot/create_shifted_cmap.py | 7 ++++---
.../test_margin_fits/test_gev/test_gev_params.py | 6 ++++++
7 files changed, 36 insertions(+), 11 deletions(-)
diff --git a/experiment/meteo_france_SCM_study/main_visualize.py b/experiment/meteo_france_SCM_study/main_visualize.py
index 9110ac96..0a96d042 100644
--- a/experiment/meteo_france_SCM_study/main_visualize.py
+++ b/experiment/meteo_france_SCM_study/main_visualize.py
@@ -19,7 +19,7 @@ def study_iterator(study_class, only_first_one=False, both_altitude=False, verbo
if verbose:
print('Loading studies....')
for nb_day in nb_days:
- for alti in AbstractStudy.ALTITUDES[::1]:
+ for alti in AbstractStudy.ALTITUDES[::-1]:
if verbose:
print('alti: {}, nb_day: {}'.format(alti, nb_day))
study = study_class(alti, nb_day) if is_safran_study else study_class(alti)
@@ -36,11 +36,12 @@ def extended_visualization():
for study_class in SCM_EXTENDED_STUDIES[1:2]:
for study in study_iterator(study_class, only_first_one=True):
study_visualizer = StudyVisualizer(study)
- study_visualizer.visualize_all_kde_graphs()
+ # study_visualizer.visualize_all_kde_graphs()
+ study_visualizer.visualize_experimental_law()
def normal_visualization():
- for study_class in SCM_STUDIES[:1]:
+ for study_class in SCM_STUDIES[1:2]:
for study in study_iterator(study_class, only_first_one=True):
study_visualizer = StudyVisualizer(study)
# study_visualizer.visualize_independent_margin_fits(threshold=[None, 20, 40, 60][0])
@@ -62,6 +63,6 @@ def complete_analysis(only_first_one=False):
if __name__ == '__main__':
- # normal_visualization()
- extended_visualization()
+ normal_visualization()
+ # extended_visualization()
# complete_analysis()
diff --git a/experiment/meteo_france_SCM_study/safran/safran_visualizer.py b/experiment/meteo_france_SCM_study/safran/safran_visualizer.py
index 11355875..b147faa9 100644
--- a/experiment/meteo_france_SCM_study/safran/safran_visualizer.py
+++ b/experiment/meteo_france_SCM_study/safran/safran_visualizer.py
@@ -43,6 +43,10 @@ class StudyVisualizer(object):
def dataset(self):
return AbstractDataset(self.observations, self.coordinates)
+ def visualize_experimental_law(self):
+ plot_name = ' experimental law'
+ self.show_or_save_to_file(plot_name)
+
def visualize_all_kde_graphs(self):
massif_names = self.study.safran_massif_names
nb_columns = 5
diff --git a/extreme_estimator/extreme_models/margin_model/margin_function/abstract_margin_function.py b/extreme_estimator/extreme_models/margin_model/margin_function/abstract_margin_function.py
index 09866a6a..cc1ac24b 100644
--- a/extreme_estimator/extreme_models/margin_model/margin_function/abstract_margin_function.py
+++ b/extreme_estimator/extreme_models/margin_model/margin_function/abstract_margin_function.py
@@ -121,6 +121,7 @@ class AbstractMarginFunction(object):
grid = []
for i, xi in enumerate(linspace):
gev_param = self.get_gev_params(np.array([xi]))
+ assert not gev_param.has_undefined_parameters, 'This case needs to be handled during display'
grid.append(gev_param.summary_dict)
grid = {gev_param: [g[gev_param] for g in grid] for gev_param in GevParams.SUMMARY_NAMES}
return grid, linspace
diff --git a/extreme_estimator/margin_fits/extreme_params.py b/extreme_estimator/margin_fits/extreme_params.py
index 01a39d3f..a20ace4f 100644
--- a/extreme_estimator/margin_fits/extreme_params.py
+++ b/extreme_estimator/margin_fits/extreme_params.py
@@ -1,4 +1,5 @@
from abc import ABC
+import numpy as np
from extreme_estimator.margin_fits.abstract_params import AbstractParams
@@ -13,4 +14,8 @@ class ExtremeParams(AbstractParams, ABC):
self.location = loc
self.scale = scale
self.shape = shape
- assert self.scale > 0
+ # Scale cannot be negative
+ # (sometimes it happens, when we want to find a quantile for every point of a 2D map
+ # then it can happen that a corner point that was not used for fitting correspond to a negative scale,
+ # in the case we set all the parameters as equal to np.nan, and we will not display those points)
+ self.has_undefined_parameters = self.scale <= 0
\ No newline at end of file
diff --git a/extreme_estimator/margin_fits/gev/gev_params.py b/extreme_estimator/margin_fits/gev/gev_params.py
index b4d9cbc2..17081517 100644
--- a/extreme_estimator/margin_fits/gev/gev_params.py
+++ b/extreme_estimator/margin_fits/gev/gev_params.py
@@ -1,5 +1,6 @@
from extreme_estimator.extreme_models.utils import r
from extreme_estimator.margin_fits.extreme_params import ExtremeParams
+import numpy as np
class GevParams(ExtremeParams):
@@ -14,8 +15,14 @@ class GevParams(ExtremeParams):
self.block_size = block_size
def quantile(self, p) -> float:
- return r.qgev(p, self.location, self.scale, self.shape)[0]
+ if self.has_undefined_parameters:
+ return np.nan
+ else:
+ return r.qgev(p, self.location, self.scale, self.shape)[0]
@property
def param_values(self):
- return [self.location, self.scale, self.shape]
\ No newline at end of file
+ if self.has_undefined_parameters:
+ return [np.nan for _ in range(3)]
+ else:
+ return [self.location, self.scale, self.shape]
\ No newline at end of file
diff --git a/extreme_estimator/margin_fits/plot/create_shifted_cmap.py b/extreme_estimator/margin_fits/plot/create_shifted_cmap.py
index b0fe67fc..db147f9c 100644
--- a/extreme_estimator/margin_fits/plot/create_shifted_cmap.py
+++ b/extreme_estimator/margin_fits/plot/create_shifted_cmap.py
@@ -24,6 +24,8 @@ def plot_extreme_param(ax, gev_param_name, values):
midpoint = 1.0
elif vmin > 0 and vmax > 0:
midpoint = 0.0
+ else:
+ raise ValueError('Unexpected values: vmin={}, vmax={}'.format(vmin, vmax))
cmap = [plt.cm.coolwarm, plt.cm.bwr, plt.cm.seismic][1]
shifted_cmap = shiftedColorMap(cmap, midpoint=midpoint, name='shifted')
norm = mpl.colors.Normalize(vmin=vmin, vmax=vmax)
@@ -47,10 +49,9 @@ def get_color_rbga_shifted(ax, gev_param_name, values):
def imshow_shifted(ax, gev_param_name, values, x, y):
- norm, shifted_cmap = plot_extreme_param(ax, gev_param_name, values)
+ masked_array = np.ma.masked_where(np.isnan(values), values)
+ norm, shifted_cmap = plot_extreme_param(ax, gev_param_name, masked_array)
shifted_cmap.set_bad(color='white')
-
- masked_array = values
if gev_param_name != ExtremeParams.SHAPE:
epsilon = 1e-2 * (np.max(values) - np.min(values))
value = np.min(values)
diff --git a/test/test_extreme_estimator/test_margin_fits/test_gev/test_gev_params.py b/test/test_extreme_estimator/test_margin_fits/test_gev/test_gev_params.py
index dab5ca55..39d3da70 100644
--- a/test/test_extreme_estimator/test_margin_fits/test_gev/test_gev_params.py
+++ b/test/test_extreme_estimator/test_margin_fits/test_gev/test_gev_params.py
@@ -14,6 +14,12 @@ class TestGevParams(unittest.TestCase):
for quantile_name, p in gev_params.quantile_name_to_p.items():
self.assertAlmostEqual(- 1 / np.log(p), quantile_dict[quantile_name])
+ def test_negative_scale(self):
+ gev_params = GevParams(loc=1.0, shape=1.0, scale=-1.0)
+ for p in [0.1, 0.5, 0.9]:
+ q = gev_params.quantile(p)
+ self.assertTrue(np.isnan(q))
+
if __name__ == '__main__':
unittest.main()
--
GitLab