From bfa9e9d82f66460b0ec8f042ec4897ef781a7a15 Mon Sep 17 00:00:00 2001
From: Le Roux Erwan <erwan.le-roux@irstea.fr>
Date: Fri, 20 Mar 2020 15:38:27 +0100
Subject: [PATCH] [quantile regression project] add test_gev_simulations.py.
 add basic version for the simulation

---
 .../quantile_estimator_from_margin.py         |   4 +-
 .../quantile_estimator_from_regression.py     |   3 +-
 .../AbstractSimulation.py                     | 104 ++++++++++++++++++
 .../GevSimulation.py                          |  18 +++
 .../test_estimator/test_quantile_estimator.py |   2 +-
 .../test_quantile_regression/__init__.py      |   0
 .../test_gev_simulations.py                   |  20 ++++
 7 files changed, 146 insertions(+), 5 deletions(-)
 create mode 100644 projects/quantile_regression_vs_evt/AbstractSimulation.py
 create mode 100644 projects/quantile_regression_vs_evt/GevSimulation.py
 create mode 100644 test/test_projects/test_quantile_regression/__init__.py
 create mode 100644 test/test_projects/test_quantile_regression/test_gev_simulations.py

diff --git a/extreme_fit/estimator/quantile_estimator/quantile_estimator_from_margin.py b/extreme_fit/estimator/quantile_estimator/quantile_estimator_from_margin.py
index fa130c1d..83763c91 100644
--- a/extreme_fit/estimator/quantile_estimator/quantile_estimator_from_margin.py
+++ b/extreme_fit/estimator/quantile_estimator/quantile_estimator_from_margin.py
@@ -22,8 +22,8 @@ from spatio_temporal_dataset.dataset.abstract_dataset import AbstractDataset
 
 class QuantileEstimatorFromMargin(LinearMarginEstimator, AbstractQuantileEstimator):
 
-    def __init__(self, dataset: AbstractDataset, margin_model: AbstractTemporalLinearMarginModel, quantile):
-        super().__init__(dataset=dataset, quantile=quantile, margin_model=margin_model)
+    def __init__(self, dataset: AbstractDataset, quantile, margin_model_class: type):
+        super().__init__(dataset=dataset, quantile=quantile, margin_model=margin_model_class(dataset.coordinates))
 
     @cached_property
     def function_from_fit(self) -> AbstractQuantileFunction:
diff --git a/extreme_fit/estimator/quantile_estimator/quantile_estimator_from_regression.py b/extreme_fit/estimator/quantile_estimator/quantile_estimator_from_regression.py
index 6cdc634c..8b7d2b0c 100644
--- a/extreme_fit/estimator/quantile_estimator/quantile_estimator_from_regression.py
+++ b/extreme_fit/estimator/quantile_estimator/quantile_estimator_from_regression.py
@@ -24,8 +24,7 @@ class QuantileRegressionEstimator(AbstractQuantileEstimator):
 
     def __init__(self, dataset: AbstractDataset, quantile: float, quantile_regression_model_class: type, **kwargs):
         super().__init__(dataset, quantile, **kwargs)
-        self.quantile_regression_model = quantile_regression_model_class(dataset,
-                                                                         quantile)  # type: AbstractQuantileRegressionModel
+        self.quantile_regression_model = quantile_regression_model_class(dataset, quantile)
 
     def _fit(self) -> AbstractResultFromModelFit:
         return self.quantile_regression_model.fit()
