diff --git a/extreme_fit/model/margin_model/polynomial_margin_model/gev_altitudinal_models.py b/extreme_fit/model/margin_model/polynomial_margin_model/gev_altitudinal_models.py
index f5612e2b6928d7fe23d7b28c2d9e59526914232e..dadc236f55c180066251165a99bfe93797bb3881 100644
--- a/extreme_fit/model/margin_model/polynomial_margin_model/gev_altitudinal_models.py
+++ b/extreme_fit/model/margin_model/polynomial_margin_model/gev_altitudinal_models.py
@@ -43,7 +43,11 @@ class AbstractAltitudinalModel(AbstractSpatioTemporalPolynomialModel):
s += '^2'
if self.dim_to_str_number(GevParams.SCALE, self.coordinates.idx_temporal_coordinates) in ['1', '2']:
s += '\\sigma_t'
- if self.dim_to_str_number(GevParams.LOC, self.coordinates.idx_temporal_coordinates) == '2':
+ if self.dim_to_str_number(GevParams.SCALE, self.coordinates.idx_temporal_coordinates) == '2':
+ s += '^2'
+ if self.dim_to_str_number(GevParams.SHAPE, self.coordinates.idx_temporal_coordinates) in ['1', '2']:
+ s += '\\zeta_t'
+ if self.dim_to_str_number(GevParams.SHAPE, self.coordinates.idx_temporal_coordinates) == '2':
s += '^2'
if len(s) == 0:
s = '0'
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 3675b84d8d20e11a685858721cbc50d11aafa693..f6e37ed3c7b31d53208b0a2ecb5a4447d19d237b 100644
--- a/extreme_fit/model/margin_model/polynomial_margin_model/utils.py
+++ b/extreme_fit/model/margin_model/polynomial_margin_model/utils.py
@@ -105,18 +105,22 @@ ALTITUDINAL_GEV_MODELS_BASED_ON_POINTWISE_ANALYSIS = [
StationaryAltitudinal,
AltitudinalShapeConstantTimeLocationLinear,
AltitudinalShapeConstantTimeScaleLinear,
- # AltitudinalShapeConstantTimeShapeLinear,
- # AltitudinalShapeConstantTimeLocShapeLinear,
+
AltitudinalShapeConstantTimeLocScaleLinear,
- # AltitudinalShapeConstantTimeScaleShapeLinear,
- # AltitudinalShapeConstantTimeLocScaleShapeLinear,
+
# With a linear shape for the altitude (8 models)
AltitudinalShapeLinearTimeStationary,
AltitudinalShapeLinearTimeLocationLinear,
AltitudinalShapeLinearTimeScaleLinear,
+
+ AltitudinalShapeLinearTimeLocScaleLinear,
+
+ # AltitudinalShapeConstantTimeShapeLinear,
+ # AltitudinalShapeConstantTimeLocShapeLinear,
+ # AltitudinalShapeConstantTimeScaleShapeLinear,
+ # AltitudinalShapeConstantTimeLocScaleShapeLinear,
# AltitudinalShapeLinearTimeShapeLinear,
# AltitudinalShapeLinearTimeLocShapeLinear,
- AltitudinalShapeLinearTimeLocScaleLinear,
# AltitudinalShapeLinearTimeScaleShapeLinear,
# AltitudinalShapeLinearTimeLocScaleShapeLinear,
diff --git a/extreme_fit/model/result_from_model_fit/result_from_extremes/abstract_result_from_extremes.py b/extreme_fit/model/result_from_model_fit/result_from_extremes/abstract_result_from_extremes.py
index ec86cb4ae69b294a789f727ced9bd03280e359c0..4870914970d9502bd0266d3d6967617ea5653f02 100644
--- a/extreme_fit/model/result_from_model_fit/result_from_extremes/abstract_result_from_extremes.py
+++ b/extreme_fit/model/result_from_model_fit/result_from_extremes/abstract_result_from_extremes.py
@@ -62,6 +62,7 @@ class AbstractResultFromExtremes(AbstractResultFromModelFit):
if self.param_name_to_dim:
d = {GevParams.greek_letter_from_param_name(param_name) + '1': r.c(transformed_temporal_covariate) for
param_name in self.param_name_to_dim.keys()}
+ print(d)
kwargs = {
"vals": r.list(**d
)
diff --git a/projects/altitude_spatial_model/altitudes_fit/altitudes_studies.py b/projects/altitude_spatial_model/altitudes_fit/altitudes_studies.py
index fff0c98828714a37e5c6693898d021b31485ace7..6221655315816d1c053a51adf8f16a65839d174c 100644
--- a/projects/altitude_spatial_model/altitudes_fit/altitudes_studies.py
+++ b/projects/altitude_spatial_model/altitudes_fit/altitudes_studies.py
@@ -139,9 +139,11 @@ class AltitudesStudies(object):
plot_name = 'Annual maxima of {} in {}'.format(SCM_STUDY_CLASS_TO_ABBREVIATION[self.study_class],
massif_name.replace('_', ' '))
ax.set_ylabel('{} ({})'.format(plot_name, self.study.variable_unit), fontsize=15)
+ plot_name = 'time series/' + plot_name
ax.set_xlabel('years', fontsize=15)
self.show_or_save_to_file(plot_name=plot_name, show=show, no_title=True)
ax.clear()
+ plt.close()
def plot_mean_maxima_against_altitude(self, massif_names=None, show=False, std=False, change=False):
ax = plt.gca()
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 2f8106d0599e0f973943c8f7ef7f0725de9fecc0..af7bb0a0fdfc2c17bb0ffc03fa5db447886a1a1e 100644
--- a/projects/altitude_spatial_model/altitudes_fit/main_altitudes_studies.py
+++ b/projects/altitude_spatial_model/altitudes_fit/main_altitudes_studies.py
@@ -47,7 +47,7 @@ def main():
# massif_names = ['Mercantour', 'Vercors', 'Ubaye']
seasons = [Season.annual, Season.winter, Season.spring, Season.automn][:1]
- main_loop(altitudes_for_groups[:1], massif_names, seasons, study_classes)
+ main_loop(altitudes_for_groups[:], massif_names, seasons, study_classes)
def main_loop(altitudes_list, massif_names, seasons, study_classes):
@@ -68,6 +68,7 @@ def main_loop(altitudes_list, massif_names, seasons, study_classes):
time.sleep(2)
+
def load_visualizer_list(season, study_class, altitudes_list, massif_names):
model_classes = ALTITUDINAL_GEV_MODELS_BASED_ON_POINTWISE_ANALYSIS
visualizer_list = []
@@ -107,7 +108,7 @@ def plots(massif_names, season, visualizer):
print('inner loop', season, type(studies.study))
# Plot time series
- # studies.plot_maxima_time_series(massif_names=massif_names)
+ studies.plot_maxima_time_series(massif_names=massif_names)
# Plot moments
# for std in [True, False][:]:
diff --git a/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/altitude_group.py b/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/altitude_group.py
index 9092fac1a0adae1351204cb87dc8962104682034..786a8684223d0bf7ec0c78c69edc405dff4508e5 100644
--- a/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/altitude_group.py
+++ b/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/altitude_group.py
@@ -5,8 +5,10 @@ altitudes_for_groups = [
[300, 600, 900],
[1200, 1500, 1800],
[2100, 2400, 2700],
- [3000, 3300, 3600]
+ [3000, 3300, 3600, 3900]
]
+
+
# altitudes_for_groups = [
# [900, 1200, 1500],
# [1800, 2100, 2400],
@@ -63,6 +65,7 @@ class HighAltitudeGroup(AbstractAltitudeGroup):
def reference_altitude(self):
return 2400
+
class VeyHighAltitudeGroup(AbstractAltitudeGroup):
@property
diff --git a/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/altitudes_studies_visualizer_for_non_stationary_models.py b/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/altitudes_studies_visualizer_for_non_stationary_models.py
index 2674d546dd1879ee9898b874dc385c802e765149..a963c82d82ee661e6d05ab1eab49715828cd6f16 100644
--- a/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/altitudes_studies_visualizer_for_non_stationary_models.py
+++ b/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/altitudes_studies_visualizer_for_non_stationary_models.py
@@ -21,7 +21,7 @@ from extreme_fit.model.margin_model.polynomial_margin_model.spatio_temporal_poly
from extreme_fit.model.margin_model.utils import MarginFitMethod
from projects.altitude_spatial_model.altitudes_fit.altitudes_studies import AltitudesStudies
from projects.altitude_spatial_model.altitudes_fit.one_fold_analysis.altitude_group import \
- get_altitude_group_from_altitudes, HighAltitudeGroup
+ get_altitude_group_from_altitudes, HighAltitudeGroup, VeyHighAltitudeGroup
from projects.altitude_spatial_model.altitudes_fit.one_fold_analysis.one_fold_fit import \
OneFoldFit
from spatio_temporal_dataset.coordinates.abstract_coordinates import AbstractCoordinates
@@ -50,12 +50,12 @@ class AltitudesStudiesVisualizerForNonStationaryModels(StudyVisualizer):
# Load one fold fit
self._massif_name_to_one_fold_fit = {}
for massif_name in self.massif_names:
- if any([massif_name in study.study_massif_names for study in self.studies.altitude_to_study.values()]):
+ if self.load_condition(massif_name):
# Load dataset
dataset = studies.spatio_temporal_dataset(massif_name=massif_name)
- # dataset_without_zero = AbstractDataset.remove_zeros(dataset.observations,
- # dataset.coordinates)
- old_fold_fit = OneFoldFit(massif_name, dataset, model_classes, self.fit_method,
+ dataset_without_zero = AbstractDataset.remove_zeros(dataset.observations,
+ dataset.coordinates)
+ old_fold_fit = OneFoldFit(massif_name, dataset_without_zero, model_classes, self.fit_method,
self.temporal_covariate_for_fit,
type(self.altitude_group),
self.display_only_model_that_pass_anderson_test)
@@ -66,7 +66,12 @@ class AltitudesStudiesVisualizerForNonStationaryModels(StudyVisualizer):
self._max_abs_for_shape = None
moment_names = ['moment', 'changes_for_moment', 'relative_changes_for_moment'][:]
- orders = [1, 2, None][:]
+ orders = [1, 2, None][2:]
+
+ def load_condition(self, massif_name):
+ altitudes_for_massif = [altitude for altitude, study in self.studies.altitude_to_study.items()
+ if massif_name in study.study_massif_names]
+ return (self.altitude_group.reference_altitude in altitudes_for_massif) and (len(altitudes_for_massif) >= 2)
@property
def massif_name_to_one_fold_fit(self) -> Dict[str, OneFoldFit]:
@@ -132,7 +137,7 @@ class AltitudesStudiesVisualizerForNonStationaryModels(StudyVisualizer):
@property
def add_colorbar(self):
- return isinstance(self.altitude_group, HighAltitudeGroup)
+ return isinstance(self.altitude_group, VeyHighAltitudeGroup)
def plot_against_years(self, method_name, order):
ax = plt.gca()
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 11cfef9587ec541da2d83309dac3194d40be604d..cb2d665122942023f1029fca96ef092f9732e7e8 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
@@ -15,6 +15,10 @@ from extreme_fit.model.margin_model.polynomial_margin_model.gumbel_altitudinal_m
from extreme_fit.model.margin_model.polynomial_margin_model.models_based_on_pariwise_analysis.gev_with_linear_shape_wrt_altitude import \
AltitudinalShapeLinearTimeStationary
from extreme_fit.model.margin_model.utils import MarginFitMethod
+from extreme_fit.model.result_from_model_fit.result_from_extremes.confidence_interval_method import \
+ ConfidenceIntervalMethodFromExtremes
+from extreme_fit.model.result_from_model_fit.result_from_extremes.eurocode_return_level_uncertainties import \
+ EurocodeConfidenceIntervalFromExtremes
from projects.altitude_spatial_model.altitudes_fit.one_fold_analysis.altitude_group import AbstractAltitudeGroup, \
DefaultAltitudeGroup
from root_utils import classproperty
@@ -110,6 +114,17 @@ class OneFoldFit(object):
sorted_estimators = sorted([estimator for estimator in estimators], key=lambda e: e.aic())
return sorted_estimators
+ @cached_property
+ def sorted_estimators_with_stationary(self):
+ if self.only_models_that_pass_anderson_test:
+ return [e for e in self.sorted_estimators if self.goodness_of_fit_test(e)]
+ else:
+ return self._sorted_estimators_without_stationary
+
+ @property
+ def has_at_least_one_valid_model(self):
+ return len(self.sorted_estimators_with_stationary) > 0
+
@cached_property
def _sorted_estimators_without_stationary(self):
return [e for e in self.sorted_estimators if not isinstance(e.margin_model, StationaryAltitudinal)]
@@ -134,8 +149,13 @@ class OneFoldFit(object):
if self.best_estimator_minimizes_total_aic and self.best_estimator_class_for_total_aic is not None:
return self.model_class_to_estimator[self.best_estimator_class_for_total_aic]
else:
- if self.has_at_least_one_valid_non_stationary_model:
- best_estimator = self.sorted_estimators_without_stationary[0]
+ # Without stationary
+ # if self.has_at_least_one_valid_non_stationary_model:
+ # best_estimator = self.sorted_estimators_without_stationary[0]
+ # return best_estimator
+ # With stationary
+ if self.has_at_least_one_valid_model:
+ best_estimator = self.sorted_estimators_with_stationary[0]
return best_estimator
else:
raise ValueError('This should not happen')
@@ -169,7 +189,10 @@ class OneFoldFit(object):
def best_name(self):
name = self.best_estimator.margin_model.name_str
latex_command = 'textbf' if self.is_significant else 'textrm'
- return '$\\' + latex_command + '{' + name + '}$'
+ best_name = '$\\' + latex_command + '{' + name + '}$'
+ if self.is_significant:
+ best_name = '\\underline{' + best_name + '}'
+ return best_name
# Significant
@@ -181,6 +204,8 @@ class OneFoldFit(object):
return self.model_class_to_estimator_with_finite_aic[StationaryAltitudinalOnlyScale]
elif isinstance(self.best_estimator.margin_model, AltitudinalShapeLinearTimeStationary):
return self.model_class_to_estimator_with_finite_aic[AltitudinalShapeLinearTimeStationary]
+ elif isinstance(self.best_estimator.margin_model, AltitudinalShapeLinearTimeStationary):
+ return self.model_class_to_estimator_with_finite_aic[AltitudinalShapeLinearTimeStationary]
else:
return self.model_class_to_estimator_with_finite_aic[StationaryAltitudinal]
@@ -194,7 +219,7 @@ class OneFoldFit(object):
@property
def is_significant(self) -> bool:
- stationary_model_classes = [StationaryAltitudinal, StationaryGumbelAltitudinal]
+ stationary_model_classes = [StationaryAltitudinal, StationaryGumbelAltitudinal, AltitudinalShapeLinearTimeStationary]
if any([isinstance(self.best_estimator.margin_model, c)
for c in stationary_model_classes]):
return False
@@ -233,3 +258,8 @@ class OneFoldFit(object):
empirical_quantiles.append(maximum_standardized)
empirical_quantiles = sorted(empirical_quantiles)
return empirical_quantiles
+
+ # def best_confidence_interval(self):
+ # EurocodeConfidenceIntervalFromExtremes.from_estimator_extremes(self.best_estimator,
+ # ci_method=ConfidenceIntervalMethodFromExtremes.ci_mle,
+ # temporal_covariate=np.array([2019, self.altitude_plot]),)
\ No newline at end of file