diff --git a/extreme_fit/model/margin_model/polynomial_margin_model/utils.py b/extreme_fit/model/margin_model/polynomial_margin_model/utils.py
index a788af48c84bc346c6da28303c18a42ac94fbd6f..44d7971e175fa8f83fa2ea327c0de1da0c9e929b 100644
--- a/extreme_fit/model/margin_model/polynomial_margin_model/utils.py
+++ b/extreme_fit/model/margin_model/polynomial_margin_model/utils.py
@@ -11,7 +11,6 @@ from extreme_fit.model.margin_model.polynomial_margin_model.spatio_temporal_poly
     NonStationaryLocationSpatioTemporalLinearityModel5, NonStationaryLocationSpatioTemporalLinearityModelAssertError1, \
     NonStationaryLocationSpatioTemporalLinearityModelAssertError2, \
     NonStationaryLocationSpatioTemporalLinearityModelAssertError3, NonStationaryLocationSpatioTemporalLinearityModel6
-
 ALTITUDINAL_MODELS = [
                          StationaryAltitudinal,
                          NonStationaryAltitudinalLocationLinear,
@@ -24,7 +23,8 @@ ALTITUDINAL_MODELS = [
                          NonStationaryAltitudinalLocationQuadraticCrossTermForLocation,
                          NonStationaryAltitudinalLocationLinearScaleLinearCrossTermForLocation,
                         NonStationaryAltitudinalLocationQuadraticScaleLinearCrossTermForLocation,
-                     ][:]
+                     ][:7]
+
 
 
 
