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 eb8ba9cf1786a152e928a19c3a37c6bc88bc2ca1..467f89068bdd07ee3beee618f4e8a12411de1a57 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 @@ -21,6 +21,7 @@ class TemporalMarginFitMethod(Enum): extremes_fevd_bayesian = 1 extremes_fevd_mle = 2 extremes_fevd_gmle = 3 + extremes_fevd_l_moments = 4 class AbstractTemporalLinearMarginModel(LinearMarginModel): @@ -53,7 +54,9 @@ class AbstractTemporalLinearMarginModel(LinearMarginModel): elif self.fit_method == TemporalMarginFitMethod.extremes_fevd_bayesian: return self.extremes_fevd_bayesian_fit(x, df_coordinates_temp) elif self.fit_method in [TemporalMarginFitMethod.extremes_fevd_mle, - TemporalMarginFitMethod.extremes_fevd_gmle]: + TemporalMarginFitMethod.extremes_fevd_gmle, + TemporalMarginFitMethod.extremes_fevd_l_moments, + ]: return self.extremes_fevd_mle_related_fit(x, df_coordinates_temp) else: raise NotImplementedError @@ -78,6 +81,9 @@ class AbstractTemporalLinearMarginModel(LinearMarginModel): method = "MLE" elif self.fit_method == TemporalMarginFitMethod.extremes_fevd_gmle: method = "GMLE" + elif self.fit_method == TemporalMarginFitMethod.extremes_fevd_l_moments: + method = "Lmoments" + assert self.margin_function_start_fit.is_a_stationary_model else: raise ValueError('wrong method') return self.run_fevd_fixed(df_coordinates_temp, method, x) diff --git a/extreme_fit/model/result_from_model_fit/result_from_extremes/result_from_mle_extremes.py b/extreme_fit/model/result_from_model_fit/result_from_extremes/result_from_mle_extremes.py index 3cf68a4673ec6eb6081b125dd02223bc1626a81e..9333ff44fc5a64efca60f61c603ba0f192dcc596 100644 --- a/extreme_fit/model/result_from_model_fit/result_from_extremes/result_from_mle_extremes.py +++ b/extreme_fit/model/result_from_model_fit/result_from_extremes/result_from_mle_extremes.py @@ -21,7 +21,10 @@ class ResultFromMleExtremes(AbstractResultFromExtremes): def margin_coef_ordered_dict(self): values = self.name_to_value['results'] d = self.get_python_dictionary(values) - values = {i: param for i, param in enumerate(np.array(d['par']))} + if 'par' in d: + values = {i: param for i, param in enumerate(np.array(d['par']))} + else: + values = {i: np.array(v)[0] for i, v in enumerate(d.values())} return get_margin_coef_ordered_dict(self.gev_param_name_to_dim, values, self.type_for_mle) def _confidence_interval_method(self, common_kwargs, ci_method, return_period): diff --git a/test/test_extreme_fit/test_distribution/test_gev/test_gev_temporal_extremes_l_moments.py b/test/test_extreme_fit/test_distribution/test_gev/test_gev_temporal_extremes_l_moments.py new file mode 100644 index 0000000000000000000000000000000000000000..c94c33ad467d7906a46c4bb6c00f30d55837b0e6 --- /dev/null +++ b/test/test_extreme_fit/test_distribution/test_gev/test_gev_temporal_extremes_l_moments.py @@ -0,0 +1,62 @@ +import unittest + +import numpy as np +import pandas as pd + +from extreme_trend.abstract_gev_trend_test import fitted_linear_margin_estimator +from extreme_fit.model.margin_model.linear_margin_model.abstract_temporal_linear_margin_model import \ + TemporalMarginFitMethod +from extreme_fit.model.margin_model.linear_margin_model.temporal_linear_margin_models import StationaryTemporalModel, \ + NonStationaryLocationTemporalModel, NonStationaryLocationAndScaleTemporalModel +from extreme_fit.model.utils import r, set_seed_r +from spatio_temporal_dataset.coordinates.abstract_coordinates import AbstractCoordinates +from spatio_temporal_dataset.coordinates.temporal_coordinates.abstract_temporal_coordinates import \ + AbstractTemporalCoordinates +from spatio_temporal_dataset.dataset.abstract_dataset import AbstractDataset +from spatio_temporal_dataset.spatio_temporal_observations.abstract_spatio_temporal_observations import \ + AbstractSpatioTemporalObservations + + +class TestGevTemporalExtremesLMoments(unittest.TestCase): + + def setUp(self) -> None: + set_seed_r() + r(""" + N <- 50 + loc = 0; scale = 1; shape <- 1 + x_gev <- rgev(N, loc = loc, scale = scale, shape = shape) + start_loc = 0; start_scale = 1; start_shape = 1 + """) + # Compute the stationary temporal margin with isMev + self.start_year = 0 + df = pd.DataFrame({AbstractCoordinates.COORDINATE_T: range(self.start_year, self.start_year + 50)}) + self.coordinates = AbstractTemporalCoordinates.from_df(df) + df2 = pd.DataFrame(data=np.array(r['x_gev']), index=df.index) + observations = AbstractSpatioTemporalObservations(df_maxima_gev=df2) + self.dataset = AbstractDataset(observations=observations, coordinates=self.coordinates) + self.fit_method = TemporalMarginFitMethod.extremes_fevd_l_moments + + ref = {'loc': 0.0813843045950251, 'scale': 1.1791830110181365, 'shape': 0.6610403806908737} + + def test_gev_temporal_margin_fit_stationary(self): + # Create estimator + estimator = fitted_linear_margin_estimator(StationaryTemporalModel, self.coordinates, self.dataset, + starting_year=0, + fit_method=self.fit_method) + ref = {'loc': 0.0813843045950251, 'scale': 1.1791830110181365, 'shape': 0.6610403806908737} + for year in range(1, 3): + mle_params_estimated = estimator.function_from_fit.get_gev_params(np.array([year])).to_dict() + for key in ref.keys(): + self.assertAlmostEqual(ref[key], mle_params_estimated[key], places=3) + + def test_gev_temporal_margin_fit_non_stationary_location(self): + # Create estimator + with self.assertRaises(AssertionError): + fitted_linear_margin_estimator(NonStationaryLocationTemporalModel, self.coordinates, + self.dataset, + starting_year=0, + fit_method=self.fit_method) + + +if __name__ == '__main__': + unittest.main()