diff --git a/projects/quantile_regression_vs_evt/AbstractSimulation.py b/projects/quantile_regression_vs_evt/AbstractSimulation.py
new file mode 100644
index 00000000..3a013a04
--- /dev/null
+++ b/projects/quantile_regression_vs_evt/AbstractSimulation.py
@@ -0,0 +1,104 @@
+from typing import Dict, List
+import matplotlib.pyplot as plt
+from collections import OrderedDict
+
+import numpy as np
+from cached_property import cached_property
+
+from extreme_fit.estimator.quantile_estimator.quantile_estimator_from_margin import QuantileEstimatorFromMargin
+from extreme_fit.estimator.quantile_estimator.quantile_estimator_from_regression import QuantileRegressionEstimator
+from extreme_fit.model.margin_model.linear_margin_model.abstract_temporal_linear_margin_model import \
+    AbstractTemporalLinearMarginModel
+from extreme_fit.model.quantile_model.quantile_regression_model import AbstractQuantileRegressionModel
+from spatio_temporal_dataset.coordinates.temporal_coordinates.generated_temporal_coordinates import \
+    ConsecutiveTemporalCoordinates
+from spatio_temporal_dataset.dataset.abstract_dataset import AbstractDataset
+from spatio_temporal_dataset.spatio_temporal_observations.abstract_spatio_temporal_observations import \
+    AbstractSpatioTemporalObservations
+
+
+class Coordinates(object):
+    pass
+
+
+class AbstractSimulation(object):
+
+    def __init__(self, nb_time_series, quantile, time_series_lengths=None, multiprocessing=False,
+                 model_classes=None):
+        self.models_classes = model_classes
+        self.multiprocessing = multiprocessing
+        self.quantile = quantile
+        self.time_series_lengths = time_series_lengths
+        self.nb_time_series = nb_time_series
+
+    def generate_all_observation(self, nb_time_series, length, coordinates) -> List[AbstractSpatioTemporalObservations]:
+        raise NotImplementedError
+
+    @cached_property
+    def time_serie_length_to_observation_list(self) -> Dict[int, List[AbstractSpatioTemporalObservations]]:
+        d = OrderedDict()
+        for length in self.time_series_lengths:
+            if self.multiprocessing:
+                raise NotImplementedError
+            else:
+                coordinates = self.time_serie_length_to_coordinates[length]
+                d[length] = self.generate_all_observation(self.nb_time_series, length, coordinates)
+        return d
+
+    @cached_property
+    def time_serie_length_to_coordinates(self) -> Dict[int, Coordinates]:
+        d = OrderedDict()
+        for length in self.time_series_lengths:
+            d[length] = ConsecutiveTemporalCoordinates.from_nb_temporal_steps(length)
+        return d
+
+    @cached_property
+    def model_class_to_time_serie_length_to_estimator_fitted(self):
+        d = OrderedDict()
+        for model_class in self.models_classes:
+            d_sub = OrderedDict()
+            for time_serie_length, observation_list in self.time_serie_length_to_observation_list.items():
+                coordinates = self.time_serie_length_to_coordinates[time_serie_length]
+                estimators_fitted = []
+                for observations in observation_list:
+                    estimators_fitted.append(self.get_fitted_quantile_estimator(model_class, observations, coordinates))
+                d_sub[time_serie_length] = estimators_fitted
+            d[model_class] = d_sub
+        return d
+
+    def get_fitted_quantile_estimator(self, model_class, observations, coordinates):
+        dataset = AbstractDataset(observations, coordinates)
+        if issubclass(model_class, AbstractTemporalLinearMarginModel):
+            estimator = QuantileEstimatorFromMargin(dataset, self.quantile, model_class)
+        elif issubclass(model_class, AbstractQuantileRegressionModel):
+            estimator = QuantileRegressionEstimator(dataset, self.quantile, model_class)
+        else:
+            raise NotImplementedError
+        estimator.fit()
+        return estimator
+
+    @cached_property
+    def model_class_to_error_last_year_quantile(self):
+        d = OrderedDict()
+        for model_class, d_sub in self.model_class_to_time_serie_length_to_estimator_fitted.items():
+            length_to_error_values = OrderedDict()
+            for length, estimators_fitted in d_sub.items():
+                errors = []
+                # Add the mean, and the quantile of the 95% confidence interval
+                error_values = [0.2, 1, 1.3]
+                length_to_error_values[length] = error_values
+            d[model_class] = length_to_error_values
+        print(d)
+        return d
+
+    def plot_error_for_last_year_quantile(self, show=True):
+        ax = plt.gca()
+        for model_class, length_to_error_values in self.model_class_to_error_last_year_quantile.items():
+            lengths = list(length_to_error_values.keys())
+            errors_values = np.array(list(length_to_error_values.values()))
+            print(errors_values.shape)
+            mean_error = errors_values[:, 1]
+            ax.plot(lengths, mean_error, label=str(model_class))
+            ax.legend()
+        if show:
+            plt.show()
diff --git a/projects/quantile_regression_vs_evt/GevSimulation.py b/projects/quantile_regression_vs_evt/GevSimulation.py
new file mode 100644
index 00000000..317c8e83
--- /dev/null
+++ b/projects/quantile_regression_vs_evt/GevSimulation.py
@@ -0,0 +1,18 @@
+from typing import List
+
+from projects.quantile_regression_vs_evt.AbstractSimulation import AbstractSimulation
+from spatio_temporal_dataset.spatio_temporal_observations.abstract_spatio_temporal_observations import \
+    AbstractSpatioTemporalObservations
+from spatio_temporal_dataset.spatio_temporal_observations.annual_maxima_observations import MarginAnnualMaxima
+
+
+class GevSimulation(AbstractSimulation):
+
+    def __init__(self, margin_model_class_for_simulation, nb_time_series, quantile, time_series_lengths=None, multiprocessing=False, model_classes=None):
+        super().__init__(nb_time_series, quantile, time_series_lengths, multiprocessing, model_classes)
+        self.margin_model_class_for_simulation = margin_model_class_for_simulation
+
+    def generate_all_observation(self, nb_time_series, length, coordinates) -> List[AbstractSpatioTemporalObservations]:
+        margin_model = self.margin_model_class_for_simulation(coordinates)
+        return [MarginAnnualMaxima.from_sampling(nb_obs=length, coordinates=coordinates, margin_model=margin_model) for _ in range(nb_time_series)]
+
diff --git a/test/test_extreme_fit/test_estimator/test_quantile_estimator.py b/test/test_extreme_fit/test_estimator/test_quantile_estimator.py
index ab8a4d60..c3fc11ae 100644
--- a/test/test_extreme_fit/test_estimator/test_quantile_estimator.py
+++ b/test/test_extreme_fit/test_estimator/test_quantile_estimator.py
@@ -24,7 +24,7 @@ class TestQuantileEstimator(unittest.TestCase):
                                                   coordinates=coordinates)
             # Load quantile estimators
             quantile_estimators = [
-                QuantileEstimatorFromMargin(dataset, constant_margin_model, quantile),
+                QuantileEstimatorFromMargin(dataset, quantile, StationaryTemporalModel),
             ]
             for quantile_model_class in load_smooth_quantile_model_classes()[:]:
                 quantile_estimator = QuantileRegressionEstimator(dataset, quantile, quantile_model_class)
diff --git a/test/test_projects/test_quantile_regression/__init__.py b/test/test_projects/test_quantile_regression/__init__.py
new file mode 100644
index 00000000..e69de29b
diff --git a/test/test_projects/test_quantile_regression/test_gev_simulations.py b/test/test_projects/test_quantile_regression/test_gev_simulations.py
new file mode 100644
index 00000000..acddcdea
--- /dev/null
+++ b/test/test_projects/test_quantile_regression/test_gev_simulations.py
@@ -0,0 +1,20 @@
+import unittest
+
+from extreme_fit.model.margin_model.linear_margin_model.temporal_linear_margin_models import StationaryTemporalModel
+from extreme_fit.model.quantile_model.quantile_regression_model import ConstantQuantileRegressionModel
+from projects.quantile_regression_vs_evt.GevSimulation import GevSimulation
+
+
+class TestGevSimulations(unittest.TestCase):
+    DISPLAY = False
+
+    def test_stationary_run(self):
+        margin_model_class_for_simulation = StationaryTemporalModel
+        simulation = GevSimulation(margin_model_class_for_simulation,
+                                   nb_time_series=1, quantile=0.5, time_series_lengths=[50, 60],
+                                   model_classes=[StationaryTemporalModel, ConstantQuantileRegressionModel])
+        simulation.plot_error_for_last_year_quantile(self.DISPLAY)
+
+
+if __name__ == '__main__':
+    unittest.main()
-- 
GitLab