From f4639557bf56e0a4eff829157eb162d60df2481b Mon Sep 17 00:00:00 2001 From: Le Roux Erwan <erwan.le-roux@irstea.fr> Date: Wed, 29 May 2019 09:40:26 +0200 Subject: [PATCH] [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) --- .../max_stable_model/max_stable_fit.R | 101 +++++++++++++++++- .../abstract_spatial_coordinates.py | 2 - .../test_max_stable_estimators.py | 48 ++++++++- test/test_utils.py | 6 +- 4 files changed, 146 insertions(+), 11 deletions(-) diff --git a/extreme_estimator/extreme_models/max_stable_model/max_stable_fit.R b/extreme_estimator/extreme_models/max_stable_model/max_stable_fit.R index 74195512..9f7635cb 100644 --- a/extreme_estimator/extreme_models/max_stable_model/max_stable_fit.R +++ b/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) diff --git a/spatio_temporal_dataset/coordinates/spatial_coordinates/abstract_spatial_coordinates.py b/spatio_temporal_dataset/coordinates/spatial_coordinates/abstract_spatial_coordinates.py index bb27d6da..f6b3daec 100644 --- a/spatio_temporal_dataset/coordinates/spatial_coordinates/abstract_spatial_coordinates.py +++ b/spatio_temporal_dataset/coordinates/spatial_coordinates/abstract_spatial_coordinates.py @@ -1,8 +1,6 @@ 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 diff --git a/test/test_extreme_estimator/test_estimator/test_max_stable_estimators.py b/test/test_extreme_estimator/test_estimator/test_max_stable_estimators.py index c8efd51c..6c7fa025 100644 --- a/test/test_extreme_estimator/test_estimator/test_max_stable_estimators.py +++ b/test/test_extreme_estimator/test_estimator/test_max_stable_estimators.py @@ -1,8 +1,11 @@ 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__': diff --git a/test/test_utils.py b/test/test_utils.py index 25596f3b..df8d68c3 100644 --- a/test/test_utils.py +++ b/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) -- GitLab