[EXTREME ESTIMATOR][MAX STABLE] add test that check the possibility to fit max...

[EXTREME ESTIMATOR][MAX STABLE] add test that check the possibility to fit max stable with 3D coordinates. This is working with all the max stable structure except Smith structure (see max_stable_fit.R for an example of the returned error)

extreme_estimator/extreme_models/max_stable_model/max_stable_fit.R

@@ -24,11 +24,103 @@ rmaxstab2D <- function (n.obs){
     scale.form = scale ~ 1
     shape.form = shape ~ 1
 
-    namedlist = list(cov11 = 1.0, cov12 = 1.2, cov22 = 2.2, locCoeff1=1.0, locCoeff2=1.0, scaleCoeff1=1.0, shapeCoeff1=1.0)
+    namedlist = list(cov11 = 1.0, cov12 = 1.2, cov22 = 2.2, locCoeff1=1.0, locCoeff2=1.0, scaleCoeff1=2.0, shapeCoeff1=0.1)
     res = fitmaxstab(data=data, coord=coord, cov.mod="gauss", start=namedlist, fit.marge=TRUE, loc.form=loc.form, scale.form=scale.form,shape.form=shape.form)
     print(res['fitted.values'])
 }
 
+# rmaxstab with 3D data
+rmaxstab3Dimprovedgauss <- function (n.obs){
+    n.site = 2
+    coord_sample <- matrix(rnorm(2*n.site, sd = sqrt(.2)), ncol = 2)
+    colnames(coord_sample) = c("E", "N")
+
+    #  Generate the data
+    data <- rmaxstab(n.obs, coord_sample, "gauss", cov11 = 100, cov12 = 25, cov22 = 220)
+    # data <- rmaxstab(n.obs, coord, "brown", range = 3, smooth = 0.5)
+    # data <- rmaxstab(n.obs, coord, "whitmat", nugget = 0.0, range = 3, smooth = 0.5)
+    #  Fit back the data
+    # print(data)n
+    # res = fitmaxstab(data, coord, "gauss", fit.marge=FALSE, )
+    # res = fitmaxstab(data, coord, "brown")
+    # res = fitmaxstab(data, coord, "whitmat", start=)
+
+    coord_fit <- matrix(rnorm(3*n.site, sd = sqrt(.2)), ncol = 3)
+    colnames(coord_fit) = c("E", "N", "A")
+
+    print(class(coord_fit))
+    print(colnames(coord_fit))
+
+    loc.form = loc ~ N
+    scale.form = scale ~ 1
+    shape.form = shape ~ 1
+
+    namedlist = list(cov11 = 1.0, cov12 = 1.2, cov22 = 2.2, cov13 = 1.0, cov23 = 1.2, cov33 = 2.2, locCoeff1=1.0, locCoeff2=1.0, scaleCoeff1=2.0, shapeCoeff1=0.1)
+    res = fitmaxstab(data=data, coord=coord_fit, cov.mod="gauss", start=namedlist, fit.marge=TRUE, loc.form=loc.form, scale.form=scale.form,shape.form=shape.form)
+    print(res['fitted.values'])
+}
+
+rmaxstab3Dimprovedbrown <- function (n.obs){
+    n.site = 2
+    coord_sample <- matrix(rnorm(2*n.site, sd = sqrt(.2)), ncol = 2)
+    colnames(coord_sample) = c("E", "N")
+
+    #  Generate the data
+    data <- rmaxstab(n.obs, coord_sample, "gauss", cov11 = 100, cov12 = 25, cov22 = 220)
+    # data <- rmaxstab(n.obs, coord, "brown", range = 3, smooth = 0.5)
+    # data <- rmaxstab(n.obs, coord, "whitmat", nugget = 0.0, range = 3, smooth = 0.5)
+    #  Fit back the data
+    # print(data)n
+    # res = fitmaxstab(data, coord, "gauss", fit.marge=FALSE, )
+    # res = fitmaxstab(data, coord, "brown")
+    # res = fitmaxstab(data, coord, "whitmat", start=)
+
+    coord_fit <- matrix(rnorm(3*n.site, sd = sqrt(.2)), ncol = 3)
+    colnames(coord_fit) = c("E", "N", "A")
+
+    print(class(coord_fit))
+    print(colnames(coord_fit))
+
+    loc.form = loc ~ N
+    scale.form = scale ~ 1
+    shape.form = shape ~ 1
+
+    # namedlist = list(cov11 = 1.0, cov12 = 1.2, cov22 = 2.2, cov13 = 1.0, cov23 = 1.2, cov33 = 2.2, locCoeff1=1.0, locCoeff2=1.0, scaleCoeff1=2.0, shapeCoeff1=0.1)
+    res = fitmaxstab(data=data, coord=coord_fit, cov.mod="brown", fit.marge=TRUE, loc.form=loc.form, scale.form=scale.form,shape.form=shape.form)
+    print(res['fitted.values'])
+}
+
+rmaxstab3Dimprovedtpowexp <- function (n.obs){
+    n.site = 2
+    coord_sample <- matrix(rnorm(2*n.site, sd = sqrt(.2)), ncol = 2)
+    colnames(coord_sample) = c("E", "N")
+
+    #  Generate the data
+    data <- rmaxstab(n.obs, coord_sample, "gauss", cov11 = 100, cov12 = 25, cov22 = 220)
+    # data <- rmaxstab(n.obs, coord, "brown", range = 3, smooth = 0.5)
+    # data <- rmaxstab(n.obs, coord, "whitmat", nugget = 0.0, range = 3, smooth = 0.5)
+    #  Fit back the data
+    # print(data)n
+    # res = fitmaxstab(data, coord, "gauss", fit.marge=FALSE, )
+    # res = fitmaxstab(data, coord, "brown")
+    # res = fitmaxstab(data, coord, "whitmat", start=)
+
+    coord_fit <- matrix(rnorm(3*n.site, sd = sqrt(.2)), ncol = 3)
+    colnames(coord_fit) = c("E", "N", "A")
+
+    print(class(coord_fit))
+    print(colnames(coord_fit))
+
+    loc.form = loc ~ N
+    scale.form = scale ~ 1
+    shape.form = shape ~ 1
+
+    # namedlist = list(cov11 = 1.0, cov12 = 1.2, cov22 = 2.2, cov13 = 1.0, cov23 = 1.2, cov33 = 2.2, locCoeff1=1.0, locCoeff2=1.0, scaleCoeff1=2.0, shapeCoeff1=0.1)
+    res = fitmaxstab(data=data, coord=coord_fit, cov.mod="tpowexp", fit.marge=TRUE, loc.form=loc.form, scale.form=scale.form,shape.form=shape.form)
+    print(res['fitted.values'])
+}
+
+
 
 # rmaxstab with 3D data
 # rmaxstab3D <- function (n.obs){
