diff --git a/extreme_fit/estimator/margin_estimator/abstract_margin_estimator.py b/extreme_fit/estimator/margin_estimator/abstract_margin_estimator.py
index cc8ae95541b8a74df94b6598a08fffc5845c6a95..d4b756e1b84bc54388ad810688c76e1b50c49a0e 100644
--- a/extreme_fit/estimator/margin_estimator/abstract_margin_estimator.py
+++ b/extreme_fit/estimator/margin_estimator/abstract_margin_estimator.py
@@ -3,11 +3,12 @@ from abc import ABC
from cached_property import cached_property
from extreme_fit.estimator.abstract_estimator import AbstractEstimator
-from extreme_fit.estimator.utils import load_margin_function
+from extreme_fit.estimator.utils import load_margin_function, compute_nllh
from extreme_fit.model.margin_model.linear_margin_model.linear_margin_model import LinearMarginModel
from extreme_fit.model.margin_model.margin_function.linear_margin_function import LinearMarginFunction
from extreme_fit.model.result_from_model_fit.abstract_result_from_model_fit import AbstractResultFromModelFit
from spatio_temporal_dataset.dataset.abstract_dataset import AbstractDataset
+from spatio_temporal_dataset.slicer.split import Split
class AbstractMarginEstimator(AbstractEstimator, ABC):
@@ -26,13 +27,24 @@ class LinearMarginEstimator(AbstractMarginEstimator):
self.margin_model = margin_model
def _fit(self) -> AbstractResultFromModelFit:
- 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_for_fit(split=self.train_split,
- starting_point=self.margin_model.starting_point)
- return self.margin_model.fitmargin_from_maxima_gev(data=maxima_gev_specialized,
+ return self.margin_model.fitmargin_from_maxima_gev(data=self.maxima_gev_train,
df_coordinates_spat=df_coordinates_spat,
- df_coordinates_temp=df_coordinates_temp)
+ df_coordinates_temp=self.coordinate_temp)
+
+ @property
+ def coordinate_temp(self):
+ return self.dataset.coordinates.df_temporal_coordinates_for_fit(split=self.train_split,
+ starting_point=self.margin_model.starting_point)
+
+ @property
+ def maxima_gev_train(self):
+ return self.dataset.maxima_gev_for_spatial_extremes_package(self.train_split)
+
+ @property
+ def nllh(self, split=Split.all):
+ assert split is Split.all
+ return compute_nllh(self, self.maxima_gev_train, self.coordinate_temp, self.margin_model)
@cached_property
def margin_function_from_fit(self) -> LinearMarginFunction:
diff --git a/extreme_fit/estimator/utils.py b/extreme_fit/estimator/utils.py
index 315b5e8a1e99f399a266a67e664040f1242d4a6c..5ffa53e8d80353e222b8b435a31ac763cd82de54 100644
--- a/extreme_fit/estimator/utils.py
+++ b/extreme_fit/estimator/utils.py
@@ -1,12 +1,26 @@
+import numpy as np
+
from extreme_fit.estimator.abstract_estimator import AbstractEstimator
from extreme_fit.model.margin_model.linear_margin_model.linear_margin_model import LinearMarginModel
from extreme_fit.model.margin_model.margin_function.linear_margin_function import LinearMarginFunction
-def load_margin_function(estimator: AbstractEstimator, margin_model: LinearMarginModel, margin_function_class=LinearMarginFunction, coef_dict=None):
+def load_margin_function(estimator: AbstractEstimator, margin_model: LinearMarginModel,
+ margin_function_class=LinearMarginFunction, coef_dict=None):
if coef_dict is None:
coef_dict = estimator.result_from_model_fit.margin_coef_ordered_dict
return margin_function_class.from_coef_dict(coordinates=estimator.dataset.coordinates,
gev_param_name_to_dims=margin_model.margin_function_start_fit.gev_param_name_to_dims,
coef_dict=coef_dict,
starting_point=margin_model.starting_point)
+
+
+def compute_nllh(estimator: AbstractEstimator, maxima, coordinate_temp, margin_model: LinearMarginModel,
+ margin_function_class=LinearMarginFunction, coef_dict=None):
+ margin_function = load_margin_function(estimator, margin_model, margin_function_class, coef_dict)
+ nllh = 0
+ for maximum, year in zip(maxima[0], coordinate_temp.values):
+ gev_params = margin_function.get_gev_params(year, is_transformed=False)
+ p = gev_params.density(maximum)
+ nllh -= np.log(p)
+ return nllh
diff --git a/test/test_extreme_fit/test_distribution/test_gev/test_gev_temporal_extremes_mle.py b/test/test_extreme_fit/test_distribution/test_gev/test_gev_temporal_extremes_mle.py
index f1dc68ff26a85c86f0dc6fa947c7a565ed521a09..94f989744077be01fb8c40f801e12f513b1a5c57 100644
--- a/test/test_extreme_fit/test_distribution/test_gev/test_gev_temporal_extremes_mle.py
+++ b/test/test_extreme_fit/test_distribution/test_gev/test_gev_temporal_extremes_mle.py
@@ -46,6 +46,7 @@ class TestGevTemporalExtremesMle(unittest.TestCase):
mle_params_estimated = estimator.margin_function_from_fit.get_gev_params(np.array([year])).to_dict()
for key in ref.keys():
self.assertAlmostEqual(ref[key], mle_params_estimated[key], places=3)
+ self.assertAlmostEqual(estimator.result_from_model_fit.nllh, estimator.nllh)
def test_gev_temporal_margin_fit_non_stationary_location(self):
# Create estimator
@@ -56,6 +57,7 @@ class TestGevTemporalExtremesMle(unittest.TestCase):
mle_params_estimated_year1 = estimator.margin_function_from_fit.get_gev_params(np.array([1])).to_dict()
mle_params_estimated_year3 = estimator.margin_function_from_fit.get_gev_params(np.array([3])).to_dict()
self.assertNotEqual(mle_params_estimated_year1, mle_params_estimated_year3)
+ self.assertAlmostEqual(estimator.result_from_model_fit.nllh, estimator.nllh)
def test_gev_temporal_margin_fit_non_stationary_location_and_scale(self):
# Create estimator
@@ -67,6 +69,7 @@ class TestGevTemporalExtremesMle(unittest.TestCase):
mle_params_estimated_year1 = estimator.margin_function_from_fit.get_gev_params(np.array([1])).to_dict()
mle_params_estimated_year3 = estimator.margin_function_from_fit.get_gev_params(np.array([3])).to_dict()
self.assertNotEqual(mle_params_estimated_year1, mle_params_estimated_year3)
+ self.assertAlmostEqual(estimator.result_from_model_fit.nllh, estimator.nllh)
if __name__ == '__main__':