From 02c210e968ab8f8ebd8bca0415505ef8a539bc16 Mon Sep 17 00:00:00 2001
From: Le Roux Erwan <erwan.le-roux@irstea.fr>
Date: Tue, 28 May 2019 18:40:13 +0200
Subject: [PATCH] [EXTREME ESTIMATOR][MAX STABLE] add argument
 use_rmaxstab_with_2_coordinates to rmaxstab that makes it possible to sample
 dataset with 3D coordinates

---
 .../abstract_max_stable_model.py              |  6 +++-
 .../dataset/simulation_dataset.py             |  6 ++--
 .../annual_maxima_observations.py             |  6 ++--
 .../test_max_stable_model.py                  | 29 +++++++++++++++++++
 .../test_dataset.py                           | 12 ++++++++
 5 files changed, 54 insertions(+), 5 deletions(-)
 create mode 100644 test/test_extreme_estimator/test_extreme_models/test_max_stable_model.py

diff --git a/extreme_estimator/extreme_models/max_stable_model/abstract_max_stable_model.py b/extreme_estimator/extreme_models/max_stable_model/abstract_max_stable_model.py
index 56add58d..029eebc9 100644
--- a/extreme_estimator/extreme_models/max_stable_model/abstract_max_stable_model.py
+++ b/extreme_estimator/extreme_models/max_stable_model/abstract_max_stable_model.py
@@ -71,13 +71,17 @@ class AbstractMaxStableModel(AbstractModel):
                                    **fit_params)
         return ResultFromSpatialExtreme(res)
 
-    def rmaxstab(self, nb_obs: int, coordinates_values: np.ndarray) -> np.ndarray:
+    def rmaxstab(self, nb_obs: int, coordinates_values: np.ndarray, use_rmaxstab_with_2_coordinates=False) -> np.ndarray:
         """
         Return an numpy of maxima. With rows being the stations and columns being the years of maxima
         """
+        if use_rmaxstab_with_2_coordinates and coordinates_values.shape[1] > 2:
+            # When we have more than 2 coordinates, then sample based on the 2 first coordinates only
+            coordinates_values = coordinates_values[:, :2]
 
         maxima_frech = np.array(
             r.rmaxstab(nb_obs, coordinates_values, *list(self.cov_mod_param.values()), **self.params_sample))
+
         return np.transpose(maxima_frech)
 
     def remove_unused_parameters(self, start_dict, coordinate_dim):
diff --git a/spatio_temporal_dataset/dataset/simulation_dataset.py b/spatio_temporal_dataset/dataset/simulation_dataset.py
index c3e90638..ebb35477 100644
--- a/spatio_temporal_dataset/dataset/simulation_dataset.py
+++ b/spatio_temporal_dataset/dataset/simulation_dataset.py
@@ -29,8 +29,10 @@ class SimulatedDataset(AbstractDataset):
 class MaxStableDataset(SimulatedDataset):
 
     @classmethod
-    def from_sampling(cls, nb_obs: int, max_stable_model: AbstractMaxStableModel, coordinates: AbstractCoordinates):
-        observations = MaxStableAnnualMaxima.from_sampling(nb_obs, max_stable_model, coordinates)
+    def from_sampling(cls, nb_obs: int, max_stable_model: AbstractMaxStableModel, coordinates: AbstractCoordinates,
+                      use_rmaxstab_with_2_coordinates=False):
+        observations = MaxStableAnnualMaxima.from_sampling(nb_obs, max_stable_model, coordinates,
+                                                           use_rmaxstab_with_2_coordinates)
         return cls(observations=observations, coordinates=coordinates, max_stable_model=max_stable_model)
 
 
diff --git a/spatio_temporal_dataset/spatio_temporal_observations/annual_maxima_observations.py b/spatio_temporal_dataset/spatio_temporal_observations/annual_maxima_observations.py
index 32817b32..162cdaa4 100644
--- a/spatio_temporal_dataset/spatio_temporal_observations/annual_maxima_observations.py
+++ b/spatio_temporal_dataset/spatio_temporal_observations/annual_maxima_observations.py
@@ -33,8 +33,10 @@ class MarginAnnualMaxima(AnnualMaxima):
 class MaxStableAnnualMaxima(AnnualMaxima):
 
     @classmethod
-    def from_sampling(cls, nb_obs: int, max_stable_model: AbstractMaxStableModel, coordinates: AbstractCoordinates):
-        maxima_frech = max_stable_model.rmaxstab(nb_obs=nb_obs, coordinates_values=coordinates.coordinates_values())
+    def from_sampling(cls, nb_obs: int, max_stable_model: AbstractMaxStableModel, coordinates: AbstractCoordinates,
+                      use_rmaxstab_with_2_coordinates=False):
+        maxima_frech = max_stable_model.rmaxstab(nb_obs=nb_obs, coordinates_values=coordinates.coordinates_values(),
+                                                 use_rmaxstab_with_2_coordinates=use_rmaxstab_with_2_coordinates)
         df_maxima_frech = pd.DataFrame(data=maxima_frech, index=coordinates.index)
         return cls(df_maxima_frech=df_maxima_frech)
 
diff --git a/test/test_extreme_estimator/test_extreme_models/test_max_stable_model.py b/test/test_extreme_estimator/test_extreme_models/test_max_stable_model.py
new file mode 100644
index 00000000..509dc059
--- /dev/null
+++ b/test/test_extreme_estimator/test_extreme_models/test_max_stable_model.py
@@ -0,0 +1,29 @@
+import unittest
+from typing import List
+
+from spatio_temporal_dataset.coordinates.abstract_coordinates import AbstractCoordinates
+from spatio_temporal_dataset.coordinates.transformed_coordinates.transformation.uniform_normalization import \
+    BetweenZeroAndOneNormalization
+from test.test_utils import load_test_3D_spatial_coordinates, load_test_1D_and_2D_spatial_coordinates, \
+    load_test_max_stable_models
+
+
+class TestMaxStableModel(unittest.TestCase):
+
+    def setUp(self) -> None:
+        self.nb_obs = 3
+        self.nb_points = 2
+        self.transformation_class = BetweenZeroAndOneNormalization
+
+    def test_rmaxstab_with_various_coordinates(self):
+        smith_process = load_test_max_stable_models()[0]
+        coordinates = load_test_1D_and_2D_spatial_coordinates(
+            nb_points=self.nb_points)  # type: List[AbstractCoordinates]
+        coordinates += load_test_3D_spatial_coordinates(nb_points=2, transformation_class=self.transformation_class)
+        for coord in coordinates:
+            res = smith_process.rmaxstab(nb_obs=self.nb_obs, coordinates_values=coord.coordinates_values(),
+                                         use_rmaxstab_with_2_coordinates=True)
+            self.assertEqual((self.nb_points, self.nb_obs), res.shape)
+
+if __name__ == '__main__':
+    unittest.main()
diff --git a/test/test_spatio_temporal_dataset/test_dataset.py b/test/test_spatio_temporal_dataset/test_dataset.py
index 357028fa..340c360f 100644
--- a/test/test_spatio_temporal_dataset/test_dataset.py
+++ b/test/test_spatio_temporal_dataset/test_dataset.py
@@ -7,6 +7,8 @@ from rpy2.rinterface import RRuntimeError
 from extreme_estimator.extreme_models.margin_model.linear_margin_model import LinearNonStationaryLocationMarginModel
 from extreme_estimator.extreme_models.utils import set_seed_for_test
 from spatio_temporal_dataset.coordinates.abstract_coordinates import AbstractCoordinates
+from spatio_temporal_dataset.coordinates.transformed_coordinates.transformation.uniform_normalization import \
+    BetweenZeroAndOneNormalization
 from spatio_temporal_dataset.dataset.simulation_dataset import MaxStableDataset, MarginDataset
 from test.test_utils import load_test_max_stable_models, load_test_3D_spatial_coordinates, \
     load_test_1D_and_2D_spatial_coordinates, load_test_spatiotemporal_coordinates
@@ -35,6 +37,16 @@ class TestDataset(unittest.TestCase):
                                            max_stable_model=smith_process,
                                            coordinates=coordinates)
 
+    def test_max_stable_dataset_R3_with_R2_spatial_structure(self):
+        """Test to create a dataset in R3, but where the spatial structure is only with respect to the 2 fisrt coordinates"""
+        coordinates = load_test_3D_spatial_coordinates(nb_points=self.nb_points, transformation_class=BetweenZeroAndOneNormalization)[0]
+        smith_process = load_test_max_stable_models()[0]
+        MaxStableDataset.from_sampling(nb_obs=self.nb_obs,
+                                       max_stable_model=smith_process,
+                                       coordinates=coordinates,
+                                       use_rmaxstab_with_2_coordinates=True)
+
+
 
 class TestSpatioTemporalDataset(unittest.TestCase):
     nb_obs = 2
-- 
GitLab