From d8718ec5c543e075dd7a4e2bc4c6f3a7c5abbacd Mon Sep 17 00:00:00 2001
From: Le Roux Erwan <erwan.le-roux@irstea.fr>
Date: Thu, 9 May 2019 15:34:32 +0200
Subject: [PATCH] [SCM][NON STATIONARITY] improve some display before first
experiments.
---
.../meteo_france_SCM_study/abstract_study.py | 11 +++----
.../main_study_visualizer.py | 32 ++++++++++++++-----
.../non_stationary_trends.py | 27 ++++++++++------
.../study_visualization/study_visualizer.py | 2 +-
.../visualization/utils.py | 17 +++++++++-
5 files changed, 63 insertions(+), 26 deletions(-)
diff --git a/experiment/meteo_france_SCM_study/abstract_study.py b/experiment/meteo_france_SCM_study/abstract_study.py
index 37c067a3..f29b6a5a 100644
--- a/experiment/meteo_france_SCM_study/abstract_study.py
+++ b/experiment/meteo_france_SCM_study/abstract_study.py
@@ -31,7 +31,7 @@ with redirect_stdout(f):
class AbstractStudy(object):
"""
- Les fichiers netcdf sont autodocumentés (on peut les comprendre avec ncdump -h notamment).
+ Les fichiers netcdf de SAFRAN et CROCUS sont autodocumentés (on peut les comprendre avec ncdump -h notamment).
"""
REANALYSIS_FOLDER = 'alp_flat/reanalysis'
@@ -164,6 +164,8 @@ class AbstractStudy(object):
df_centroid = self.load_df_centroid()
for coord_column in [AbstractCoordinates.COORDINATE_X, AbstractCoordinates.COORDINATE_Y]:
df_centroid.loc[:, coord_column] = df_centroid[coord_column].str.replace(',', '.').astype(float)
+ # Filter, keep massifs present at the altitude of interest
+ df_centroid = df_centroid.loc[self.study_massif_names]
# Build coordinate object from df_centroid
return AbstractSpatialCoordinates.from_df(df_centroid)
@@ -311,7 +313,8 @@ class AbstractStudy(object):
@property
def title(self):
- return "{} at altitude {}m".format(self.variable_name, self.altitude)
+ return "{} at altitude {}m ({} mountain chains)".format(self.variable_name, self.altitude,
+ len(self.study_massif_names))
@property
def variable_name(self):
@@ -344,7 +347,3 @@ class AbstractStudy(object):
assert self.model_name in ['Safran', 'Crocus']
study_folder = 'meteo' if self.model_name is 'Safran' else 'pro'
return op.join(self.full_path, self.REANALYSIS_FOLDER, study_folder)
-
-
-
-
diff --git a/experiment/meteo_france_SCM_study/visualization/study_visualization/main_study_visualizer.py b/experiment/meteo_france_SCM_study/visualization/study_visualization/main_study_visualizer.py
index fbd7d185..358ff9ef 100644
--- a/experiment/meteo_france_SCM_study/visualization/study_visualization/main_study_visualizer.py
+++ b/experiment/meteo_france_SCM_study/visualization/study_visualization/main_study_visualizer.py
@@ -12,6 +12,14 @@ SCM_EXTENDED_STUDIES = [ExtendedSafranSnowfall, ExtendedCrocusSwe, ExtendedCrocu
SCM_STUDY_TO_EXTENDED_STUDY = OrderedDict(zip(SCM_STUDIES, SCM_EXTENDED_STUDIES))
+def study_iterator_global(study_classes, only_first_one=False, both_altitude=False, verbose=True, altitudes=None):
+ for study_class in study_classes:
+ for study in study_iterator(study_class, only_first_one, both_altitude, verbose, altitudes):
+ yield study
+ if only_first_one:
+ break
+
+
def study_iterator(study_class, only_first_one=False, both_altitude=False, verbose=True, altitudes=None):
all_studies = []
is_safran_study = study_class in [SafranSnowfall, ExtendedSafranSnowfall]
@@ -20,11 +28,16 @@ def study_iterator(study_class, only_first_one=False, both_altitude=False, verbo
print('Loading studies....')
for nb_day in nb_days:
altis = [1800] if altitudes is None else altitudes
+
for alti in altis:
+
if verbose:
print('alti: {}, nb_day: {}'.format(alti, nb_day))
study = study_class(altitude=alti, nb_consecutive_days=nb_day) if is_safran_study \
else study_class(altitude=alti)
+ massifs = study.altitude_to_massif_names[alti]
+ if verbose:
+ print('{} massifs: {}'.format(len(massifs), massifs))
yield study
if only_first_one and not both_altitude:
break
@@ -75,7 +88,8 @@ def normal_visualization(temporal_non_stationarity=False):
# for study_class in SCM_STUDIES[:1]:
for study_class in [CrocusDepth, SafranSnowfall, SafranRainfall, SafranTemperature][1:2]:
for study in study_iterator(study_class, only_first_one=only_first_one):
- study_visualizer = StudyVisualizer(study, save_to_file=save_to_file, temporal_non_stationarity=temporal_non_stationarity)
+ 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)
@@ -92,17 +106,19 @@ def complete_analysis(only_first_one=False):
study_visualizer.visualize_all_experimental_law()
print('Study normal')
for study in study_iterator(study_class, only_first_one=only_first_one):
- study_visualizer = StudyVisualizer(study, save_to_file=True)
- study_visualizer.visualize_linear_margin_fit()
+ study_visualizer = StudyVisualizer(study, save_to_file=True)
+ study_visualizer.visualize_linear_margin_fit()
def trend_analysis():
- save_to_file = True
+ save_to_file = False
only_first_one = True
- for study_class in [CrocusSwe, CrocusDepth, SafranSnowfall, SafranRainfall, SafranTemperature][:3]:
- for study in study_iterator(study_class, only_first_one=only_first_one, altitudes=[1800, 2100, 2400, 2700]):
- study_visualizer = StudyVisualizer(study, save_to_file=save_to_file)
- study_visualizer.visualize_temporal_trend_relevance(complete_analysis=not only_first_one)
+ altitudes = [300, 1200, 2100, 3000][3:]
+ study_classes = [CrocusSwe, CrocusDepth, SafranSnowfall, SafranRainfall, SafranTemperature]
+ for study in study_iterator_global(study_classes, only_first_one=only_first_one, altitudes=altitudes):
+ study_visualizer = StudyVisualizer(study, save_to_file=save_to_file)
+ study_visualizer.visualize_temporal_trend_relevance(complete_analysis=False)
+
if __name__ == '__main__':
# annual_mean_vizu_compare_durand_study(safran=True, take_mean_value=True, altitude=2100)
diff --git a/experiment/meteo_france_SCM_study/visualization/study_visualization/non_stationary_trends.py b/experiment/meteo_france_SCM_study/visualization/study_visualization/non_stationary_trends.py
index 6d30e662..4300cf2f 100644
--- a/experiment/meteo_france_SCM_study/visualization/study_visualization/non_stationary_trends.py
+++ b/experiment/meteo_france_SCM_study/visualization/study_visualization/non_stationary_trends.py
@@ -1,5 +1,6 @@
from typing import Union
+from experiment.meteo_france_SCM_study.visualization.utils import align_yaxis_on_zero
from extreme_estimator.estimator.abstract_estimator import AbstractEstimator
from scipy.stats import chi2
from extreme_estimator.estimator.full_estimator.abstract_full_estimator import \
@@ -68,6 +69,14 @@ class AbstractNonStationaryTrendTest(object):
ax.plot(years, difference, color_difference + 'o-')
ax.set_ylabel(self.RESULT_ATTRIBUTE_METRIC + ' difference', color=color_difference)
+ # Plot zero line
+ years_line = [years[0] -10, years[-1] + 10]
+ ax.plot(years_line, [0 for _ in years_line], 'k-', label='zero line')
+ # Plot significative line corresponding to 0.05 relevance
+ alpha = 0.05
+ significative_deviance = chi2.ppf(q=1 - alpha, df=1)
+ ax.plot(years_line, [significative_deviance for _ in years_line], 'g-', label='significative line')
+
# Plot the mu1 parameter
mu1_trends = [self.get_mu1(starting_point=year) for year in years]
ax2 = ax.twinx()
@@ -75,16 +84,14 @@ class AbstractNonStationaryTrendTest(object):
ax.plot(years, mu1_trends, color_mu1 + 'o-')
ax2.set_ylabel('mu1 parameter', color=color_mu1)
- # Plot zero line
- ax.plot(years, [0 for _ in years], 'k-', label='zero line')
- # Plot significative line corresponding to 0.05 relevance
- alpha = 0.05
- significative_deviance = chi2.ppf(q=1 - alpha, df=1)
- ax.plot(years, [significative_deviance for _ in years], 'g-', label='significative line')
- # Add some informations about the graph
-
- ax.set_xlabel('year')
- ax.set_title(self.display_name)
+ ax.set_xlabel('starting year for the linear trend of mu1')
+ align_yaxis_on_zero(ax, ax2)
+ title = self.display_name
+ first = mu1_trends[0]
+ last = mu1_trends[-1]
+ title += '\n mu1: first {}\nlast {}'.format(first, last)
+ title += '; equals? {}'.format(first == last)
+ ax.set_title(title)
ax.legend()
@property
diff --git a/experiment/meteo_france_SCM_study/visualization/study_visualization/study_visualizer.py b/experiment/meteo_france_SCM_study/visualization/study_visualization/study_visualizer.py
index ea9d5b5c..799ad66a 100644
--- a/experiment/meteo_france_SCM_study/visualization/study_visualization/study_visualizer.py
+++ b/experiment/meteo_france_SCM_study/visualization/study_visualization/study_visualizer.py
@@ -141,7 +141,7 @@ class StudyVisualizer(object):
trend_tests = [ConditionalIndedendenceLocationTrendTest(self.dataset)]
max_stable_models = load_test_max_stable_models(default_covariance_function=self.default_covariance_function)
- for max_stable_model in max_stable_models[:1]:
+ for max_stable_model in [max_stable_models[1], max_stable_models[-2]]:
trend_tests.append(MaxStableLocationTrendTest(self.dataset, max_stable_model))
nb_trend_tests = len(trend_tests)
diff --git a/experiment/meteo_france_SCM_study/visualization/utils.py b/experiment/meteo_france_SCM_study/visualization/utils.py
index 19101351..48ac8b30 100644
--- a/experiment/meteo_france_SCM_study/visualization/utils.py
+++ b/experiment/meteo_france_SCM_study/visualization/utils.py
@@ -35,11 +35,26 @@ def get_km_formatter():
return tkr.FuncFormatter(numfmt) # create your custom formatter function
-def create_adjusted_axes(nb_rows, nb_columns, figsize=(16,10), subplot_space=0.5):
+def create_adjusted_axes(nb_rows, nb_columns, figsize=(16, 10), subplot_space=0.5):
fig, axes = plt.subplots(nb_rows, nb_columns, figsize=figsize)
fig.subplots_adjust(hspace=subplot_space, wspace=subplot_space)
return axes
+def align_yaxis(ax1, v1, ax2, v2):
+ # From https://stackoverflow.com/questions/10481990/matplotlib-axis-with-two-scales-shared-origin
+ """adjust ax2 ylimit so that v2 in ax2 is aligned to v1 in ax1"""
+ _, y1 = ax1.transData.transform((0, v1))
+ _, y2 = ax2.transData.transform((0, v2))
+ inv = ax2.transData.inverted()
+ _, dy = inv.transform((0, 0)) - inv.transform((0, y1 - y2))
+ miny, maxy = ax2.get_ylim()
+ ax2.set_ylim(miny + dy, maxy + dy)
+
+
+def align_yaxis_on_zero(ax1, ax2):
+ align_yaxis(ax1, 0, ax2, 0)
+
+
if __name__ == '__main__':
example_plot_df()
--
GitLab