diff --git a/experiment/meteo_france_data/stations_data/comparison_analysis.py b/experiment/meteo_france_data/stations_data/comparison_analysis.py
index ee5160279ff0ab29cadfe7a968ae0ca9bb98dbdc..9b1a9d802742171f38058449f03533ee4a869078 100644
--- a/experiment/meteo_france_data/stations_data/comparison_analysis.py
+++ b/experiment/meteo_france_data/stations_data/comparison_analysis.py
@@ -280,11 +280,11 @@ class ComparisonAnalysis(object):
                                                                        margin_model=margin_model,
                                                                        max_stable_model=max_stable_model)
                 estimator.fit()
-                print(estimator.result_from_model_fit.margin_coef_dict)
+                print(estimator.result_from_model_fit.margin_coef_ordered_dict)
                 print(estimator.result_from_model_fit.other_coef_dict)
                 estimators.append(estimator)
             # Compare the sign of them margin coefficient for the estimators
-            coefs = [{k: v for k, v in e.result_from_model_fit.margin_coef_dict.items() if 'Coeff1' not in k} for e in
+            coefs = [{k: v for k, v in e.result_from_model_fit.margin_coef_ordered_dict.items() if 'Coeff1' not in k} for e in
                      estimators]
             different_sign = [k for k, v in coefs[0].items() if np.sign(coefs[1][k]) != np.sign(v)]
             print('All linear coefficient have the same sign: {}, different_signs for: {}'.format(
diff --git a/extreme_fit/estimator/utils.py b/extreme_fit/estimator/utils.py
index dec6231828cbc6e60f2aaff9ca4fee0db0ae2f4c..315b5e8a1e99f399a266a67e664040f1242d4a6c 100644
--- a/extreme_fit/estimator/utils.py
+++ b/extreme_fit/estimator/utils.py
@@ -5,7 +5,7 @@ from extreme_fit.model.margin_model.margin_function.linear_margin_function impor
 
 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_dict
+        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,
diff --git a/extreme_fit/model/result_from_model_fit/abstract_result_from_model_fit.py b/extreme_fit/model/result_from_model_fit/abstract_result_from_model_fit.py
index 5560e1c862877450498c1f9807fb03cc71465456..3a91b1dfa6497cf9be9d8ee1a476062784a0f04f 100644
--- a/extreme_fit/model/result_from_model_fit/abstract_result_from_model_fit.py
+++ b/extreme_fit/model/result_from_model_fit/abstract_result_from_model_fit.py
@@ -20,9 +20,13 @@ class AbstractResultFromModelFit(object):
         raise NotImplementedError
 
     @property
-    def margin_coef_dict(self):
+    def margin_coef_ordered_dict(self):
         raise NotImplementedError
 
+    @property
+    def margin_coef_ordered_names(self):
+        return list(self.margin_coef_ordered_dict.keys())
+
     @property
     def other_coef_dict(self):
         raise NotImplementedError
@@ -39,6 +43,10 @@ class AbstractResultFromModelFit(object):
     def convergence(self) -> str:
         raise NotImplementedError
 
+    @property
+    def covariance(self):
+        raise NotImplementedError
+
 
 
 
diff --git a/extreme_fit/model/result_from_model_fit/result_from_extremes.py b/extreme_fit/model/result_from_model_fit/result_from_extremes.py
index 9bbb1e22f9aa425bf3c890e8f9df2f95e2fd0977..dce2219d7624d4e4a893581f7d74418fbf9e31b7 100644
--- a/extreme_fit/model/result_from_model_fit/result_from_extremes.py
+++ b/extreme_fit/model/result_from_model_fit/result_from_extremes.py
@@ -5,7 +5,7 @@ from rpy2 import robjects
 from extreme_fit.distribution.gev.gev_params import GevParams
 from extreme_fit.model.result_from_model_fit.abstract_result_from_model_fit import \
     AbstractResultFromModelFit
-from extreme_fit.model.result_from_model_fit.utils import get_margin_coef_dict
+from extreme_fit.model.result_from_model_fit.utils import get_margin_coef_ordered_dict
 from extreme_fit.model.utils import r
 
 
@@ -40,10 +40,10 @@ class ResultFromExtremes(AbstractResultFromModelFit):
     def get_coef_dict_from_posterior_sample(self, s: pd.Series):
         assert len(s) >= 3
         values = {i: v for i, v in enumerate(s)}
-        return get_margin_coef_dict(self.gev_param_name_to_dim, values)
+        return get_margin_coef_ordered_dict(self.gev_param_name_to_dim, values)
 
     @property
-    def margin_coef_dict(self):
+    def margin_coef_ordered_dict(self):
         """ It is the coef for the margin function corresponding to the mean posterior parameters """
         mean_posterior_parameters = self.df_posterior_samples.iloc[:, :-2].mean(axis=0)
         return self.get_coef_dict_from_posterior_sample(mean_posterior_parameters)
diff --git a/extreme_fit/model/result_from_model_fit/result_from_ismev.py b/extreme_fit/model/result_from_model_fit/result_from_ismev.py
index d83dacecbc9b1610fd35e04092feecc3ceade29e..2332a6608d920736b671aa8bf1aa9f73645403b4 100644
--- a/extreme_fit/model/result_from_model_fit/result_from_ismev.py
+++ b/extreme_fit/model/result_from_model_fit/result_from_ismev.py
@@ -1,13 +1,11 @@
-from rpy2 import robjects
-
-from extreme_fit.model.result_from_model_fit.abstract_result_from_model_fit import \
-    AbstractResultFromModelFit
-from extreme_fit.model.result_from_model_fit.utils import convertFloatVector_to_float, get_margin_coef_dict
+import numpy as np
+import pandas as pd
 from rpy2 import robjects
 
 from extreme_fit.model.result_from_model_fit.abstract_result_from_model_fit import \
     AbstractResultFromModelFit
 from extreme_fit.model.result_from_model_fit.utils import convertFloatVector_to_float
+from extreme_fit.model.result_from_model_fit.utils import get_margin_coef_ordered_dict
 
 
 class ResultFromIsmev(AbstractResultFromModelFit):
@@ -17,12 +15,12 @@ class ResultFromIsmev(AbstractResultFromModelFit):
         self.gev_param_name_to_dim = gev_param_name_to_dim
 
     @property
-    def margin_coef_dict(self):
-        return get_margin_coef_dict(self.gev_param_name_to_dim, self.name_to_value['mle'])
+    def margin_coef_ordered_dict(self):
+        return get_margin_coef_ordered_dict(self.gev_param_name_to_dim, self.name_to_value['mle'])
 
     @property
     def all_parameters(self):
-        return self.margin_coef_dict
+        return self.margin_coef_ordered_dict
 
     @property
     def nllh(self):
@@ -35,3 +33,9 @@ class ResultFromIsmev(AbstractResultFromModelFit):
     @property
     def convergence(self) -> str:
         return convertFloatVector_to_float(self.name_to_value['conv']) == 0
+
+    @property
+    def covariance(self):
+        return pd.DataFrame(np.array(self.name_to_value['cov']),
+                            index=self.margin_coef_ordered_names,
+                            columns=self.margin_coef_ordered_names)
diff --git a/extreme_fit/model/result_from_model_fit/result_from_spatial_extreme.py b/extreme_fit/model/result_from_model_fit/result_from_spatial_extreme.py
index 8e46fde4f763df41beba4a53e84b684784e0ecca..ed7156a455e580af9ac035885b15746610a54383 100644
--- a/extreme_fit/model/result_from_model_fit/result_from_spatial_extreme.py
+++ b/extreme_fit/model/result_from_model_fit/result_from_spatial_extreme.py
@@ -33,7 +33,7 @@ class ResultFromSpatialExtreme(AbstractResultFromModelFit):
         return {key: fitted_values.rx2(key)[0] for key in fitted_values.names}
 
     @property
-    def margin_coef_dict(self):
+    def margin_coef_ordered_dict(self):
         return {k: v for k, v in self.all_parameters.items() if LinearCoef.COEFF_STR in k}
 
     @property
diff --git a/extreme_fit/model/result_from_model_fit/utils.py b/extreme_fit/model/result_from_model_fit/utils.py
index 6d559b4f79abc69925565f3f5cb54a2a02dd4941..7dabcc92a3b20fa409488dbefa93d871563d168e 100644
--- a/extreme_fit/model/result_from_model_fit/utils.py
+++ b/extreme_fit/model/result_from_model_fit/utils.py
@@ -1,3 +1,5 @@
+from collections import OrderedDict
+
 import numpy as np
 
 from extreme_fit.distribution.gev.gev_params import GevParams
@@ -9,10 +11,10 @@ def convertFloatVector_to_float(f):
     return np.array(f)[0]
 
 
-def get_margin_coef_dict(gev_param_name_to_dim, mle_values):
+def get_margin_coef_ordered_dict(gev_param_name_to_dim, mle_values):
     assert gev_param_name_to_dim is not None
     # Build the Coeff dict from gev_param_name_to_dim
-    coef_dict = {}
+    coef_dict = OrderedDict()
     i = 0
     for gev_param_name in GevParams.PARAM_NAMES:
         # Add intercept
diff --git a/test/test_extreme_fit/test_distribution/test_gev/test_gev_temporal.py b/test/test_extreme_fit/test_distribution/test_gev/test_gev_temporal.py
index 8e96aa051ad031d70a32c4b44b74d8884294079d..116c041e60357d80307e66787879c6e06e228080 100644
--- a/test/test_extreme_fit/test_distribution/test_gev/test_gev_temporal.py
+++ b/test/test_extreme_fit/test_distribution/test_gev/test_gev_temporal.py
@@ -22,7 +22,7 @@ class TestGevTemporal(unittest.TestCase):
         set_seed_r()
         r("""
         N <- 50
-        loc = 0; scale = 1; shape <- 1
+        loc = 0; scale = 2; shape <- 1
         x_gev <- rgev(N, loc = loc, scale = scale, shape = shape)
         start_loc = 0; start_scale = 1; start_shape = 1
         """)
@@ -39,7 +39,7 @@ class TestGevTemporal(unittest.TestCase):
         margin_model = StationaryTemporalModel(self.coordinates)
         estimator = LinearMarginEstimator(self.dataset, margin_model)
         estimator.fit()
-        ref = {'loc': 0.0219, 'scale': 1.0347, 'shape': 0.8295}
+        ref = {'loc': 0.04309190816463247, 'scale': 2.0688696961628437, 'shape': 0.8291528207825063}
         for year in range(1, 3):
             mle_params_estimated = estimator.margin_function_from_fit.get_gev_params(np.array([year])).to_dict()
             for key in ref.keys():