[SCM] add spatio_temporal max stable fit

experiment/meteo_france_SCM_study/visualization/study_visualization/main_study_visualizer.py

@@ -77,7 +77,7 @@ def normal_visualization(temporal_non_stationarity=False):
             study_visualizer = StudyVisualizer(study, save_to_file=save_to_file, temporal_non_stationarity=temporal_non_stationarity)
             # study_visualizer.visualize_independent_margin_fits(threshold=[None, 20, 40, 60][0])
             # study_visualizer.visualize_annual_mean_values()
-            study_visualizer.visualize_linear_margin_fit(only_first_max_stable=None)
+            study_visualizer.visualize_linear_margin_fit(only_first_max_stable=True)
 
 
 def complete_analysis(only_first_one=False):

extreme_estimator/estimator/full_estimator/abstract_full_estimator.py

@@ -56,8 +56,9 @@ class FullEstimatorInASingleStepWithSmoothMargin(AbstractFullEstimator):
     def _fit(self):
         # Estimate both the margin and the max-stable structure
         self._result_from_fit = self.max_stable_model.fitmaxstab(
-            maxima_gev=self.dataset.maxima_gev(split=self.train_split),
-            df_coordinates=self.dataset.df_coordinates(split=self.train_split),
+            data_gev=self.dataset.maxima_gev_for_spatial_extremes_package(self.train_split),
+            df_coordinates_spat=self.dataset.coordinates.df_spatial_coordinates(self.train_split),
+            df_coordinates_temp=self.dataset.coordinates.df_temporal_coordinates(self.train_split),
             fit_marge=True,
             fit_marge_form_dict=self.linear_margin_function_to_fit.form_dict,
             margin_start_dict=self.linear_margin_function_to_fit.coef_dict

extreme_estimator/estimator/margin_estimator/abstract_margin_estimator.py

@@ -38,10 +38,10 @@ class SmoothMarginEstimator(AbstractMarginEstimator):
         self.margin_model = margin_model
 
     def _fit(self):
-        maxima_gev = self.dataset.maxima_gev(self.train_split)
+        maxima_gev_specialized = self.dataset.maxima_gev_for_spatial_extremes_package(self.train_split)
         df_coordinates_spat = self.dataset.coordinates.df_spatial_coordinates(self.train_split)
         df_coordinates_temp = self.dataset.coordinates.df_temporal_coordinates(self.train_split)
-        self._result_from_fit = self.margin_model.fitmargin_from_maxima_gev(maxima_gev=maxima_gev,
+        self._result_from_fit = self.margin_model.fitmargin_from_maxima_gev(data=maxima_gev_specialized,
                                                                             df_coordinates_spat=df_coordinates_spat,
                                                                             df_coordinates_temp=df_coordinates_temp)
         self.extract_fitted_models_from_fitted_params(self.margin_model.margin_function_start_fit, self.fitted_values)

extreme_estimator/estimator/max_stable_estimator/abstract_max_stable_estimator.py

@@ -19,8 +19,8 @@ class MaxStableEstimator(AbstractMaxStableEstimator):
     def _fit(self):
         assert self.dataset.maxima_frech(split=self.train_split) is not None
         self._result_from_fit = self.max_stable_model.fitmaxstab(
-            maxima_frech=self.dataset.maxima_frech(split=self.train_split),
-            df_coordinates=self.dataset.df_coordinates(split=self.train_split))
+            data_frech=self.dataset.maxima_frech_for_spatial_extremes_package(split=self.train_split),
+            df_coordinates_spat=self.dataset.df_coordinates(split=self.train_split))
         self.max_stable_params_fitted = self.fitted_values
 
     def scalars(self, true_max_stable_params: dict):

extreme_estimator/extreme_models/margin_model/parametric_margin_model.py

@@ -20,17 +20,12 @@ class ParametricMarginModel(AbstractMarginModel, ABC):
         self.margin_function_start_fit = None  # type: ParametricMarginFunction
         super().__init__(coordinates, use_start_value, params_start_fit, params_sample)
 
-    def fitmargin_from_maxima_gev(self, maxima_gev: np.ndarray, df_coordinates_spat: pd.DataFrame,
+    def fitmargin_from_maxima_gev(self, data: np.ndarray, df_coordinates_spat: pd.DataFrame,
                                   df_coordinates_temp: pd.DataFrame) -> ResultFromFit:
-        # The reshaping on the line below is only valid if we have a single observation per spatio-temporal point
-        if maxima_gev.shape[1] == 1:
-            maxima_gev = maxima_gev.reshape([len(df_coordinates_spat), len(df_coordinates_temp)])
-        data = np.transpose(maxima_gev)
         assert data.shape[1] == len(df_coordinates_spat)
 
         fit_params = get_margin_formula(self.margin_function_start_fit.form_dict)
 
-
         # Covariables
         covariables = get_coord(df_coordinates=df_coordinates_spat)
         fit_params['temp.cov'] = get_coord(df_coordinates=df_coordinates_temp)

extreme_estimator/extreme_models/max_stable_model/abstract_max_stable_model.py

@@ -20,32 +20,30 @@ class AbstractMaxStableModel(AbstractModel):
     def cov_mod_param(self):
         return {'cov.mod': self.cov_mod}
 
-    def fitmaxstab(self, df_coordinates: pd.DataFrame, maxima_frech: np.ndarray = None, maxima_gev: np.ndarray = None,
+    def fitmaxstab(self, df_coordinates_spat: pd.DataFrame, df_coordinates_temp: pd.DataFrame = None,
+                   data_frech: np.ndarray = None, data_gev: np.ndarray = None,
                    fit_marge=False, fit_marge_form_dict=None, margin_start_dict=None) -> ResultFromFit:
-        assert isinstance(df_coordinates, pd.DataFrame)
+        assert isinstance(df_coordinates_spat, pd.DataFrame)
         if fit_marge:
             assert fit_marge_form_dict is not None
             assert margin_start_dict is not None
 
         # Prepare the data
-        maxima = maxima_gev if fit_marge else maxima_frech
-        assert isinstance(maxima, np.ndarray)
-        assert len(df_coordinates) == len(maxima), 'Coordinates and observations sizes should match,' \
-                                                   'check that the same split was used for both objects, \n' \
-                                                   'df_coordinates size: {}, data size {}'.format(len(df_coordinates),
-                                                                                                  len(maxima))
-        data = np.transpose(maxima)
+        data = data_gev if fit_marge else data_frech
+        assert isinstance(data, np.ndarray)
+        assert len(df_coordinates_spat) == data.shape[1]
 
         # Prepare the coord
-        df_coordinates = df_coordinates.copy()
+        df_coordinates_spat = df_coordinates_spat.copy()
         # In the one dimensional case, fitmaxstab isn't working
         # therefore, we treat our 1D coordinate as 2D coordinate on the line y=x, and enforce iso=TRUE
-        fitmaxstab_with_one_dimensional_data = len(df_coordinates.columns) == 1
+        fitmaxstab_with_one_dimensional_data = len(df_coordinates_spat.columns) == 1
         if fitmaxstab_with_one_dimensional_data:
-            assert AbstractCoordinates.COORDINATE_X in df_coordinates.columns
-            df_coordinates[AbstractCoordinates.COORDINATE_Y] = df_coordinates[AbstractCoordinates.COORDINATE_X]
+            assert AbstractCoordinates.COORDINATE_X in df_coordinates_spat.columns
+            df_coordinates_spat[AbstractCoordinates.COORDINATE_Y] = df_coordinates_spat[
+                AbstractCoordinates.COORDINATE_X]
         # Give names to columns to enable a specification of the shape of each marginal parameter
-        coord = get_coord(df_coordinates)
+        coord = get_coord(df_coordinates_spat)
 
         #  Prepare the fit_params (a dictionary containing all additional parameters)
         fit_params = self.cov_mod_param.copy()
@@ -61,8 +59,16 @@ class AbstractMaxStableModel(AbstractModel):
         fit_params['start'] = r.list(**start_dict)
         fit_params['fit.marge'] = fit_marge
 
+        # Add some temporal covariates
+        # Check the shape of the data
+        has_temp_cov = df_coordinates_temp is not None and len(df_coordinates_temp) > 0
+        if has_temp_cov:
+            assert len(df_coordinates_temp) == len(data)
+            fit_params['temp.cov'] = get_coord(df_coordinates_temp)
+
         # Run the fitmaxstab in R
-        return safe_run_r_estimator(function=r.fitmaxstab, use_start=self.use_start_value, data=data, coord=coord, **fit_params)
+        return safe_run_r_estimator(function=r.fitmaxstab, use_start=self.use_start_value, data=data, coord=coord,
+                                    **fit_params)
 
     def rmaxstab(self, nb_obs: int, coordinates_values: np.ndarray) -> np.ndarray:
         """
@@ -86,7 +92,8 @@ class CovarianceFunction(Enum):
 
 class AbstractMaxStableModelWithCovarianceFunction(AbstractMaxStableModel):
 
-    def __init__(self, use_start_value=False, params_start_fit=None, params_sample=None, covariance_function: CovarianceFunction = None):
+    def __init__(self, use_start_value=False, params_start_fit=None, params_sample=None,
+                 covariance_function: CovarianceFunction = None):
         super().__init__(use_start_value, params_start_fit, params_sample)
         assert covariance_function is not None
         self.covariance_function = covariance_function

spatio_temporal_dataset/coordinates/abstract_coordinates.py

@@ -166,7 +166,7 @@ class AbstractCoordinates(object):
         return len(self.coordinates_spatial_names)
 
     @property
-    def has_spatial_coordinates(self):
+    def has_spatial_coordinates(self) -> bool:
         return self.nb_coordinates_spatial > 0
 
     def df_spatial_coordinates(self, split: Split = Split.all) -> pd.DataFrame:
@@ -189,13 +189,9 @@ class AbstractCoordinates(object):
         return len(self.coordinates_temporal_names)
 
     @property
-    def has_temporal_coordinates(self):
+    def has_temporal_coordinates(self) -> bool:
         return self.nb_coordinates_temporal > 0
 
-    @property
-    def has_spatio_temporal_coordinates(self):
-        return self.has_spatial_coordinates and self.has_temporal_coordinates
-
     def df_temporal_coordinates(self, split: Split = Split.all) -> pd.DataFrame:
         if self.nb_coordinates_temporal == 0:
             return pd.DataFrame()
@@ -206,6 +202,15 @@ class AbstractCoordinates(object):
         df_temporal_coordinates = self.df_temporal_coordinates(split)
         return int(df_temporal_coordinates.min()), int(df_temporal_coordinates.max()),
 
+    # Spatio temporal attributes
+
+    @property
+    def has_spatio_temporal_coordinates(self) -> bool:
+        return self.has_spatial_coordinates and self.has_temporal_coordinates
+
+    def spatio_temporal_shape(self, split: Split.all) -> Tuple[int, int]:
+        return len(self.df_spatial_coordinates(split)), len(self.df_temporal_coordinates(split))
+
     #  Visualization
 
     @property

spatio_temporal_dataset/dataset/abstract_dataset.py

@@ -51,6 +51,20 @@ class AbstractDataset(object):
     def set_maxima_frech(self, maxima_frech_values: np.ndarray, split: Split = Split.all):
         self.observations.set_maxima_frech(maxima_frech_values, split, self.slicer)
 
+    # Observation wrapper for fit function
+
+    def transform_maxima_for_spatial_extreme_package(self, maxima_function, split) -> np.ndarray:
+        array = maxima_function(split)
+        if self.coordinates.has_spatio_temporal_coordinates:
+            array = array.reshape(self.coordinates.spatio_temporal_shape(split))
+        return np.transpose(array)
+
+    def maxima_gev_for_spatial_extremes_package(self, split: Split = Split.all) -> np.ndarray:
+        return self.transform_maxima_for_spatial_extreme_package(self.maxima_gev, split)
+
+    def maxima_frech_for_spatial_extremes_package(self, split: Split = Split.all) -> np.ndarray:
+        return self.transform_maxima_for_spatial_extreme_package(self.maxima_frech, split)
+
     # Coordinates wrapper
 
     def df_coordinates(self, split: Split = Split.all) -> pd.DataFrame: