From 4a01d9573b27546c8a95abf4c8c67844a4c84aa9 Mon Sep 17 00:00:00 2001
From: Le Roux Erwan <erwan.le-roux@irstea.fr>
Date: Tue, 30 Jun 2020 15:34:48 +0200
Subject: [PATCH] [contrasting] add Gumbel altitudinal models. But so far only
the first three models are passing tests (models with quadratic & cross terms
are crashing)
---
...al_models.py => gev_altitudinal_models.py} | 20 ++++-
.../gumbel_altitudinal_models.py | 89 +++++++++++++++++++
.../polynomial_margin_model/utils.py | 54 ++++++++---
.../altitudes_fit/main_altitudes_studies.py | 14 +--
.../one_fold_analysis/one_fold_fit.py | 15 ++--
.../test_gev_spatio_temporal_extremes_mle.py | 15 ++--
6 files changed, 176 insertions(+), 31 deletions(-)
rename extreme_fit/model/margin_model/polynomial_margin_model/{altitudinal_models.py => gev_altitudinal_models.py} (87%)
create mode 100644 extreme_fit/model/margin_model/polynomial_margin_model/gumbel_altitudinal_models.py
diff --git a/extreme_fit/model/margin_model/polynomial_margin_model/altitudinal_models.py b/extreme_fit/model/margin_model/polynomial_margin_model/gev_altitudinal_models.py
similarity index 87%
rename from extreme_fit/model/margin_model/polynomial_margin_model/altitudinal_models.py
rename to extreme_fit/model/margin_model/polynomial_margin_model/gev_altitudinal_models.py
index fceadfce..5b67925c 100644
--- a/extreme_fit/model/margin_model/polynomial_margin_model/altitudinal_models.py
+++ b/extreme_fit/model/margin_model/polynomial_margin_model/gev_altitudinal_models.py
@@ -5,6 +5,7 @@ from extreme_fit.model.margin_model.polynomial_margin_model.spatio_temporal_poly
class AbstractAltitudinalModel(AbstractSpatioTemporalPolynomialModel):
+ DISTRIBUTION_STR = 'Gev'
def load_margin_function(self, param_name_to_dims=None):
return super().load_margin_function(self.param_name_to_list_dim_and_degree_for_margin_function)
@@ -22,7 +23,7 @@ class AbstractAltitudinalModel(AbstractSpatioTemporalPolynomialModel):
def dim_to_str_number(self, param_name, dim):
list_dim_and_degree = self.param_name_to_list_dim_and_degree[param_name]
- dims = [d for d, _ in list_dim_and_degree]
+ dims = [d for d, _ in list_dim_and_degree]
if dim not in dims:
return '0'
else:
@@ -31,7 +32,7 @@ class AbstractAltitudinalModel(AbstractSpatioTemporalPolynomialModel):
@property
def name_str(self):
- name = 'Gev'
+ name = self.DISTRIBUTION_STR
name += self.dim_to_str_number(GevParams.LOC, self.coordinates.idx_temporal_coordinates)
name += self.dim_to_str_number(GevParams.SCALE, self.coordinates.idx_temporal_coordinates)
if isinstance(self, AbstractAddCrossTermForLocation):
@@ -49,6 +50,16 @@ class StationaryAltitudinal(AbstractAltitudinalModel):
}
+class NonStationaryAltitudinalScaleLinear(AbstractAltitudinalModel):
+
+ @property
+ def param_name_to_list_dim_and_degree(self):
+ return {
+ GevParams.LOC: [(self.coordinates.idx_x_coordinates, 1)],
+ GevParams.SCALE: [(self.coordinates.idx_x_coordinates, 1), (self.coordinates.idx_temporal_coordinates, 1)]
+ }
+
+
class NonStationaryAltitudinalLocationLinear(AbstractAltitudinalModel):
@property
@@ -124,3 +135,8 @@ class NonStationaryAltitudinalLocationLinearScaleLinearCrossTermForLocation(Abst
class NonStationaryAltitudinalLocationQuadraticScaleLinearCrossTermForLocation(AbstractAddCrossTermForLocation,
NonStationaryAltitudinalLocationQuadraticScaleLinear):
pass
+
+
+class NonStationaryAltitudinalScaleLinearCrossTermForLocation(AbstractAddCrossTermForLocation,
+ NonStationaryAltitudinalScaleLinear):
+ pass
diff --git a/extreme_fit/model/margin_model/polynomial_margin_model/gumbel_altitudinal_models.py b/extreme_fit/model/margin_model/polynomial_margin_model/gumbel_altitudinal_models.py
new file mode 100644
index 00000000..35da1544
--- /dev/null
+++ b/extreme_fit/model/margin_model/polynomial_margin_model/gumbel_altitudinal_models.py
@@ -0,0 +1,89 @@
+from extreme_fit.distribution.gev.gev_params import GevParams
+from extreme_fit.model.margin_model.linear_margin_model.temporal_linear_margin_models import GumbelTemporalModel
+from extreme_fit.model.margin_model.polynomial_margin_model.gev_altitudinal_models import AbstractAltitudinalModel, \
+ StationaryAltitudinal, NonStationaryAltitudinalScaleLinear, NonStationaryAltitudinalLocationLinear, \
+ NonStationaryAltitudinalLocationQuadratic, NonStationaryAltitudinalLocationLinearScaleLinear, \
+ NonStationaryAltitudinalLocationQuadraticScaleLinear, NonStationaryAltitudinalScaleLinearCrossTermForLocation, \
+ NonStationaryCrossTermForLocation, NonStationaryAltitudinalLocationLinearCrossTermForLocation, \
+ NonStationaryAltitudinalLocationQuadraticCrossTermForLocation, \
+ NonStationaryAltitudinalLocationLinearScaleLinearCrossTermForLocation, \
+ NonStationaryAltitudinalLocationQuadraticScaleLinearCrossTermForLocation, AbstractAddCrossTermForLocation
+from extreme_fit.model.margin_model.polynomial_margin_model.polynomial_margin_model import PolynomialMarginModel
+from extreme_fit.model.margin_model.polynomial_margin_model.spatio_temporal_polynomial_model import \
+ AbstractSpatioTemporalPolynomialModel
+from extreme_fit.model.margin_model.utils import MarginFitMethod
+from spatio_temporal_dataset.coordinates.abstract_coordinates import AbstractCoordinates
+
+
+class AbstractGumbelAltitudinalModel(AbstractAltitudinalModel):
+ DISTRIBUTION_STR = 'Gum'
+
+ def __init__(self, coordinates: AbstractCoordinates, params_user=None, starting_point=None,
+ fit_method=MarginFitMethod.extremes_fevd_mle, nb_iterations_for_bayesian_fit=5000,
+ params_initial_fit_bayesian=None, params_class=GevParams, max_degree=2):
+ super().__init__(coordinates, params_user, starting_point, fit_method, nb_iterations_for_bayesian_fit,
+ params_initial_fit_bayesian, "Gumbel", params_class, max_degree)
+
+ @property
+ def nb_params(self):
+ return super().nb_params - 1
+
+
+class StationaryGumbelAltitudinal(AbstractGumbelAltitudinalModel, StationaryAltitudinal):
+ pass
+
+
+class NonStationaryGumbelAltitudinalScaleLinear(AbstractGumbelAltitudinalModel, NonStationaryAltitudinalScaleLinear):
+ pass
+
+
+class NonStationaryGumbelAltitudinalLocationLinear(AbstractGumbelAltitudinalModel,
+ NonStationaryAltitudinalLocationLinear):
+ pass
+
+
+class NonStationaryGumbelAltitudinalLocationQuadratic(AbstractGumbelAltitudinalModel,
+ NonStationaryAltitudinalLocationQuadratic):
+ pass
+
+
+class NonStationaryGumbelAltitudinalLocationLinearScaleLinear(AbstractGumbelAltitudinalModel,
+ NonStationaryAltitudinalLocationLinearScaleLinear):
+ pass
+
+
+class NonStationaryGumbelAltitudinalLocationQuadraticScaleLinear(AbstractGumbelAltitudinalModel,
+ NonStationaryAltitudinalLocationQuadraticScaleLinear):
+ pass
+
+
+# Add cross terms
+
+
+class NonStationaryGumbelCrossTermForLocation(AbstractAddCrossTermForLocation, AbstractGumbelAltitudinalModel, StationaryAltitudinal):
+ pass
+
+
+class NonStationaryGumbelAltitudinalLocationLinearCrossTermForLocation(AbstractAddCrossTermForLocation, AbstractGumbelAltitudinalModel,
+ NonStationaryAltitudinalLocationLinear):
+ pass
+
+
+class NonStationaryGumbelAltitudinalLocationQuadraticCrossTermForLocation(AbstractGumbelAltitudinalModel,
+ NonStationaryAltitudinalLocationQuadraticCrossTermForLocation):
+ pass
+
+
+class NonStationaryGumbelAltitudinalLocationLinearScaleLinearCrossTermForLocation(AbstractGumbelAltitudinalModel,
+ NonStationaryAltitudinalLocationLinearScaleLinearCrossTermForLocation):
+ pass
+
+
+class NonStationaryGumbelAltitudinalLocationQuadraticScaleLinearCrossTermForLocation(AbstractGumbelAltitudinalModel,
+ NonStationaryAltitudinalLocationQuadraticScaleLinearCrossTermForLocation):
+ pass
+
+
+class NonStationaryGumbelAltitudinalScaleLinearCrossTermForLocation(AbstractGumbelAltitudinalModel,
+ NonStationaryAltitudinalScaleLinearCrossTermForLocation):
+ pass
diff --git a/extreme_fit/model/margin_model/polynomial_margin_model/utils.py b/extreme_fit/model/margin_model/polynomial_margin_model/utils.py
index 615edc68..0ff5aaad 100644
--- a/extreme_fit/model/margin_model/polynomial_margin_model/utils.py
+++ b/extreme_fit/model/margin_model/polynomial_margin_model/utils.py
@@ -1,32 +1,58 @@
-from extreme_fit.model.margin_model.polynomial_margin_model.altitudinal_models import StationaryAltitudinal, \
+from extreme_fit.model.margin_model.polynomial_margin_model.gev_altitudinal_models import StationaryAltitudinal, \
NonStationaryAltitudinalLocationLinear, NonStationaryAltitudinalLocationQuadratic, \
NonStationaryAltitudinalLocationLinearScaleLinear, NonStationaryAltitudinalLocationQuadraticScaleLinear, \
NonStationaryCrossTermForLocation, NonStationaryAltitudinalLocationLinearCrossTermForLocation, \
NonStationaryAltitudinalLocationQuadraticCrossTermForLocation, \
NonStationaryAltitudinalLocationLinearScaleLinearCrossTermForLocation, \
- NonStationaryAltitudinalLocationQuadraticScaleLinearCrossTermForLocation
+ NonStationaryAltitudinalLocationQuadraticScaleLinearCrossTermForLocation, NonStationaryAltitudinalScaleLinear, \
+ NonStationaryAltitudinalScaleLinearCrossTermForLocation
+from extreme_fit.model.margin_model.polynomial_margin_model.gumbel_altitudinal_models import \
+ StationaryGumbelAltitudinal, NonStationaryGumbelAltitudinalScaleLinear, \
+ NonStationaryGumbelAltitudinalLocationLinear, NonStationaryGumbelAltitudinalLocationQuadratic, \
+ NonStationaryGumbelAltitudinalLocationLinearScaleLinear, NonStationaryGumbelAltitudinalLocationQuadraticScaleLinear, \
+ NonStationaryGumbelAltitudinalLocationLinearCrossTermForLocation, \
+ NonStationaryGumbelAltitudinalLocationQuadraticCrossTermForLocation, \
+ NonStationaryGumbelAltitudinalLocationLinearScaleLinearCrossTermForLocation, \
+ NonStationaryGumbelAltitudinalLocationQuadraticScaleLinearCrossTermForLocation, \
+ NonStationaryGumbelAltitudinalScaleLinearCrossTermForLocation, NonStationaryGumbelCrossTermForLocation
from extreme_fit.model.margin_model.polynomial_margin_model.spatio_temporal_polynomial_model import \
NonStationaryLocationSpatioTemporalLinearityModel1, NonStationaryLocationSpatioTemporalLinearityModel2, \
NonStationaryLocationSpatioTemporalLinearityModel3, NonStationaryLocationSpatioTemporalLinearityModel4, \
NonStationaryLocationSpatioTemporalLinearityModel5, NonStationaryLocationSpatioTemporalLinearityModelAssertError1, \
NonStationaryLocationSpatioTemporalLinearityModelAssertError2, \
NonStationaryLocationSpatioTemporalLinearityModelAssertError3, NonStationaryLocationSpatioTemporalLinearityModel6
-ALTITUDINAL_MODELS = [
- StationaryAltitudinal,
- NonStationaryAltitudinalLocationLinear,
- NonStationaryAltitudinalLocationQuadratic,
- NonStationaryAltitudinalLocationLinearScaleLinear,
- NonStationaryAltitudinalLocationQuadraticScaleLinear,
- NonStationaryCrossTermForLocation,
- NonStationaryAltitudinalLocationLinearCrossTermForLocation,
- NonStationaryAltitudinalLocationQuadraticCrossTermForLocation,
- NonStationaryAltitudinalLocationLinearScaleLinearCrossTermForLocation,
- NonStationaryAltitudinalLocationQuadraticScaleLinearCrossTermForLocation,
- ][:]
+ALTITUDINAL_GEV_MODELS = [
+ StationaryAltitudinal,
+ NonStationaryAltitudinalScaleLinear,
+ NonStationaryAltitudinalLocationLinear,
+ NonStationaryAltitudinalLocationQuadratic,
+ NonStationaryAltitudinalLocationLinearScaleLinear,
+ NonStationaryAltitudinalLocationQuadraticScaleLinear,
+ NonStationaryCrossTermForLocation,
+ NonStationaryAltitudinalScaleLinearCrossTermForLocation,
+ NonStationaryAltitudinalLocationLinearCrossTermForLocation,
+ NonStationaryAltitudinalLocationQuadraticCrossTermForLocation,
+ NonStationaryAltitudinalLocationLinearScaleLinearCrossTermForLocation,
+ NonStationaryAltitudinalLocationQuadraticScaleLinearCrossTermForLocation,
+ ][:]
+ALTITUDINAL_GUMBEL_MODELS = [
+ StationaryGumbelAltitudinal,
+ NonStationaryGumbelAltitudinalScaleLinear,
+ NonStationaryGumbelAltitudinalLocationLinear,
+ NonStationaryGumbelAltitudinalLocationQuadratic,
+ NonStationaryGumbelAltitudinalLocationLinearScaleLinear,
+ NonStationaryGumbelAltitudinalLocationQuadraticScaleLinear,
+ NonStationaryGumbelCrossTermForLocation,
+ NonStationaryGumbelAltitudinalLocationLinearCrossTermForLocation,
+ NonStationaryGumbelAltitudinalLocationQuadraticCrossTermForLocation,
+ NonStationaryGumbelAltitudinalLocationLinearScaleLinearCrossTermForLocation,
+ NonStationaryGumbelAltitudinalLocationQuadraticScaleLinearCrossTermForLocation,
+ NonStationaryGumbelAltitudinalScaleLinearCrossTermForLocation,
+ ][:]
MODELS_THAT_SHOULD_RAISE_AN_ASSERTION_ERROR = [NonStationaryLocationSpatioTemporalLinearityModelAssertError1,
NonStationaryLocationSpatioTemporalLinearityModelAssertError2,
diff --git a/projects/altitude_spatial_model/altitudes_fit/main_altitudes_studies.py b/projects/altitude_spatial_model/altitudes_fit/main_altitudes_studies.py
index 99e903a5..576b2ca2 100644
--- a/projects/altitude_spatial_model/altitudes_fit/main_altitudes_studies.py
+++ b/projects/altitude_spatial_model/altitudes_fit/main_altitudes_studies.py
@@ -10,19 +10,22 @@ from extreme_data.meteo_france_data.scm_models_data.safran.safran import SafranS
SafranSnowfall5Days, SafranSnowfall7Days, SafranPrecipitation1Day, SafranPrecipitation3Days, \
SafranPrecipitation5Days, SafranPrecipitation7Days
from extreme_data.meteo_france_data.scm_models_data.utils import Season
-from extreme_fit.model.margin_model.polynomial_margin_model.utils import ALTITUDINAL_MODELS
+from extreme_fit.model.margin_model.polynomial_margin_model.utils import ALTITUDINAL_GEV_MODELS, \
+ ALTITUDINAL_GUMBEL_MODELS
from projects.altitude_spatial_model.altitudes_fit.altitudes_studies import AltitudesStudies
from projects.altitude_spatial_model.altitudes_fit.one_fold_analysis.altitudes_studies_visualizer_for_non_stationary_models import \
AltitudesStudiesVisualizerForNonStationaryModels
def plot_altitudinal_fit(studies, massif_names=None):
- model_classes = ALTITUDINAL_MODELS
+ model_classes = ALTITUDINAL_GEV_MODELS + ALTITUDINAL_GUMBEL_MODELS
+ # model_classes = ALTITUDINAL_GUMBEL_MODELS
visualizer = AltitudesStudiesVisualizerForNonStationaryModels(studies=studies,
model_classes=model_classes,
massif_names=massif_names,
show=False)
# Plot the results for the model that minimizes the individual aic
+ OneFoldFit.best_estimator_minimizes_mean_aic = False
visualizer.plot_moments()
visualizer.plot_shape_map()
@@ -31,13 +34,14 @@ def plot_altitudinal_fit(studies, massif_names=None):
model_class_to_aic_scores = {model_class: [] for model_class in model_classes}
for one_fold_fit in visualizer.massif_name_to_one_fold_fit.values():
for model_class, estimator in one_fold_fit.model_class_to_estimator_with_finite_aic.items():
- print(estimator.n())
aic_score_normalized = estimator.aic() / estimator.n()
model_class_to_aic_scores[model_class].append(aic_score_normalized)
model_class_to_mean_aic_score = {model_class: np.array(aic_scores).mean()
for model_class, aic_scores in model_class_to_aic_scores.items()}
+ print(model_class_to_mean_aic_score)
sorted_model_class = sorted(model_classes, key=lambda m: model_class_to_mean_aic_score[m])
best_model_class_for_mean_aic = sorted_model_class[0]
+ print(best_model_class_for_mean_aic)
for one_fold_fit in visualizer.massif_name_to_one_fold_fit.values():
one_fold_fit.best_estimator_class_for_mean_aic = best_model_class_for_mean_aic
@@ -58,8 +62,8 @@ def plot_moments(studies, massif_names=None):
def main():
altitudes = [900, 1200, 1500, 1800, 2100, 2400, 2700, 3000][4:7]
- altitudes = [900, 1200, 1500, 1800, 2100, 2400, 2700, 3000][:]
- altitudes = [900, 1200, 1500, 1800, 2100, 2400, 2700, 3000, 3300, 3600, 3900]
+ # altitudes = [900, 1200, 1500, 1800, 2100, 2400, 2700, 3000][:]
+ # altitudes = [900, 1200, 1500, 1800, 2100, 2400, 2700, 3000, 3300, 3600, 3900]
study_classes = [SafranSnowfall1Day, SafranSnowfall3Days, SafranSnowfall5Days, SafranSnowfall7Days][:2]
study_classes = [SafranPrecipitation1Day, SafranPrecipitation3Days, SafranPrecipitation5Days,
SafranPrecipitation7Days][:]
diff --git a/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/one_fold_fit.py b/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/one_fold_fit.py
index 113e1329..cc348b81 100644
--- a/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/one_fold_fit.py
+++ b/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/one_fold_fit.py
@@ -5,7 +5,9 @@ from cached_property import cached_property
from scipy.stats import chi2
from extreme_fit.estimator.margin_estimator.utils import fitted_linear_margin_estimator_short
-from extreme_fit.model.margin_model.polynomial_margin_model.altitudinal_models import StationaryAltitudinal
+from extreme_fit.model.margin_model.polynomial_margin_model.gev_altitudinal_models import StationaryAltitudinal
+from extreme_fit.model.margin_model.polynomial_margin_model.gumbel_altitudinal_models import \
+ StationaryGumbelAltitudinal, AbstractGumbelAltitudinalModel
from extreme_fit.model.margin_model.utils import MarginFitMethod
from root_utils import classproperty
@@ -137,9 +139,10 @@ class OneFoldFit(object):
@property
def stationary_estimator(self):
- for estimator in self.sorted_estimators_with_finite_aic:
- if isinstance(estimator.margin_model, StationaryAltitudinal):
- return estimator
+ if isinstance(self.best_estimator.margin_model, AbstractGumbelAltitudinalModel):
+ return self.model_class_to_estimator_with_finite_aic[StationaryGumbelAltitudinal]
+ else:
+ return self.model_class_to_estimator_with_finite_aic[StationaryAltitudinal]
@property
def likelihood_ratio(self):
@@ -151,7 +154,9 @@ class OneFoldFit(object):
@property
def is_significant(self) -> bool:
- if isinstance(self.best_estimator.margin_model, StationaryAltitudinal):
+ stationary_model_classes = [StationaryAltitudinal, StationaryGumbelAltitudinal]
+ if any([isinstance(self.best_estimator.margin_model, c)
+ for c in stationary_model_classes]):
return False
else:
return self.likelihood_ratio > chi2.ppf(q=1 - self.SIGNIFICANCE_LEVEL, df=self.degree_freedom_chi2)
diff --git a/test/test_extreme_fit/test_estimator/test_gev_spatio_temporal_extremes_mle.py b/test/test_extreme_fit/test_estimator/test_gev_spatio_temporal_extremes_mle.py
index 387cc3e7..6e147115 100644
--- a/test/test_extreme_fit/test_estimator/test_gev_spatio_temporal_extremes_mle.py
+++ b/test/test_extreme_fit/test_estimator/test_gev_spatio_temporal_extremes_mle.py
@@ -2,12 +2,12 @@ import unittest
from random import sample
from extreme_data.meteo_france_data.scm_models_data.safran.safran import SafranSnowfall1Day, SafranPrecipitation1Day
-from extreme_fit.model.margin_model.polynomial_margin_model.altitudinal_models import \
+from extreme_fit.model.margin_model.polynomial_margin_model.gev_altitudinal_models import \
NonStationaryAltitudinalLocationQuadraticScaleLinearCrossTermForLocation, \
NonStationaryAltitudinalLocationQuadraticCrossTermForLocation, NonStationaryAltitudinalLocationLinear, \
NonStationaryAltitudinalLocationLinearCrossTermForLocation
-from extreme_fit.model.margin_model.polynomial_margin_model.utils import ALTITUDINAL_MODELS, \
- MODELS_THAT_SHOULD_RAISE_AN_ASSERTION_ERROR, VARIOUS_SPATIO_TEMPORAL_MODELS
+from extreme_fit.model.margin_model.polynomial_margin_model.utils import ALTITUDINAL_GEV_MODELS, \
+ MODELS_THAT_SHOULD_RAISE_AN_ASSERTION_ERROR, VARIOUS_SPATIO_TEMPORAL_MODELS, ALTITUDINAL_GUMBEL_MODELS
from extreme_fit.model.margin_model.utils import \
MarginFitMethod
from projects.altitude_spatial_model.altitudes_fit.altitudes_studies import AltitudesStudies
@@ -52,8 +52,13 @@ class TestGevTemporalQuadraticExtremesMle(unittest.TestCase):
for model_class in VARIOUS_SPATIO_TEMPORAL_MODELS:
self.common_test(model_class)
- def test_altitudinal_models(self):
- for model_class in ALTITUDINAL_MODELS:
+ def test_altitudinal_gev_models(self):
+ for model_class in ALTITUDINAL_GEV_MODELS:
+ self.common_test(model_class)
+
+ def test_altitudinal_gumbel_models(self):
+ for model_class in ALTITUDINAL_GUMBEL_MODELS[:3]:
+ # print(model_class)
self.common_test(model_class)
--
GitLab