diff --git a/extreme_fit/model/utils.py b/extreme_fit/model/utils.py
index 1b6e82d34622c7252dfeb29ece1ffdc386f1da88..3d183036f879c00fb097772c78a66f085e8374be 100644
--- a/extreme_fit/model/utils.py
+++ b/extreme_fit/model/utils.py
@@ -182,7 +182,7 @@ def new_coef_name_to_old_coef_names():
 def get_margin_formula_extremes(fit_marge_form_dict) -> Dict:
     v_to_str = lambda v: ' '.join(v.split()[2:]) if v != 'NULL' else ' 1'
     form_dict = {
-        k: '~' + ' + '.join(
+        k: '~ ' + ' + '.join(
             [v_to_str(fit_marge_form_dict[e]) for e in l if e in fit_marge_form_dict])
         for k, l in new_coef_name_to_old_coef_names().items()
     }
diff --git a/projects/altitude_spatial_model/altitudes_fit/altitudes_studies.py b/projects/altitude_spatial_model/altitudes_fit/altitudes_studies.py
index acd7f9d5c010670fe8a4b752ca330bb9b89fc8ad..a3e071e0968a9bc148ada06fdb9bb45691ae39a2 100644
--- a/projects/altitude_spatial_model/altitudes_fit/altitudes_studies.py
+++ b/projects/altitude_spatial_model/altitudes_fit/altitudes_studies.py
@@ -139,7 +139,7 @@ class AltitudesStudies(object):
             moment += ' Relative'
         if change is True or change is None:
             moment += ' change (between two block of 30 years) for'
-        moment = 'mean' if not std else 'std'
+        moment += 'mean' if not std else 'std'
         plot_name = '{} of annual maxima of {}'.format(moment, SCM_STUDY_CLASS_TO_ABBREVIATION[self.study_class])
         ax.set_ylabel('{} ({})'.format(plot_name, self.study.variable_unit), fontsize=15)
         ax.set_xlabel('altitudes', fontsize=15)
diff --git a/projects/altitude_spatial_model/altitudes_fit/main_altitudes_studies.py b/projects/altitude_spatial_model/altitudes_fit/main_altitudes_studies.py
index 8c800b5dd67949bd31c9b85a8a363eb0ffc6d764..b2f2d47ee234ec20539950ed8b0d6192a81f0807 100644
--- a/projects/altitude_spatial_model/altitudes_fit/main_altitudes_studies.py
+++ b/projects/altitude_spatial_model/altitudes_fit/main_altitudes_studies.py
@@ -34,9 +34,9 @@ def main():
     study_classes = [SafranSnowfall1Day, SafranSnowfall3Days, SafranSnowfall5Days, SafranSnowfall7Days][:2]
     study_classes = [SafranPrecipitation1Day, SafranPrecipitation3Days, SafranPrecipitation5Days,
                      SafranPrecipitation7Days][:]
-    study_classes = [SafranPrecipitation1Day, SafranSnowfall1Day, SafranSnowfall3Days, SafranPrecipitation3Days][:1]
+    study_classes = [SafranSnowfall1Day, SafranPrecipitation1Day, SafranSnowfall3Days, SafranPrecipitation3Days][:1]
     massif_names = None
-    massif_names = ['Aravis']
+    # massif_names = ['Aravis']
     # massif_names = ['Chartreuse', 'Belledonne']
 
     for study_class in study_classes:
diff --git a/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/one_fold_fit.py b/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/one_fold_fit.py
index 1ed8d6795523d4ec762fed46b0f3b142d16c02c7..302e9a0914b2a9fdb9f8b137d1557c85d2f6160e 100644
--- a/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/one_fold_fit.py
+++ b/projects/altitude_spatial_model/altitudes_fit/one_fold_analysis/one_fold_fit.py
@@ -1,4 +1,6 @@
+import numpy.testing as npt
 import numpy as np
+import rpy2
 from cached_property import cached_property
 
 from extreme_fit.estimator.margin_estimator.utils import fitted_linear_margin_estimator_short
@@ -49,10 +51,12 @@ class OneFoldFit(object):
         for estimator in estimators:
             try:
                 aic = estimator.aic()
+                npt.assert_almost_equal(estimator.result_from_model_fit.aic, aic, decimal=5)
                 print(self.massif_name, estimator.margin_model.name_str, aic)
                 estimators_with_finite_aic.append(estimator)
-            except AssertionError:
+            except (AssertionError, rpy2.rinterface.RRuntimeError):
                 print(self.massif_name, estimator.margin_model.name_str, 'infinite aic')
+        print('Summary {}: {}/{} fitted'.format(self.massif_name, len(estimators_with_finite_aic), len(estimators)))
         sorted_estimators_with_finite_aic = sorted([estimator for estimator in estimators_with_finite_aic],
                                                    key=lambda e: e.aic())
         return sorted_estimators_with_finite_aic
diff --git a/test/test_extreme_fit/test_estimator/test_gev_spatio_temporal_extremes_mle.py b/test/test_extreme_fit/test_estimator/test_gev_spatio_temporal_extremes_mle.py
index c68d57f139dac4a7769105f8388408e813b90dec..53d35149a3bb8cbb939e629403040c900c820ffe 100644
--- a/test/test_extreme_fit/test_estimator/test_gev_spatio_temporal_extremes_mle.py
+++ b/test/test_extreme_fit/test_estimator/test_gev_spatio_temporal_extremes_mle.py
@@ -1,7 +1,10 @@
 import unittest
 from random import sample
 
-from extreme_data.meteo_france_data.scm_models_data.safran.safran import SafranSnowfall1Day
+from extreme_data.meteo_france_data.scm_models_data.safran.safran import SafranSnowfall1Day, SafranPrecipitation1Day
+from extreme_fit.model.margin_model.polynomial_margin_model.altitudinal_models import \
+    NonStationaryAltitudinalLocationQuadraticScaleLinearCrossTermForLocation, \
+    NonStationaryAltitudinalLocationQuadraticCrossTermForLocation
 from extreme_fit.model.margin_model.polynomial_margin_model.utils import ALTITUDINAL_MODELS, \
     MODELS_THAT_SHOULD_RAISE_AN_ASSERTION_ERROR, VARIOUS_SPATIO_TEMPORAL_MODELS
 from extreme_fit.model.margin_model.utils import \
@@ -14,16 +17,19 @@ from projects.altitude_spatial_model.altitudes_fit.two_fold_analysis.two_fold_fi
 
 class TestGevTemporalQuadraticExtremesMle(unittest.TestCase):
 
+    def setUp(self) -> None:
+        self.altitudes = [900, 1200]
+        self.massif_name = 'Vercors'
+        self.study_class = SafranSnowfall1Day
+
     def get_estimator_fitted(self, model_class):
-        altitudes = [900, 1200]
-        study_class = SafranSnowfall1Day
-        studies = AltitudesStudies(study_class, altitudes, year_max=2019)
-        two_fold_datasets_generator = TwoFoldDatasetsGenerator(studies, nb_samples=1, massif_names=['Vercors'])
+        studies = AltitudesStudies(self.study_class, self.altitudes, year_max=2019)
+        two_fold_datasets_generator = TwoFoldDatasetsGenerator(studies, nb_samples=1, massif_names=[self.massif_name])
         model_family_name_to_model_class = {'Non stationary': [model_class]}
         two_fold_fit = TwoFoldFit(two_fold_datasets_generator=two_fold_datasets_generator,
                                   model_family_name_to_model_classes=model_family_name_to_model_class,
                                   fit_method=MarginFitMethod.extremes_fevd_mle)
-        massif_fit = two_fold_fit.massif_name_to_massif_fit['Vercors']
+        massif_fit = two_fold_fit.massif_name_to_massif_fit[self.massif_name]
         sample_fit = massif_fit.sample_id_to_sample_fit[0]
         model_fit = sample_fit.model_class_to_model_fit[model_class]  # type: TwoFoldModelFit
         estimator = model_fit.estimator_fold_1
@@ -37,18 +43,55 @@ class TestGevTemporalQuadraticExtremesMle(unittest.TestCase):
         self.assertAlmostEqual(estimator.result_from_model_fit.bic, estimator.bic(split=estimator.train_split))
 
     def test_assert_error(self):
-        for model_class in sample(MODELS_THAT_SHOULD_RAISE_AN_ASSERTION_ERROR, 1):
+        for model_class in MODELS_THAT_SHOULD_RAISE_AN_ASSERTION_ERROR:
             with self.assertRaises(AssertionError):
                 self.common_test(model_class)
 
     def test_location_spatio_temporal_models(self):
-        for model_class in sample(VARIOUS_SPATIO_TEMPORAL_MODELS, 3):
+        for model_class in VARIOUS_SPATIO_TEMPORAL_MODELS:
             self.common_test(model_class)
 
     def test_altitudinal_models(self):
-        for model_class in sample(ALTITUDINAL_MODELS, 3):
+        for model_class in ALTITUDINAL_MODELS:
+            print(model_class)
             self.common_test(model_class)
 
 
+# class MyTest(unittest.TestCase):
+#
+#     def setUp(self) -> None:
+#         self.study_class = SafranPrecipitation1Day
+#         self.altitudes = [900, 1200, 1500, 1800, 2100, 2400, 2700, 3000]
+#         self.massif_name = 'Aravis'
+#
+#     def get_estimator_fitted(self, model_class):
+#         studies = AltitudesStudies(self.study_class, self.altitudes, year_max=2019)
+#         two_fold_datasets_generator = TwoFoldDatasetsGenerator(studies, nb_samples=1, massif_names=[self.massif_name])
+#         model_family_name_to_model_class = {'Non stationary': [model_class]}
+#         two_fold_fit = TwoFoldFit(two_fold_datasets_generator=two_fold_datasets_generator,
+#                                   model_family_name_to_model_classes=model_family_name_to_model_class,
+#                                   fit_method=MarginFitMethod.extremes_fevd_mle)
+#         massif_fit = two_fold_fit.massif_name_to_massif_fit[self.massif_name]
+#         sample_fit = massif_fit.sample_id_to_sample_fit[0]
+#         model_fit = sample_fit.model_class_to_model_fit[model_class]  # type: TwoFoldModelFit
+#         estimator = model_fit.estimator_fold_1
+#         return estimator
+#
+#     def common_test(self, model_class):
+#         estimator = self.get_estimator_fitted(model_class)
+#         # Assert that indicators are correctly computed
+#         self.assertAlmostEqual(estimator.result_from_model_fit.nllh, estimator.nllh(split=estimator.train_split))
+#         self.assertAlmostEqual(estimator.result_from_model_fit.aic, estimator.aic(split=estimator.train_split))
+#         self.assertAlmostEqual(estimator.result_from_model_fit.bic, estimator.bic(split=estimator.train_split))
+#
+#     # def test_altitudinal_models(self):
+#     #     for model_class in ALTITUDINAL_MODELS:
+#     #         self.common_test(model_class)
+#
+#     def test_wrong(self):
+#         self.common_test(NonStationaryAltitudinalLocationQuadraticCrossTermForLocation)
+#         # self.common_test(NonStationaryAltitudinalLocationQuadraticScaleLinearCrossTermForLocation)
+
+
 if __name__ == '__main__':
     unittest.main()