@@ -106,10 +198,13 @@ rmaxstab1D <- function (n.obs){
 call_main = !exists("python_wrapping")
 if (call_main) {
     set.seed(42)
-    n.obs = 500
+    n.obs = 50
+    rmaxstab3Dimprovedgauss(n.obs)
+    # rmaxstab3Dimprovedbrown(n.obs)
+    # rmaxstab3Dimprovedtpowexp(n.obs)
     # rmaxstab2D(n.obs)
     # rmaxstab3D(n.obs)
-    rmaxstab1D(n.obs)
+    # rmaxstab1D(n.obs)
 
     # namedlist = list(cov11 = 1.0, cov12 = 1.2, cov22 = 2.2)
     # res = fitmaxstab(data=data, coord=coord, cov.mod="gauss", start=namedlist)

spatio_temporal_dataset/coordinates/spatial_coordinates/abstract_spatial_coordinates.py

 import pandas as pd
 
 from spatio_temporal_dataset.coordinates.abstract_coordinates import AbstractCoordinates
-from spatio_temporal_dataset.coordinates.transformed_coordinates.transformation.abstract_transformation import \
-    AbstractTransformation
 from spatio_temporal_dataset.slicer.spatial_slicer import SpatialSlicer
 
 

test/test_extreme_estimator/test_estimator/test_max_stable_estimators.py

 import unittest
 
+from extreme_estimator.extreme_models.utils import SafeRunException
+from spatio_temporal_dataset.coordinates.transformed_coordinates.transformation.uniform_normalization import \
+    BetweenZeroAndOneNormalization
 from spatio_temporal_dataset.dataset.simulation_dataset import MaxStableDataset
 from test.test_utils import load_test_max_stable_models, load_test_1D_and_2D_spatial_coordinates, \
-    load_test_max_stable_estimators
+    load_test_max_stable_estimators, load_test_3D_spatial_coordinates
 
 
 class TestMaxStableEstimators(unittest.TestCase):
@@ -13,22 +16,59 @@ class TestMaxStableEstimators(unittest.TestCase):
     def setUp(self):
         super().setUp()
         self.coordinates = load_test_1D_and_2D_spatial_coordinates(nb_points=self.nb_points)
+
         self.max_stable_models = load_test_max_stable_models()
 
     def test_max_stable_estimators(self):
+        self.fit_max_stable_estimator_for_all_coordinates()
+        self.assertTrue(True)
+
+    def fit_max_stable_estimator_for_all_coordinates(self):
         for coordinates in self.coordinates:
             for max_stable_model in self.max_stable_models:
+                use_rmaxstab_with_2_coordinates = coordinates.nb_coordinates_spatial > 2
                 dataset = MaxStableDataset.from_sampling(nb_obs=self.nb_obs,
                                                          max_stable_model=max_stable_model,
-                                                         coordinates=coordinates)
-
+                                                         coordinates=coordinates,
+                                                         use_rmaxstab_with_2_coordinates=use_rmaxstab_with_2_coordinates)
                 for max_stable_estimator in load_test_max_stable_estimators(dataset, max_stable_model):
                     max_stable_estimator.fit()
                     if self.DISPLAY:
                         print(type(max_stable_model))
                         print(dataset.df_dataset.head())
                         print(max_stable_estimator.additional_information)
-        self.assertTrue(True)
+
+
+class TestMaxStableEstimatorWorkingFor3DCoordinates(TestMaxStableEstimators):
+
+    def setUp(self):
+        super().setUp()
+        self.coordinates = load_test_3D_spatial_coordinates(nb_points=self.nb_points,
+                                                            transformation_class=BetweenZeroAndOneNormalization)
+        # Select only the max stable structure that work with 3D coordinates
+        self.max_stable_models = load_test_max_stable_models()[1:]
+
+
+class TestMaxStableEstimatorGaussFor3DCoordinates(TestMaxStableEstimators):
+    """
+    See the fhe function rmaxstab3Dimprovedgauss in my file max_stable_fit.R
+    it returns the following error when we try to fit a 3D smipth process:
+
+    Error in nplk(p[1], p[2], p[3], p[4], p[5], p[6], p[7], p[8], p[9], p[10],  :
+    objet 'C_smithdsgnmat3d' introuvable
+    Calls: source ... rmaxstab3Dimprovedgauss -> fitmaxstab -> smithform -> do.call -> nllh -> nplk
+    Exécution arrêtée
+    """
+
+    def setUp(self):
+        super().setUp()
+        self.coordinates = load_test_3D_spatial_coordinates(nb_points=self.nb_points,
+                                                            transformation_class=BetweenZeroAndOneNormalization)
+        self.max_stable_models = load_test_max_stable_models()[:1]
+
+    def test_max_stable_estimators(self):
+        with self.assertRaises(SafeRunException):
+            self.fit_max_stable_estimator_for_all_coordinates()
 
 
 if __name__ == '__main__':

test/test_utils.py

@@ -16,6 +16,8 @@ from extreme_estimator.extreme_models.max_stable_model.max_stable_models import
     Geometric, ExtremalT, ISchlather
 from experiment.meteo_france_data.scm_models_data.safran.safran import SafranSnowfall, Safran, SafranRainfall, \
     SafranTemperature, SafranTotalPrecip
+from spatio_temporal_dataset.coordinates.spatial_coordinates.abstract_spatial_coordinates import \
+    AbstractSpatialCoordinates
 from spatio_temporal_dataset.coordinates.spatial_coordinates.alps_station_3D_coordinates import \
     AlpsStation3DCoordinatesWithAnisotropy
 from spatio_temporal_dataset.coordinates.spatial_coordinates.generated_spatial_coordinates import \
@@ -82,13 +84,13 @@ def load_test_spatial_coordinates(nb_points, coordinate_types, train_split_ratio
             for coordinate_class in coordinate_types]
 
 
-def load_test_1D_and_2D_spatial_coordinates(nb_points, train_split_ratio=None, transformation_class=None):
+def load_test_1D_and_2D_spatial_coordinates(nb_points, train_split_ratio=None, transformation_class=None) -> List[AbstractSpatialCoordinates]:
     return load_test_spatial_coordinates(nb_points, TEST_1D_AND_2D_SPATIAL_COORDINATES,
                                          train_split_ratio=train_split_ratio,
                                          transformation_class=transformation_class)
 
 
-def load_test_3D_spatial_coordinates(nb_points, transformation_class=None):
+def load_test_3D_spatial_coordinates(nb_points, transformation_class=None) -> List[AbstractSpatialCoordinates]:
     return load_test_spatial_coordinates(nb_points, TEST_3D_SPATIAL_COORDINATES, transformation_class=transformation_class)