Commit 4b277185 authored by Le Roux Erwan's avatar Le Roux Erwan
Browse files

[ESTIMATOR] add margin estimator fit with fitspatgev

parent bf835ec0
No related merge requests found
Showing with 199 additions and 139 deletions
+199 -139
...@@ -73,10 +73,12 @@ class SplitCurve(object): ...@@ -73,10 +73,12 @@ class SplitCurve(object):
# Create bins of data, each with an associated color corresponding to its error # Create bins of data, each with an associated color corresponding to its error
data = self.mean_error_dict[gev_value_name].values data = self.mean_error_dict[gev_value_name].values
data_min, data_max = data.min(), data.max()
nb_bins = 10 nb_bins = 10
limits = np.linspace(data.min(), data.max(), num=nb_bins + 1) limits = np.linspace(data_min, data_max, num=nb_bins + 1)
limits[-1] += 0.01 limits[-1] += 0.01
colors = cm.binary(limits) # Binary color should
colors = cm.binary((limits - data_min/ (data_max - data_min)))
# Display train/test points # Display train/test points
for split, marker in [(self.dataset.train_split, 'o'), (self.dataset.test_split, 'x')]: for split, marker in [(self.dataset.train_split, 'o'), (self.dataset.test_split, 'x')]:
......
from experiment.fit_diagnosis.split_curve import SplitCurve
from extreme_estimator.estimator.full_estimator import FullEstimatorInASingleStepWithSmoothMargin
from extreme_estimator.estimator.margin_estimator import SmoothMarginEstimator
from extreme_estimator.extreme_models.margin_model.smooth_margin_model import ConstantMarginModel, \
LinearAllParametersAllDimsMarginModel
from extreme_estimator.extreme_models.max_stable_model.max_stable_models import Smith
from extreme_estimator.gev_params import GevParams
from spatio_temporal_dataset.coordinates.spatial_coordinates.alps_station_2D_coordinates import \
AlpsStation2DCoordinates, AlpsStation2DCoordinatesBetweenZeroAndOne
from spatio_temporal_dataset.coordinates.spatial_coordinates.coordinates_1D import LinSpaceSpatialCoordinates
from spatio_temporal_dataset.dataset.simulation_dataset import FullSimulatedDataset
class SplitCurveAlps(SplitCurve):
def __init__(self, nb_fit: int = 1):
super().__init__(nb_fit)
self.nb_obs = 2
# nb_points = len(AlpsStation2DCoordinates.from_csv())
nb_points = 10
self.coordinates = AlpsStation2DCoordinatesBetweenZeroAndOne.from_nb_points(nb_points, train_split_ratio=0.8)
# MarginModel Linear with respect to the shape (from 0.01 to 0.02)
params_sample = {
(GevParams.GEV_LOC, 0): 10,
(GevParams.GEV_SHAPE, 0): 1.0,
(GevParams.GEV_SCALE, 0): 10,
}
self.margin_model = ConstantMarginModel(coordinates=self.coordinates, params_sample=params_sample)
self.max_stable_model = Smith()
def load_dataset(self):
return FullSimulatedDataset.from_double_sampling(nb_obs=self.nb_obs, margin_model=self.margin_model,
coordinates=self.coordinates,
max_stable_model=self.max_stable_model)
def load_estimator(self, dataset):
max_stable_model = Smith()
margin_model_for_estimator = LinearAllParametersAllDimsMarginModel(dataset.coordinates)
# estimator = FullEstimatorInASingleStepWithSmoothMargin(dataset, margin_model_for_estimator, max_stable_model)
estimator = SmoothMarginEstimator(dataset, margin_model_for_estimator)
return estimator
if __name__ == '__main__':
curve = SplitCurveAlps(nb_fit=2)
curve.fit()
from experiment.fit_diagnosis.split_curve import SplitCurve from experiment.fit_diagnosis.split_curve import SplitCurve
from extreme_estimator.estimator.full_estimator import FullEstimatorInASingleStepWithSmoothMargin from extreme_estimator.estimator.full_estimator import FullEstimatorInASingleStepWithSmoothMargin
from extreme_estimator.estimator.margin_estimator import SmoothMarginEstimator
from extreme_estimator.extreme_models.margin_model.smooth_margin_model import ConstantMarginModel, \ from extreme_estimator.extreme_models.margin_model.smooth_margin_model import ConstantMarginModel, \
LinearAllParametersAllDimsMarginModel LinearAllParametersAllDimsMarginModel
from extreme_estimator.extreme_models.max_stable_model.max_stable_models import Smith from extreme_estimator.extreme_models.max_stable_model.max_stable_models import Smith
...@@ -32,14 +33,11 @@ class SplitCurveExample(SplitCurve): ...@@ -32,14 +33,11 @@ class SplitCurveExample(SplitCurve):
def load_estimator(self, dataset): def load_estimator(self, dataset):
max_stable_model = Smith() max_stable_model = Smith()
margin_model_for_estimator = LinearAllParametersAllDimsMarginModel(dataset.coordinates) margin_model_for_estimator = LinearAllParametersAllDimsMarginModel(dataset.coordinates)
estimator = FullEstimatorInASingleStepWithSmoothMargin(dataset, margin_model_for_estimator, max_stable_model) # estimator = FullEstimatorInASingleStepWithSmoothMargin(dataset, margin_model_for_estimator, max_stable_model)
# estimator = SmoothMarginEstimator(dataset, margin_model_for_estimator) estimator = SmoothMarginEstimator(dataset, margin_model_for_estimator)
return estimator return estimator
if __name__ == '__main__': if __name__ == '__main__':
curve = SplitCurveExample(nb_fit=2) curve = SplitCurveExample(nb_fit=2)
curve.fit() curve.fit()
import time import time
from extreme_estimator.extreme_models.margin_model.margin_function.abstract_margin_function import \
AbstractMarginFunction
from extreme_estimator.extreme_models.margin_model.margin_function.linear_margin_function import LinearMarginFunction
from spatio_temporal_dataset.dataset.abstract_dataset import AbstractDataset from spatio_temporal_dataset.dataset.abstract_dataset import AbstractDataset
from spatio_temporal_dataset.slicer.split import Split from spatio_temporal_dataset.slicer.split import Split
...@@ -15,10 +18,9 @@ class AbstractEstimator(object): ...@@ -15,10 +18,9 @@ class AbstractEstimator(object):
def __init__(self, dataset: AbstractDataset): def __init__(self, dataset: AbstractDataset):
self.dataset = dataset # type: AbstractDataset self.dataset = dataset # type: AbstractDataset
self.additional_information = dict() self.additional_information = dict()
self._params_fitted = None
@property self._margin_function_fitted = None
def train_split(self): self._max_stable_model_fitted = None
return self.dataset.train_split
def fit(self): def fit(self):
ts = time.time() ts = time.time()
...@@ -26,6 +28,31 @@ class AbstractEstimator(object): ...@@ -26,6 +28,31 @@ class AbstractEstimator(object):
te = time.time() te = time.time()
self.additional_information[self.DURATION] = int((te - ts) * 1000) self.additional_information[self.DURATION] = int((te - ts) * 1000)
@property
def params_fitted(self):
assert self.is_fitted
return self._params_fitted
@property
def margin_function_fitted(self) -> AbstractMarginFunction:
assert self.is_fitted
assert self._margin_function_fitted is not None, 'No margin function has been fitted'
return self._margin_function_fitted
def extract_fitted_models_from_fitted_params(self, margin_function_to_fit, full_params_fitted):
coef_dict = {k: v for k, v in full_params_fitted.items() if 'Coeff' in k}
self._margin_function_fitted = LinearMarginFunction.from_coef_dict(coordinates=self.dataset.coordinates,
gev_param_name_to_linear_dims=margin_function_to_fit.gev_param_name_to_linear_dims,
coef_dict=coef_dict)
@property
def is_fitted(self):
return self._params_fitted is not None
@property
def train_split(self):
return self.dataset.train_split
def scalars(self, true_max_stable_params: dict): def scalars(self, true_max_stable_params: dict):
error = self._error(true_max_stable_params) error = self._error(true_max_stable_params)
return {**error, **self.additional_information} return {**error, **self.additional_information}
......
...@@ -14,18 +14,6 @@ class AbstractFullEstimator(AbstractEstimator): ...@@ -14,18 +14,6 @@ class AbstractFullEstimator(AbstractEstimator):
def __init__(self, dataset: AbstractDataset): def __init__(self, dataset: AbstractDataset):
super().__init__(dataset) super().__init__(dataset)
self._margin_function_fitted = None
self._max_stable_model_fitted = None
@property
def margin_function_fitted(self) -> AbstractMarginFunction:
assert self._margin_function_fitted is not None, 'Error: estimator has not been fitted'
return self._margin_function_fitted
# @property
# def max_stable_fitted(self) -> AbstractMarginFunction:
# assert self._margin_function_fitted is not None, 'Error: estimator has not been fitted'
# return self._margin_function_fitted
class SmoothMarginalsThenUnitaryMsp(AbstractFullEstimator): class SmoothMarginalsThenUnitaryMsp(AbstractFullEstimator):
...@@ -43,7 +31,8 @@ class SmoothMarginalsThenUnitaryMsp(AbstractFullEstimator): ...@@ -43,7 +31,8 @@ class SmoothMarginalsThenUnitaryMsp(AbstractFullEstimator):
# Compute the maxima_frech # Compute the maxima_frech
maxima_gev_train = self.dataset.maxima_gev(split=self.train_split) maxima_gev_train = self.dataset.maxima_gev(split=self.train_split)
maxima_frech = AbstractMarginModel.gev2frech(maxima_gev=maxima_gev_train, maxima_frech = AbstractMarginModel.gev2frech(maxima_gev=maxima_gev_train,
coordinates_values=self.dataset.coordinates_values(self.train_split), coordinates_values=self.dataset.coordinates_values(
self.train_split),
margin_function=self.margin_estimator.margin_function_fitted) margin_function=self.margin_estimator.margin_function_fitted)
# Update maxima frech field through the dataset object # Update maxima frech field through the dataset object
self.dataset.set_maxima_frech(maxima_frech, split=self.train_split) self.dataset.set_maxima_frech(maxima_frech, split=self.train_split)
...@@ -68,7 +57,7 @@ class FullEstimatorInASingleStepWithSmoothMargin(AbstractFullEstimator): ...@@ -68,7 +57,7 @@ class FullEstimatorInASingleStepWithSmoothMargin(AbstractFullEstimator):
def _fit(self): def _fit(self):
# Estimate both the margin and the max-stable structure # Estimate both the margin and the max-stable structure
self.full_params_fitted = self.max_stable_model.fitmaxstab( self._params_fitted = self.max_stable_model.fitmaxstab(
maxima_gev=self.dataset.maxima_gev(split=self.train_split), maxima_gev=self.dataset.maxima_gev(split=self.train_split),
df_coordinates=self.dataset.df_coordinates(split=self.train_split), df_coordinates=self.dataset.df_coordinates(split=self.train_split),
fit_marge=True, fit_marge=True,
...@@ -76,10 +65,7 @@ class FullEstimatorInASingleStepWithSmoothMargin(AbstractFullEstimator): ...@@ -76,10 +65,7 @@ class FullEstimatorInASingleStepWithSmoothMargin(AbstractFullEstimator):
margin_start_dict=self.linear_margin_function_to_fit.coef_dict margin_start_dict=self.linear_margin_function_to_fit.coef_dict
) )
# Create the fitted margin function # Create the fitted margin function
coef_dict = {k: v for k, v in self.full_params_fitted.items() if 'Coeff' in k} self.extract_fitted_models_from_fitted_params(self.linear_margin_function_to_fit, self._params_fitted)
self._margin_function_fitted = LinearMarginFunction.from_coef_dict(coordinates=self.dataset.coordinates,
gev_param_name_to_linear_dims=self.linear_margin_function_to_fit.gev_param_name_to_linear_dims,
coef_dict=coef_dict)
class PointwiseAndThenUnitaryMsp(AbstractFullEstimator): class PointwiseAndThenUnitaryMsp(AbstractFullEstimator):
......
...@@ -43,6 +43,8 @@ class SmoothMarginEstimator(AbstractMarginEstimator): ...@@ -43,6 +43,8 @@ class SmoothMarginEstimator(AbstractMarginEstimator):
def _fit(self): def _fit(self):
maxima_gev = self.dataset.maxima_gev(split=self.train_split) maxima_gev = self.dataset.maxima_gev(split=self.train_split)
coordinate_values = self.dataset.coordinates_values(self.train_split) df_coordinates = self.dataset.df_coordinates(self.train_split)
self._margin_function_fitted = self.margin_model.fitmargin_from_maxima_gev(maxima_gev=maxima_gev, self._params_fitted = self.margin_model.fitmargin_from_maxima_gev(maxima_gev=maxima_gev,
coordinates_values=coordinate_values) df_coordinates=df_coordinates)
self.extract_fitted_models_from_fitted_params(self.margin_model.margin_function_start_fit, self._params_fitted)
assert isinstance(self.margin_function_fitted, AbstractMarginFunction)
class AbstractModel(object): class AbstractModel(object):
def __init__(self, params_start_fit=None, params_sample=None): def __init__(self, use_start_value=True, params_start_fit=None, params_sample=None):
self.default_params_start_fit = None self.default_params_start_fit = None
self.default_params_sample = None self.default_params_sample = None
self.use_start_value = use_start_value
self.user_params_start_fit = params_start_fit self.user_params_start_fit = params_start_fit
self.user_params_sample = params_sample self.user_params_sample = params_sample
......
import numpy as np import numpy as np
import pandas as pd
from extreme_estimator.extreme_models.abstract_model import AbstractModel from extreme_estimator.extreme_models.abstract_model import AbstractModel
from extreme_estimator.extreme_models.margin_model.margin_function.abstract_margin_function \ from extreme_estimator.extreme_models.margin_model.margin_function.abstract_margin_function \
...@@ -10,8 +11,8 @@ from spatio_temporal_dataset.coordinates.abstract_coordinates import AbstractCoo ...@@ -10,8 +11,8 @@ from spatio_temporal_dataset.coordinates.abstract_coordinates import AbstractCoo
class AbstractMarginModel(AbstractModel): class AbstractMarginModel(AbstractModel):
def __init__(self, coordinates: AbstractCoordinates, params_start_fit=None, params_sample=None): def __init__(self, coordinates: AbstractCoordinates, use_start_value=True, params_start_fit=None, params_sample=None):
super().__init__(params_start_fit, params_sample) super().__init__(use_start_value, params_start_fit, params_sample)
assert isinstance(coordinates, AbstractCoordinates), type(coordinates) assert isinstance(coordinates, AbstractCoordinates), type(coordinates)
self.coordinates = coordinates self.coordinates = coordinates
self.margin_function_sample = None # type: AbstractMarginFunction self.margin_function_sample = None # type: AbstractMarginFunction
...@@ -68,7 +69,7 @@ class AbstractMarginModel(AbstractModel): ...@@ -68,7 +69,7 @@ class AbstractMarginModel(AbstractModel):
# Fitting methods needs to be defined in child classes # Fitting methods needs to be defined in child classes
def fitmargin_from_maxima_gev(self, maxima_gev: np.ndarray, coordinates_values: np.ndarray) \ def fitmargin_from_maxima_gev(self, maxima_gev: np.ndarray, df_coordinates: pd.DataFrame) \
-> AbstractMarginFunction: -> AbstractMarginFunction:
pass pass
......
...@@ -24,6 +24,9 @@ class AbstractMarginFunction(object): ...@@ -24,6 +24,9 @@ class AbstractMarginFunction(object):
self.color = 'skyblue' self.color = 'skyblue'
self.filter = None self.filter = None
self._grid_2D = None
self._grid_1D = None
def get_gev_params(self, coordinate: np.ndarray) -> GevParams: def get_gev_params(self, coordinate: np.ndarray) -> GevParams:
"""Main method that maps each coordinate to its GEV parameters""" """Main method that maps each coordinate to its GEV parameters"""
pass pass
...@@ -76,7 +79,7 @@ class AbstractMarginFunction(object): ...@@ -76,7 +79,7 @@ class AbstractMarginFunction(object):
def visualize_1D(self, gev_value_name=GevParams.GEV_LOC, ax=None, show=True): def visualize_1D(self, gev_value_name=GevParams.GEV_LOC, ax=None, show=True):
x = self.coordinates.x_coordinates x = self.coordinates.x_coordinates
grid, linspace = self.get_grid_values_1D(x) grid, linspace = self.grid_1D(x)
if ax is None: if ax is None:
ax = plt.gca() ax = plt.gca()
if self.datapoint_display: if self.datapoint_display:
...@@ -91,6 +94,11 @@ class AbstractMarginFunction(object): ...@@ -91,6 +94,11 @@ class AbstractMarginFunction(object):
if show: if show:
plt.show() plt.show()
def grid_1D(self, x):
if self._grid_1D is None:
self._grid_1D = self.get_grid_values_1D(x)
return self._grid_1D
def get_grid_values_1D(self, x): def get_grid_values_1D(self, x):
# TODO: to avoid getting the value several times, I could cache the results # TODO: to avoid getting the value several times, I could cache the results
if self.datapoint_display: if self.datapoint_display:
...@@ -115,7 +123,7 @@ class AbstractMarginFunction(object): ...@@ -115,7 +123,7 @@ class AbstractMarginFunction(object):
def visualize_2D(self, gev_value_name=GevParams.GEV_LOC, ax=None, show=True): def visualize_2D(self, gev_value_name=GevParams.GEV_LOC, ax=None, show=True):
x = self.coordinates.x_coordinates x = self.coordinates.x_coordinates
y = self.coordinates.y_coordinates y = self.coordinates.y_coordinates
grid = self.get_grid_2D(x, y) grid = self.grid_2D(x, y)
if ax is None: if ax is None:
ax = plt.gca() ax = plt.gca()
imshow_method = ax.imshow imshow_method = ax.imshow
...@@ -133,6 +141,11 @@ class AbstractMarginFunction(object): ...@@ -133,6 +141,11 @@ class AbstractMarginFunction(object):
if show: if show:
plt.show() plt.show()
def grid_2D(self, x, y):
if self._grid_2D is None:
self._grid_2D = self.get_grid_2D(x, y)
return self._grid_2D
def get_grid_2D(self, x, y): def get_grid_2D(self, x, y):
resolution = 100 resolution = 100
grid = [] grid = []
......
...@@ -22,6 +22,6 @@ def error_dict_between_margin_functions(reference: AbstractMarginFunction, fitte ...@@ -22,6 +22,6 @@ def error_dict_between_margin_functions(reference: AbstractMarginFunction, fitte
fitted_values = fitted.gev_value_name_to_serie fitted_values = fitted.gev_value_name_to_serie
gev_param_name_to_error_serie = {} gev_param_name_to_error_serie = {}
for gev_value_name in GevParams.GEV_VALUE_NAMES: for gev_value_name in GevParams.GEV_VALUE_NAMES:
error = relative_abs_error(reference_values[gev_value_name], fitted_values[gev_value_name]) error = 100 * relative_abs_error(reference_values[gev_value_name], fitted_values[gev_value_name])
gev_param_name_to_error_serie[gev_value_name] = error gev_param_name_to_error_serie[gev_value_name] = error
return gev_param_name_to_error_serie return gev_param_name_to_error_serie
from typing import Dict
import numpy as np import numpy as np
import pandas as pd
from extreme_estimator.extreme_models.margin_model.margin_function.abstract_margin_function import \ from extreme_estimator.extreme_models.margin_model.margin_function.abstract_margin_function import \
AbstractMarginFunction AbstractMarginFunction
from extreme_estimator.extreme_models.margin_model.abstract_margin_model import AbstractMarginModel from extreme_estimator.extreme_models.margin_model.abstract_margin_model import AbstractMarginModel
from extreme_estimator.extreme_models.margin_model.margin_function.linear_margin_function import LinearMarginFunction from extreme_estimator.extreme_models.margin_model.margin_function.linear_margin_function import LinearMarginFunction
from extreme_estimator.extreme_models.margin_model.param_function.linear_coef import LinearCoef from extreme_estimator.extreme_models.margin_model.param_function.linear_coef import LinearCoef
from extreme_estimator.extreme_models.utils import safe_run_r_estimator, r, retrieve_fitted_values, get_coord, \
get_margin_formula
from extreme_estimator.gev_params import GevParams from extreme_estimator.gev_params import GevParams
...@@ -47,10 +52,6 @@ class LinearMarginModel(AbstractMarginModel): ...@@ -47,10 +52,6 @@ class LinearMarginModel(AbstractMarginModel):
gev_param_name_to_linear_coef[gev_param_name] = linear_coef gev_param_name_to_linear_coef[gev_param_name] = linear_coef
return gev_param_name_to_linear_coef return gev_param_name_to_linear_coef
def fitmargin_from_maxima_gev(self, maxima_gev: np.ndarray,
coordinates_values: np.ndarray) -> AbstractMarginFunction:
return self.margin_function_start_fit
@classmethod @classmethod
def from_coef_list(cls, coordinates, gev_param_name_to_coef_list): def from_coef_list(cls, coordinates, gev_param_name_to_coef_list):
params = {} params = {}
...@@ -59,6 +60,16 @@ class LinearMarginModel(AbstractMarginModel): ...@@ -59,6 +60,16 @@ class LinearMarginModel(AbstractMarginModel):
params[(gev_param_name, dim)] = coef params[(gev_param_name, dim)] = coef
return cls(coordinates, params_sample=params, params_start_fit=params) return cls(coordinates, params_sample=params, params_start_fit=params)
def fitmargin_from_maxima_gev(self, maxima_gev: np.ndarray,
df_coordinates: pd.DataFrame) -> Dict[str, float]:
data = np.transpose(maxima_gev)
covariables = get_coord(df_coordinates)
fit_params = get_margin_formula(self.margin_function_start_fit.form_dict)
if self.use_start_value:
fit_params['start'] = r.list(**self.margin_function_start_fit.coef_dict)
res = safe_run_r_estimator(function=r.fitspatgev, data=data, covariables=covariables, **fit_params)
return retrieve_fitted_values(res)
class ConstantMarginModel(LinearMarginModel): class ConstantMarginModel(LinearMarginModel):
......
import pandas as pd
from enum import Enum from enum import Enum
import numpy as np import numpy as np
import rpy2 import pandas as pd
from rpy2.rinterface import RRuntimeError, RRuntimeWarning
import rpy2.robjects as robjects import rpy2.robjects as robjects
from extreme_estimator.extreme_models.abstract_model import AbstractModel from extreme_estimator.extreme_models.abstract_model import AbstractModel
from extreme_estimator.extreme_models.utils import r, set_seed_r from extreme_estimator.extreme_models.utils import r, safe_run_r_estimator, retrieve_fitted_values, get_coord, \
get_margin_formula
from spatio_temporal_dataset.coordinates.abstract_coordinates import AbstractCoordinates from spatio_temporal_dataset.coordinates.abstract_coordinates import AbstractCoordinates
class AbstractMaxStableModel(AbstractModel): class AbstractMaxStableModel(AbstractModel):
def __init__(self, params_start_fit=None, params_sample=None): def __init__(self, use_start_value=True, params_start_fit=None, params_sample=None):
super().__init__(params_start_fit, params_sample) super().__init__(use_start_value, params_start_fit, params_sample)
self.cov_mod = None self.cov_mod = None
self.start_value_for_fitmaxstab = True
@property @property
def cov_mod_param(self): def cov_mod_param(self):
...@@ -48,8 +46,7 @@ class AbstractMaxStableModel(AbstractModel): ...@@ -48,8 +46,7 @@ class AbstractMaxStableModel(AbstractModel):
assert AbstractCoordinates.COORDINATE_X in df_coordinates.columns assert AbstractCoordinates.COORDINATE_X in df_coordinates.columns
df_coordinates[AbstractCoordinates.COORDINATE_Y] = df_coordinates[AbstractCoordinates.COORDINATE_X] df_coordinates[AbstractCoordinates.COORDINATE_Y] = df_coordinates[AbstractCoordinates.COORDINATE_X]
# Give names to columns to enable a specification of the shape of each marginal parameter # Give names to columns to enable a specification of the shape of each marginal parameter
coord = robjects.vectors.Matrix(df_coordinates.values) coord = get_coord(df_coordinates)
coord.colnames = robjects.StrVector(list(df_coordinates.columns))
# Prepare the fit_params (a dictionary containing all additional parameters) # Prepare the fit_params (a dictionary containing all additional parameters)
fit_params = self.cov_mod_param.copy() fit_params = self.cov_mod_param.copy()
...@@ -58,27 +55,17 @@ class AbstractMaxStableModel(AbstractModel): ...@@ -58,27 +55,17 @@ class AbstractMaxStableModel(AbstractModel):
start_dict = self.remove_unused_parameters(start_dict, fitmaxstab_with_one_dimensional_data) start_dict = self.remove_unused_parameters(start_dict, fitmaxstab_with_one_dimensional_data)
if fit_marge: if fit_marge:
start_dict.update(margin_start_dict) start_dict.update(margin_start_dict)
fit_params.update({k: robjects.Formula(v) for k, v in fit_marge_form_dict.items()}) margin_formulas = get_margin_formula(fit_marge_form_dict)
fit_params.update(margin_formulas)
if fitmaxstab_with_one_dimensional_data: if fitmaxstab_with_one_dimensional_data:
fit_params['iso'] = True fit_params['iso'] = True
if self.start_value_for_fitmaxstab: if self.use_start_value:
fit_params['start'] = r.list(**start_dict) fit_params['start'] = r.list(**start_dict)
fit_params['fit.marge'] = fit_marge fit_params['fit.marge'] = fit_marge
# Run the fitmaxstab in R # Run the fitmaxstab in R
try: res = safe_run_r_estimator(function=r.fitmaxstab, data=data, coord=coord, **fit_params)
res = r.fitmaxstab(data=data, coord=coord, **fit_params) # type: robjects.ListVector return retrieve_fitted_values(res)
except (RRuntimeError, RRuntimeWarning) as e:
if isinstance(e, RRuntimeError):
raise Exception('Some R exception have been launched at RunTime: \n {}'.format(e.__repr__()))
if isinstance(e, RRuntimeWarning):
print(e.__repr__())
print('WARNING')
# todo: maybe if the convergence was not successful I could try other starting point several times
# Retrieve the resulting fitted values
fitted_values = res.rx2('fitted.values')
fitted_values = {key: fitted_values.rx2(key)[0] for key in fitted_values.names}
return fitted_values
def rmaxstab(self, nb_obs: int, coordinates_values: np.ndarray) -> np.ndarray: def rmaxstab(self, nb_obs: int, coordinates_values: np.ndarray) -> np.ndarray:
""" """
...@@ -102,8 +89,8 @@ class CovarianceFunction(Enum): ...@@ -102,8 +89,8 @@ class CovarianceFunction(Enum):
class AbstractMaxStableModelWithCovarianceFunction(AbstractMaxStableModel): class AbstractMaxStableModelWithCovarianceFunction(AbstractMaxStableModel):
def __init__(self, params_start_fit=None, params_sample=None, covariance_function: CovarianceFunction = None): def __init__(self, use_start_value=True, params_start_fit=None, params_sample=None, covariance_function: CovarianceFunction = None):
super().__init__(params_start_fit, params_sample) super().__init__(use_start_value, params_start_fit, params_sample)
assert covariance_function is not None assert covariance_function is not None
self.covariance_function = covariance_function self.covariance_function = covariance_function
self.default_params_sample = { self.default_params_sample = {
......
# test how to call the max stable method
import pandas as pd
import numpy as np
from extreme_estimator.extreme_models.utils import R
from extreme_estimator.gev.gev_mle_fit import GevMleFit
import rpy2.robjects.numpy2ri as rpyn
import rpy2.robjects as robjects
def max_stable_fit():
r = R().r
r("""
set.seed(42)
n.obs = 50
n.site = 2
coord <- matrix(rnorm(2*n.site, sd = sqrt(.2)), ncol = 2)
colnames(coord) = c("E", "N")
# Generate the data
data <- rmaxstab(n.obs, coord, "gauss", cov11 = 100, cov12 = 25, cov22 = 220)
loc.form = loc ~ 1
scale.form = scale ~ 1
shape.form = shape ~ 1
namedlist = list(cov11 = 1.0, cov12 = 1.2, cov22 = 2.2, locCoeff1=1.0, scaleCoeff1=1.0, shapeCoeff1=1.0)
# 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)
""")
# coord = np.array(r['coord'])
data = r['data']
coord = pd.DataFrame(r['coord'])
coord.colnames = robjects.StrVector(['E', 'N'])
print(r['loc.form'])
print(type(r['loc.form']))
# x_gev = rpyn.numpy2ri(x_gev)
print(coord)
print(coord.colnames)
# res = r.fitmaxstab(data=data, coord=coord, covmod="gauss", start=namedlist, fit.marge=TRUE, loc.form=loc.form, scale.form=scale.form,shape.form=shape.form)
# m2.colnames = R.StrVector(['x', 'y'])
# res = r.fitmaxstab()
# mle_params = dict(r.attr(res, 'coef').items())
# print(mle_params)
if __name__ == '__main__':
max_stable_fit()
...@@ -6,8 +6,8 @@ from extreme_estimator.extreme_models.max_stable_model.abstract_max_stable_model ...@@ -6,8 +6,8 @@ from extreme_estimator.extreme_models.max_stable_model.abstract_max_stable_model
class Smith(AbstractMaxStableModel): class Smith(AbstractMaxStableModel):
def __init__(self, params_start_fit=None, params_sample=None): def __init__(self, *args, **kwargs):
super().__init__(params_start_fit=params_start_fit, params_sample=params_sample) super().__init__(*args, **kwargs)
self.cov_mod = 'gauss' self.cov_mod = 'gauss'
self.default_params_start_fit = { self.default_params_start_fit = {
'var': 1, 'var': 1,
...@@ -27,8 +27,8 @@ class Smith(AbstractMaxStableModel): ...@@ -27,8 +27,8 @@ class Smith(AbstractMaxStableModel):
class BrownResnick(AbstractMaxStableModel): class BrownResnick(AbstractMaxStableModel):
def __init__(self, params_start_fit=None, params_sample=None): def __init__(self, *args, **kwargs):
super().__init__(params_start_fit=params_start_fit, params_sample=params_sample) super().__init__(*args, **kwargs)
self.cov_mod = 'brown' self.cov_mod = 'brown'
self.default_params_start_fit = { self.default_params_start_fit = {
'range': 3, 'range': 3,
...@@ -42,8 +42,8 @@ class BrownResnick(AbstractMaxStableModel): ...@@ -42,8 +42,8 @@ class BrownResnick(AbstractMaxStableModel):
class Schlather(AbstractMaxStableModelWithCovarianceFunction): class Schlather(AbstractMaxStableModelWithCovarianceFunction):
def __init__(self, params_start_fit=None, params_sample=None, covariance_function: CovarianceFunction = None): def __init__(self, *args, **kwargs):
super().__init__(params_start_fit, params_sample, covariance_function) super().__init__(*args, **kwargs)
self.cov_mod = self.covariance_function.name self.cov_mod = self.covariance_function.name
self.default_params_sample.update({}) self.default_params_sample.update({})
self.default_params_start_fit = self.default_params_sample.copy() self.default_params_start_fit = self.default_params_sample.copy()
...@@ -51,8 +51,8 @@ class Schlather(AbstractMaxStableModelWithCovarianceFunction): ...@@ -51,8 +51,8 @@ class Schlather(AbstractMaxStableModelWithCovarianceFunction):
class Geometric(AbstractMaxStableModelWithCovarianceFunction): class Geometric(AbstractMaxStableModelWithCovarianceFunction):
def __init__(self, params_start_fit=None, params_sample=None, covariance_function: CovarianceFunction = None): def __init__(self, *args, **kwargs):
super().__init__(params_start_fit, params_sample, covariance_function) super().__init__(*args, **kwargs)
self.cov_mod = 'g' + self.covariance_function.name self.cov_mod = 'g' + self.covariance_function.name
self.default_params_sample.update({'sigma2': 0.5}) self.default_params_sample.update({'sigma2': 0.5})
self.default_params_start_fit = self.default_params_sample.copy() self.default_params_start_fit = self.default_params_sample.copy()
...@@ -60,8 +60,8 @@ class Geometric(AbstractMaxStableModelWithCovarianceFunction): ...@@ -60,8 +60,8 @@ class Geometric(AbstractMaxStableModelWithCovarianceFunction):
class ExtremalT(AbstractMaxStableModelWithCovarianceFunction): class ExtremalT(AbstractMaxStableModelWithCovarianceFunction):
def __init__(self, params_start_fit=None, params_sample=None, covariance_function: CovarianceFunction = None): def __init__(self, *args, **kwargs):
super().__init__(params_start_fit, params_sample, covariance_function) super().__init__(*args, **kwargs)
self.cov_mod = 't' + self.covariance_function.name self.cov_mod = 't' + self.covariance_function.name
self.default_params_sample.update({'DoF': 2}) self.default_params_sample.update({'DoF': 2})
self.default_params_start_fit = self.default_params_sample.copy() self.default_params_start_fit = self.default_params_sample.copy()
...@@ -69,8 +69,8 @@ class ExtremalT(AbstractMaxStableModelWithCovarianceFunction): ...@@ -69,8 +69,8 @@ class ExtremalT(AbstractMaxStableModelWithCovarianceFunction):
class ISchlather(AbstractMaxStableModelWithCovarianceFunction): class ISchlather(AbstractMaxStableModelWithCovarianceFunction):
def __init__(self, params_start_fit=None, params_sample=None, covariance_function: CovarianceFunction = None): def __init__(self, *args, **kwargs):
super().__init__(params_start_fit, params_sample, covariance_function) super().__init__(*args, **kwargs)
self.cov_mod = 'i' + self.covariance_function.name self.cov_mod = 'i' + self.covariance_function.name
self.default_params_sample.update({'alpha': 0.5}) self.default_params_sample.update({'alpha': 0.5})
self.default_params_start_fit = self.default_params_sample.copy() self.default_params_start_fit = self.default_params_sample.copy()
import os.path as op import os.path as op
import random import random
import sys import sys
from typing import Dict
import pandas as pd
import rpy2.robjects as ro import rpy2.robjects as ro
from rpy2 import robjects
from rpy2.rinterface import RRuntimeWarning
from rpy2.rinterface._rinterface import RRuntimeError
from rpy2.robjects import numpy2ri from rpy2.robjects import numpy2ri
from rpy2.robjects import pandas2ri from rpy2.robjects import pandas2ri
...@@ -23,6 +28,36 @@ def get_associated_r_file(python_filepath: str) -> str: ...@@ -23,6 +28,36 @@ def get_associated_r_file(python_filepath: str) -> str:
assert op.exists(r_filepath), r_filepath assert op.exists(r_filepath), r_filepath
return r_filepath return r_filepath
def safe_run_r_estimator(function, **parameters):
try:
res = function(**parameters) # type:
except (RRuntimeError, RRuntimeWarning) as e:
if isinstance(e, RRuntimeError):
raise Exception('Some R exception have been launched at RunTime: \n {}'.format(e.__repr__()))
if isinstance(e, RRuntimeWarning):
print(e.__repr__())
print('WARNING')
return res
def retrieve_fitted_values(res: robjects.ListVector):
# todo: maybe if the convergence was not successful I could try other starting point several times
# Retrieve the resulting fitted values
fitted_values = res.rx2('fitted.values')
fitted_values = {key: fitted_values.rx2(key)[0] for key in fitted_values.names}
return fitted_values
def get_coord(df_coordinates: pd.DataFrame):
coord = robjects.vectors.Matrix(df_coordinates.values)
coord.colnames = robjects.StrVector(list(df_coordinates.columns))
return coord
def get_margin_formula(fit_marge_form_dict) -> Dict:
return {k: robjects.Formula(v) for k, v in fit_marge_form_dict.items()}
# def conversion_to_FloatVector(data): # def conversion_to_FloatVector(data):
# """Convert DataFrame or numpy array into FloatVector for r""" # """Convert DataFrame or numpy array into FloatVector for r"""
# if isinstance(data, pd.DataFrame): # if isinstance(data, pd.DataFrame):
......
File moved
File moved
...@@ -18,6 +18,11 @@ class AbstractDataset(object): ...@@ -18,6 +18,11 @@ class AbstractDataset(object):
self.observations = observations self.observations = observations
self.coordinates = coordinates self.coordinates = coordinates
# @property
# def max_stable_fitted(self) -> AbstractMarginFunction:
# assert self._margin_function_fitted is not None, 'Error: estimator has not been fitted'
# return self._margin_function_fitted
@property @property
def slicer(self): def slicer(self):
return self.coordinates.slicer return self.coordinates.slicer
...@@ -42,9 +47,6 @@ class AbstractDataset(object): ...@@ -42,9 +47,6 @@ class AbstractDataset(object):
# todo: maybe I should add the split from the temporal observations # todo: maybe I should add the split from the temporal observations
return self.observations.df_maxima_gev.join(self.coordinates.df_merged) return self.observations.df_maxima_gev.join(self.coordinates.df_merged)
def df_coordinates(self, split: Split = Split.all) -> pd.DataFrame:
return self.coordinates.df_coordinates(split=split)
# Observation wrapper # Observation wrapper
def maxima_gev(self, split: Split = Split.all) -> np.ndarray: def maxima_gev(self, split: Split = Split.all) -> np.ndarray:
...@@ -58,6 +60,9 @@ class AbstractDataset(object): ...@@ -58,6 +60,9 @@ class AbstractDataset(object):
# Coordinates wrapper # Coordinates wrapper
def df_coordinates(self, split: Split = Split.all) -> pd.DataFrame:
return self.coordinates.df_coordinates(split=split)
def coordinates_values(self, split: Split = Split.all) -> np.ndarray: def coordinates_values(self, split: Split = Split.all) -> np.ndarray:
return self.coordinates.coordinates_values(split=split) return self.coordinates.coordinates_values(split=split)
......
...@@ -28,6 +28,7 @@ class TestSmoothMarginEstimator(unittest.TestCase): ...@@ -28,6 +28,7 @@ class TestSmoothMarginEstimator(unittest.TestCase):
# Plot # Plot
if self.DISPLAY: if self.DISPLAY:
margin_model.margin_function_sample.visualize(show=True) margin_model.margin_function_sample.visualize(show=True)
estimator.margin_function_fitted.visualize(show=True)
self.assertTrue(True) self.assertTrue(True)
......
Supports Markdown
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment