diff --git a/extreme_fit/model/margin_model/linear_margin_model/temporal_linear_margin_models.py b/extreme_fit/model/margin_model/linear_margin_model/temporal_linear_margin_models.py index 989db4712c9aa57ece5b431523912ac615137ddb..c9626a486720972dcd780d276da4b44dcf2f0dec 100644 --- a/extreme_fit/model/margin_model/linear_margin_model/temporal_linear_margin_models.py +++ b/extreme_fit/model/margin_model/linear_margin_model/temporal_linear_margin_models.py @@ -57,7 +57,7 @@ class NonStationaryLocationAndScaleTemporalModel(AbstractTemporalLinearMarginMod def load_margin_function(self, param_name_to_dims=None): return super().load_margin_function({GevParams.LOC: [self.coordinates.idx_temporal_coordinates], - GevParams.SCALE: [self.coordinates.idx_temporal_coordinates]}) + GevParams.SCALE: [self.coordinates.idx_temporal_coordinates]}) @property def mul(self): @@ -68,6 +68,56 @@ class NonStationaryLocationAndScaleTemporalModel(AbstractTemporalLinearMarginMod return 1 +class NonStationaryLocationAndShapeTemporalModel(AbstractTemporalLinearMarginModel): + + def load_margin_function(self, param_name_to_dims=None): + return super().load_margin_function({GevParams.LOC: [self.coordinates.idx_temporal_coordinates], + GevParams.SHAPE: [self.coordinates.idx_temporal_coordinates]}) + + @property + def mul(self): + return 1 + + @property + def shl(self): + return 1 + + +class NonStationaryScaleAndShapeTemporalModel(AbstractTemporalLinearMarginModel): + + def load_margin_function(self, param_name_to_dims=None): + return super().load_margin_function({GevParams.SHAPE: [self.coordinates.idx_temporal_coordinates], + GevParams.SCALE: [self.coordinates.idx_temporal_coordinates]}) + + @property + def shl(self): + return 1 + + @property + def sigl(self): + return 1 + + +class NonStationaryLocationAndScaleAndShapeTemporalModel(AbstractTemporalLinearMarginModel): + + def load_margin_function(self, param_name_to_dims=None): + return super().load_margin_function({GevParams.LOC: [self.coordinates.idx_temporal_coordinates], + GevParams.SCALE: [self.coordinates.idx_temporal_coordinates], + GevParams.SHAPE: [self.coordinates.idx_temporal_coordinates]}) + + @property + def mul(self): + return 1 + + @property + def sigl(self): + return 1 + + @property + def shl(self): + return 1 + + class GumbelTemporalModel(StationaryTemporalModel): def __init__(self, coordinates: AbstractCoordinates, diff --git a/extreme_trend/trend_test_one_parameter/gev_trend_test_one_parameter.py b/extreme_trend/trend_test_one_parameter/gev_trend_test_one_parameter.py index e23df33795f03adee5ecf17181cfe9fa5b5fe31c..e0f942a18800346c9b742e072bbc20c8ca04a16c 100644 --- a/extreme_trend/trend_test_one_parameter/gev_trend_test_one_parameter.py +++ b/extreme_trend/trend_test_one_parameter/gev_trend_test_one_parameter.py @@ -6,6 +6,7 @@ from extreme_fit.model.margin_model.linear_margin_model.temporal_linear_margin_m NonStationaryLocationTemporalModel, NonStationaryScaleTemporalModel, NonStationaryShapeTemporalModel, \ StationaryTemporalModel from extreme_fit.distribution.gev.gev_params import GevParams +from root_utils import classproperty class GevTrendTestOneParameter(AbstractGevTrendTest): @@ -15,6 +16,28 @@ class GevTrendTestOneParameter(AbstractGevTrendTest): return 1 +class GevVersusGev(GevTrendTestOneParameter): + + def __init__(self, years, maxima, starting_year, quantile_level=EUROCODE_QUANTILE, + fit_method=MarginFitMethod.extremes_fevd_mle): + super().__init__(years, maxima, starting_year, + unconstrained_model_class=StationaryTemporalModel, + constrained_model_class=StationaryTemporalModel, + quantile_level=quantile_level, + fit_method=fit_method) + + @property + def is_significant(self) -> bool: + return False + + @classproperty + def total_number_of_parameters_for_unconstrained_model(cls) -> int: + return 3 + + def _slope_strength(self): + return 0.0 + + class GevTrendTestOneParameterAgainstStationary(GevTrendTestOneParameter): def __init__(self, years, maxima, starting_year, unconstrained_model_class, param_name, @@ -32,10 +55,15 @@ class GevTrendTestOneParameterAgainstStationary(GevTrendTestOneParameter): def non_stationary_linear_coef(self): return self.get_non_stationary_linear_coef(param_name=self.param_name) + @classproperty + def total_number_of_parameters_for_unconstrained_model(cls) -> int: + return 4 + class GevLocationTrendTest(GevTrendTestOneParameterAgainstStationary): - def __init__(self, years, maxima, starting_year, quantile_level=EUROCODE_QUANTILE, constrained_model_class=StationaryTemporalModel, + def __init__(self, years, maxima, starting_year, quantile_level=EUROCODE_QUANTILE, + constrained_model_class=StationaryTemporalModel, fit_method=MarginFitMethod.extremes_fevd_mle): super().__init__(years, maxima, starting_year, unconstrained_model_class=NonStationaryLocationTemporalModel, @@ -58,9 +86,11 @@ class GevLocationTrendTest(GevTrendTestOneParameterAgainstStationary): def variance_difference_same_sign_as_slope_strenght(self) -> bool: return False + class GevScaleTrendTest(GevTrendTestOneParameterAgainstStationary): - def __init__(self, years, maxima, starting_year, quantile_level=EUROCODE_QUANTILE, constrained_model_class=StationaryTemporalModel, + def __init__(self, years, maxima, starting_year, quantile_level=EUROCODE_QUANTILE, + constrained_model_class=StationaryTemporalModel, fit_method=MarginFitMethod.extremes_fevd_mle): super().__init__(years, maxima, starting_year, unconstrained_model_class=NonStationaryScaleTemporalModel, @@ -88,11 +118,10 @@ class GevScaleTrendTest(GevTrendTestOneParameterAgainstStationary): class GevShapeTrendTest(GevTrendTestOneParameterAgainstStationary): - def __init__(self, years, maxima, starting_year, quantile_level=EUROCODE_QUANTILE, fit_method=MarginFitMethod.extremes_fevd_mle): + def __init__(self, years, maxima, starting_year, quantile_level=EUROCODE_QUANTILE, + fit_method=MarginFitMethod.extremes_fevd_mle): super().__init__(years, maxima, starting_year, unconstrained_model_class=NonStationaryShapeTemporalModel, param_name=GevParams.SHAPE, quantile_level=quantile_level, fit_method=fit_method) - - diff --git a/extreme_trend/trend_test_three_parameters/gev_trend_test_three_parameters.py b/extreme_trend/trend_test_three_parameters/gev_trend_test_three_parameters.py index e69fd8a43036bcf7d05d2f3550981bf11b9de975..caaed063997b9531984cf46771bb90b7d3971440 100644 --- a/extreme_trend/trend_test_three_parameters/gev_trend_test_three_parameters.py +++ b/extreme_trend/trend_test_three_parameters/gev_trend_test_three_parameters.py @@ -3,7 +3,8 @@ from extreme_fit.model.margin_model.utils import \ from extreme_data.eurocode_data.utils import EUROCODE_QUANTILE from extreme_trend.abstract_gev_trend_test import AbstractGevTrendTest from extreme_fit.model.margin_model.linear_margin_model.temporal_linear_margin_models import \ - NonStationaryLocationAndScaleTemporalModel, GumbelTemporalModel + NonStationaryLocationAndScaleTemporalModel, GumbelTemporalModel, StationaryTemporalModel, \ + NonStationaryLocationAndScaleAndShapeTemporalModel from extreme_fit.distribution.gev.gev_params import GevParams from root_utils import classproperty @@ -15,9 +16,25 @@ class GevTrendTestThreeParameters(AbstractGevTrendTest): return 3 +class GevLocationAndScaleAndShapeTrendTest(GevTrendTestThreeParameters): + + def __init__(self, years, maxima, starting_year, quantile_level=EUROCODE_QUANTILE, + fit_method=MarginFitMethod.extremes_fevd_mle): + super().__init__(years, maxima, starting_year, + unconstrained_model_class=StationaryTemporalModel, + constrained_model_class=NonStationaryLocationAndScaleAndShapeTemporalModel, + quantile_level=quantile_level, + fit_method=fit_method) + + @classproperty + def total_number_of_parameters_for_unconstrained_model(cls) -> int: + return 6 + + class GevLocationAndScaleTrendTestAgainstGumbel(GevTrendTestThreeParameters): - def __init__(self, years, maxima, starting_year, quantile_level=EUROCODE_QUANTILE, fit_method=MarginFitMethod.extremes_fevd_mle): + def __init__(self, years, maxima, starting_year, quantile_level=EUROCODE_QUANTILE, + fit_method=MarginFitMethod.extremes_fevd_mle): super().__init__(years, maxima, starting_year, unconstrained_model_class=NonStationaryLocationAndScaleTemporalModel, constrained_model_class=GumbelTemporalModel, @@ -47,4 +64,4 @@ class GevLocationAndScaleTrendTestAgainstGumbel(GevTrendTestThreeParameters): @classproperty def total_number_of_parameters_for_unconstrained_model(cls) -> int: - return 5 \ No newline at end of file + return 5 diff --git a/extreme_trend/trend_test_two_parameters/gev_trend_test_two_parameters.py b/extreme_trend/trend_test_two_parameters/gev_trend_test_two_parameters.py index ae6b8f54f04afcfb2bfcaabfb1b3e2f9407e03ce..20b6edef969d212c49290a2bf1bb696a53f471cb 100644 --- a/extreme_trend/trend_test_two_parameters/gev_trend_test_two_parameters.py +++ b/extreme_trend/trend_test_two_parameters/gev_trend_test_two_parameters.py @@ -5,7 +5,8 @@ from extreme_trend.trend_test_one_parameter.gev_trend_test_one_parameter import from extreme_fit.model.margin_model.utils import \ MarginFitMethod from extreme_fit.model.margin_model.linear_margin_model.temporal_linear_margin_models import \ - NonStationaryLocationAndScaleTemporalModel, StationaryTemporalModel, GumbelTemporalModel + NonStationaryLocationAndScaleTemporalModel, StationaryTemporalModel, GumbelTemporalModel, \ + NonStationaryScaleAndShapeTemporalModel, NonStationaryLocationAndShapeTemporalModel from extreme_fit.distribution.gev.gev_params import GevParams from root_utils import classproperty @@ -16,8 +17,37 @@ class GevTrendTestTwoParameters(AbstractGevTrendTest): def degree_freedom_chi2(self) -> int: return 2 +class GevTrendTestTwoParametersAgainstGev(AbstractGevTrendTest): -class GevLocationAndScaleTrendTest(GevTrendTestTwoParameters): + @classproperty + def total_number_of_parameters_for_unconstrained_model(cls) -> int: + return 5 + + +class GevLocationAndShapeTrendTest(GevTrendTestTwoParametersAgainstGev): + + def __init__(self, years, maxima, starting_year, constrained_model_class=StationaryTemporalModel, + quantile_level=EUROCODE_QUANTILE, fit_method=MarginFitMethod.extremes_fevd_mle): + super().__init__(years, maxima, starting_year, + unconstrained_model_class=NonStationaryLocationAndShapeTemporalModel, + constrained_model_class=constrained_model_class, + quantile_level=quantile_level, + fit_method=fit_method) + + +class GevScaleAndShapeTrendTest(GevTrendTestTwoParametersAgainstGev): + + + def __init__(self, years, maxima, starting_year, constrained_model_class=StationaryTemporalModel, + quantile_level=EUROCODE_QUANTILE, fit_method=MarginFitMethod.extremes_fevd_mle): + super().__init__(years, maxima, starting_year, + unconstrained_model_class=NonStationaryScaleAndShapeTemporalModel, + constrained_model_class=constrained_model_class, + quantile_level=quantile_level, + fit_method=fit_method) + + +class GevLocationAndScaleTrendTest(GevTrendTestTwoParametersAgainstGev): def __init__(self, years, maxima, starting_year, constrained_model_class=StationaryTemporalModel, quantile_level=EUROCODE_QUANTILE, fit_method=MarginFitMethod.extremes_fevd_mle): diff --git a/extreme_trend/visualizers/study_visualizer_for_non_stationary_trends.py b/extreme_trend/visualizers/study_visualizer_for_non_stationary_trends.py index b91631b3304b77e614d4733f6c9a84b7e553fb84..1376bb39104e9c447a39b335eb503035f3e221b3 100644 --- a/extreme_trend/visualizers/study_visualizer_for_non_stationary_trends.py +++ b/extreme_trend/visualizers/study_visualizer_for_non_stationary_trends.py @@ -17,7 +17,7 @@ from extreme_data.meteo_france_data.scm_models_data.abstract_extended_study impo from extreme_data.meteo_france_data.scm_models_data.abstract_study import AbstractStudy from extreme_data.meteo_france_data.scm_models_data.visualization.study_visualizer import \ StudyVisualizer -from projects.exceeding_snow_loads.utils import NON_STATIONARY_TREND_TEST_PAPER +from projects.exceeding_snow_loads.utils import NON_STATIONARY_TREND_TEST_PAPER_1 from extreme_trend.abstract_gev_trend_test import AbstractGevTrendTest from extreme_trend.trend_test_one_parameter.gumbel_trend_test_one_parameter import \ GumbelLocationTrendTest, GevStationaryVersusGumbel, GumbelScaleTrendTest, GumbelVersusGumbel @@ -87,7 +87,7 @@ class StudyVisualizerForNonStationaryTrends(StudyVisualizer): assert set(self.uncertainty_massif_names).issubset(set(self.study.study_massif_names)) if self.non_stationary_trend_test_to_marker is None: # Assign default argument for the non stationary trends - self.non_stationary_trend_test = NON_STATIONARY_TREND_TEST_PAPER + self.non_stationary_trend_test = NON_STATIONARY_TREND_TEST_PAPER_1 self.non_stationary_trend_test_to_marker = {t: t.marker for t in self.non_stationary_trend_test} else: self.non_stationary_trend_test = list(self.non_stationary_trend_test_to_marker.keys()) diff --git a/projects/contrasting_trends_in_snow_loads/article2_snowfall_versus_time_and_altitude/data.py b/projects/contrasting_trends_in_snow_loads/article2_snowfall_versus_time_and_altitude/data.py index af194e84cf7fec16eb42b944d1b70d92f434e478..aff3b1848ffe9c1918ed2dd546360539a18b72fc 100644 --- a/projects/contrasting_trends_in_snow_loads/article2_snowfall_versus_time_and_altitude/data.py +++ b/projects/contrasting_trends_in_snow_loads/article2_snowfall_versus_time_and_altitude/data.py @@ -63,7 +63,7 @@ def max_graph_annual_maxima_poster(): tight_pad = {'h_pad': 0.2} snow_abbreviation = 'max ' + snow_abbreviation for color, _, massif_name in marker_altitude_massif_name_for_paper1[::-1]: - last_plot = color == "magenta" and altitude == second_altitude + last_plot = massif_name == 'Ubaye' and altitude == second_altitude label = '{} massif at {}m'.format(massif_name, altitude) study_visualizer.visualize_max_graphs_poster(massif_name, altitude, snow_abbreviation, color, label, last_plot, ax, tight_pad=tight_pad, diff --git a/projects/contrasting_trends_in_snow_loads/article2_snowfall_versus_time_and_altitude/main_snowfall_article.py b/projects/contrasting_trends_in_snow_loads/article2_snowfall_versus_time_and_altitude/main_snowfall_article.py index 5a9ccf61c32ae9b8232d0d7de9ff0acb9876cd22..06759a89d404ad0c400bdd248074d6e9eee0cbc7 100644 --- a/projects/contrasting_trends_in_snow_loads/article2_snowfall_versus_time_and_altitude/main_snowfall_article.py +++ b/projects/contrasting_trends_in_snow_loads/article2_snowfall_versus_time_and_altitude/main_snowfall_article.py @@ -72,10 +72,10 @@ def intermediate_result(altitudes, massif_names=None, _ = compute_minimized_aic(visualizer) # Plots - validation_plot(altitude_to_visualizer, order_derivative=0) - validation_plot(altitude_to_visualizer, order_derivative=1) + # validation_plot(altitude_to_visualizer, order_derivative=0) + # validation_plot(altitude_to_visualizer, order_derivative=1) plot_snowfall_mean(altitude_to_visualizer) - plot_snowfall_time_derivative_mean(altitude_to_visualizer) + # plot_snowfall_time_derivative_mean(altitude_to_visualizer) def major_result(): @@ -86,7 +86,7 @@ def major_result(): model_subsets_for_uncertainty = None altitudes = paper_altitudes altitudes = [300, 600, 900, 1200, 1500, 1800, 2100, 2400, 2700, 3000, 3300, 3600, 3900] - # altitudes = [900, 1200] + altitudes = [900, 1200] for study_class in study_classes: intermediate_result(altitudes, massif_names, model_subsets_for_uncertainty, diff --git a/projects/exceeding_snow_loads/utils.py b/projects/exceeding_snow_loads/utils.py index 0847e667c23c2882260f81c2ea4562185081c046..0faed07afe2d6b2ff0232d726e32fea91d9adb4d 100644 --- a/projects/exceeding_snow_loads/utils.py +++ b/projects/exceeding_snow_loads/utils.py @@ -4,12 +4,15 @@ from extreme_data.meteo_france_data.scm_models_data.crocus.crocus import CrocusS CrocusSnowLoad3Days from extreme_data.meteo_france_data.scm_models_data.visualization.main_study_visualizer import \ ALL_ALTITUDES_WITHOUT_NAN +from extreme_trend.trend_test_one_parameter.gev_trend_test_one_parameter import GevVersusGev, GevScaleTrendTest, \ + GevLocationTrendTest, GevShapeTrendTest from extreme_trend.trend_test_one_parameter.gumbel_trend_test_one_parameter import \ GumbelVersusGumbel, GumbelLocationTrendTest, GumbelScaleTrendTest, GevStationaryVersusGumbel from extreme_trend.trend_test_three_parameters.gev_trend_test_three_parameters import \ - GevLocationAndScaleTrendTestAgainstGumbel + GevLocationAndScaleTrendTestAgainstGumbel, GevLocationAndScaleAndShapeTrendTest from extreme_trend.trend_test_two_parameters.gev_trend_test_two_parameters import \ - GevLocationAgainstGumbel, GevScaleAgainstGumbel + GevLocationAgainstGumbel, GevScaleAgainstGumbel, GevLocationAndScaleTrendTest, GevScaleAndShapeTrendTest, \ + GevLocationAndShapeTrendTest from extreme_trend.trend_test_two_parameters.gumbel_test_two_parameters import \ GumbelLocationAndScaleTrendTest @@ -17,12 +20,20 @@ paper_altitudes = ALL_ALTITUDES_WITHOUT_NAN paper_study_classes = [CrocusSnowLoadTotal, CrocusSnowLoadEurocode, CrocusSnowLoad3Days][:2] # dpi_paper1_figure = 700 dpi_paper1_figure = None -NON_STATIONARY_TREND_TEST_PAPER = [GumbelVersusGumbel, - GumbelLocationTrendTest, GumbelScaleTrendTest, - GumbelLocationAndScaleTrendTest, - GevStationaryVersusGumbel, - GevLocationAgainstGumbel, GevScaleAgainstGumbel, - GevLocationAndScaleTrendTestAgainstGumbel] +NON_STATIONARY_TREND_TEST_PAPER_1 = [GumbelVersusGumbel, + GumbelLocationTrendTest, GumbelScaleTrendTest, + GumbelLocationAndScaleTrendTest, + GevStationaryVersusGumbel, + GevLocationAgainstGumbel, GevScaleAgainstGumbel, + GevLocationAndScaleTrendTestAgainstGumbel] +NON_STATIONARY_TREND_TEST_PAPER_2 = [ + # Gumbel models + GumbelVersusGumbel, GumbelLocationTrendTest, GumbelScaleTrendTest, GumbelLocationAndScaleTrendTest, + # GEV models with constant shape + GevVersusGev, GevLocationTrendTest, GevScaleTrendTest, GevLocationAndScaleTrendTest, + # GEV models with linear shape + GevShapeTrendTest, GevLocationAndShapeTrendTest, GevScaleAndShapeTrendTest, GevLocationAndScaleAndShapeTrendTest +] diff --git a/test/test_extreme_trend/test_extreme_trend.py b/test/test_extreme_trend/test_extreme_trend.py index e5bf0416f074da142635a9a6f66333358d80b1a3..efd2d0a821a96055a522368aa114e17277ca954e 100644 --- a/test/test_extreme_trend/test_extreme_trend.py +++ b/test/test_extreme_trend/test_extreme_trend.py @@ -1,16 +1,24 @@ import unittest -from projects.exceeding_snow_loads.utils import NON_STATIONARY_TREND_TEST_PAPER +from projects.exceeding_snow_loads.utils import NON_STATIONARY_TREND_TEST_PAPER_1, NON_STATIONARY_TREND_TEST_PAPER_2 class TestTrendAnalysis(unittest.TestCase): - def test_nb_parameters(self): - trend_test_classes = NON_STATIONARY_TREND_TEST_PAPER + def test_nb_parameters_paper1(self): + trend_test_classes = NON_STATIONARY_TREND_TEST_PAPER_1 nb_expected = [2, 3, 3, 4, 3, 4, 4, 5] for trend_test_class, nb in zip(trend_test_classes, nb_expected): self.assertEqual(trend_test_class.total_number_of_parameters_for_unconstrained_model, nb) + def test_nb_parameters_paper2(self): + trend_test_classes = NON_STATIONARY_TREND_TEST_PAPER_2 + nb_expected = [2, 3, 3, 4, + 3, 4, 4, 5, + 4, 5, 5, 6] + for trend_test_class, nb in zip(trend_test_classes, nb_expected): + self.assertEqual(trend_test_class.total_number_of_parameters_for_unconstrained_model, nb) + if __name__ == '__main__': unittest.main()