Aller au contenu

Risques climatiques (kadi.weather.risk)

Le module RiskIndicators analyse les données météo pour produire des alertes opérationnelles : probabilité de pluie imminente et indices de sécheresse multi-méthodes.

Probabilité de pluie

La probabilité de pluie combine deux sources complémentaires :

  1. Prévisions Open-Meteo : proportion de jours avec pluie dans la fenêtre future (information directe sur les prochains jours).
  2. Chaînes de Markov : probabilités de transition calculées sur l'historique (si il a plu hier, quelle est la probabilité d'une pluie aujourd'hui ?).

Combinaison pondérée : 

P_combinée = 0.7 × P_prevision + 0.3 × P_markov

Ce ratio donne plus de poids aux prévisions récentes tout en intégrant la tendance historique.

Indice de sécheresse (SPI)

Le Standardized Precipitation Index (SPI) mesure l'écart des précipitations d'une période par rapport à la normale historique, exprimé en nombre d'écarts-types.

Formule (approximation Z-Score) : 

SPI = (P_période - P_moyenne) / σ_précipitations

Un SPI < -1.0 indique une période sèche, un SPI < -1.5 une sécheresse marquée.

Échelle de sévérité

SPI Sévérité Implication agronomique
> -0.5 no_drought Conditions normales
-0.5 à -1.0 mild Surveillance conseillée
-1.0 à -1.5 moderate Stress hydrique modéré, irrigation possible
< -1.5 severe Risque de pertes de récolte, alerte

Méthodes d'analyse de sécheresse

Méthode Description Usage recommandé
'spi' Z-Score sur les précipitations Suivi courant, simple et rapide
'markov' Probabilité de persistance d'un état sec Prévision à court terme (J+7)
'hurst' Exposant de Hurst (mémoire longue) Tendances multi-annuelles
'combined' Combinaison pondérée des 3 méthodes Vue d'ensemble complète

Utilisation

Probabilité de pluie

from kadi.weather import WeatherSession

session = WeatherSession(latitude=9.3333, longitude=2.6333, name="Parakou")

# Probabilité sur les 3 prochains jours (seuil : 1 mm minimum)
prob = session.rain_probability(days_ahead=3, min_rainfall_mm=1.0)

print(f"Probabilité demain : {prob['tomorrow'] * 100:.0f}%")
print(f"Message            : {prob['message']}")
print(f"Recommandation     : {prob['recommendation']}")

Retour de rain_probability :

Clé Type Description
tomorrow float Probabilité J+1 (0.0 à 1.0)
message str Résumé lisible de la situation
recommendation str Recommandation opérationnelle

Exemples de recommandations : - "Évitez les traitements phytosanitaires — pluie probable demain." - "Bon moment pour les semis — 3 jours secs prévus." - "Risque de lessivage des intrants — attendez 48h."

Indice de sécheresse

# SPI sur 3 mois glissants (suivi courant)
spi3 = session.drought_index(method="spi", window_months=3)
print(f"SPI 3 mois : {spi3['spi_3month']:.2f}")
print(f"Sévérité   : {spi3['drought_severity']}")

# Analyse combinée pour une vue complète
combined = session.drought_index(method="combined", window_months=6)
print(f"Indice combiné : {combined.get('combined_score', 'N/A'):.2f}")
print(f"Sévérité       : {combined['drought_severity']}")

Retour de drought_index :

Clé Type Description
spi_3month float Valeur SPI calculée
drought_severity str Niveau de sécheresse
window_months int Fenêtre temporelle utilisée
method str Méthode appliquée

Intégration avec kadi.market (Phase 4)

Le module de risque est utilisé directement par kadi.market.logistics pour ajuster les coûts de transport selon les conditions climatiques.

from kadi.weather import WeatherSession
from kadi.market import Market

ws = WeatherSession(latitude=9.30, longitude=2.08, name="Parakou")
marche = Market(lat=9.30, lon=2.08, location="Parakou", weather_session=ws)

# Le risque climatique entre dans le calcul logistique automatiquement
cout = marche.logistics.calculate_transfer_cost("Parakou", "Cotonou", crop="tomato")
print(f"Probabilité de pluie utilisée : {cout['prob_pluie'] * 100:.0f}%")
print(f"Gamma route effectif          : {cout['gamma_effectif']:.3f}")

# Vue synthétique du risque climatique
risque = marche.assess_climate_risk(days_ahead=7)
print(risque["recommendation"])

RiskIndicators

Évalue les risques climatiques (sécheresse, probabilité de pluie).

Source code in kadi/weather/risk.py 
 16
 17
 18
 19
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
class RiskIndicators:
    """
    Évalue les risques climatiques (sécheresse, probabilité de pluie).
    """

    def __init__(self, location: Location, rainfall_historical: pd.Series, forecast_data: pd.DataFrame):
        """
        Initialise les indicateurs de risque.

        :param location: Instance de Location.
        :param rainfall_historical: Série historique de pluie.
        :param forecast_data: DataFrame de prévisions météorologiques.
        """
        self.location = location
        self.rainfall_historical = rainfall_historical
        self.forecast_data = forecast_data

    def drought_index(self, method: str = 'spi', window_months: int = 3) -> dict:
        """
        Calcule l'indice de sécheresse avec la méthode spécifiée.

        :param method: Méthode ('spi', 'markov', 'hurst', 'combined').
        :param window_months: Fenêtre temporelle pour le calcul du SPI.
        :return: Dictionnaire avec les résultats de sécheresse.
        """
        results = {}

        if method in ['spi', 'combined']:
            spi_val = self.spi(window_months)
            results[f'spi_{window_months}month'] = spi_val
            results['drought_severity'] = self._get_severity_level(spi_val)

        if method in ['markov', 'combined']:
            markov_res = self.markov_transition()
            results['markov_p_dry'] = markov_res.get('p_dry_dry', 0.0)

        if method in ['hurst', 'combined']:
            # L'exposant de Hurst nécessite une série assez longue (au moins qq années)
            hurst = self.hurst_exponent()
            results['hurst_exponent'] = round(hurst, 2)

        if method not in ['spi', 'markov', 'hurst', 'combined']:
            raise ValidationError(f"Méthode {method} non supportée pour l'indice de sécheresse.")

        return results

    def spi(self, window_months: int) -> float:
        """
        Calcule le Standardized Precipitation Index (SPI) pour la période récente.
        Implémentation simplifiée pour le MVP.

        :param window_months: Fenêtre d'accumulation en mois.
        :return: Valeur du SPI.
        """
        if self.rainfall_historical.empty:
            raise InsufficientData("Aucune donnée historique pour le calcul du SPI.")

        # Cumul de précipitation sur la fenêtre temporelle
        days = window_months * 30
        rolling_sum = self.rainfall_historical.rolling(window=days, min_periods=days//2).sum()
        rolling_sum = rolling_sum.dropna()

        if len(rolling_sum) < 30: # Pas assez de données pour faire une loi Gamma
            raise InsufficientData("Pas assez de jours de données pour ajuster le modèle SPI (minimum 30 requis).")

        # Ajustement d'une distribution Gamma simplifié
        # On évite les zéros pour la distribution Gamma
        valid_data = rolling_sum[rolling_sum > 0]
        if len(valid_data) == 0:
            return -3.0 # Extrême sécheresse

        # Dans un environnement de production, on utiliserait scipy.stats.gamma.fit
        # et une correction de probabilité mixte pour p(x=0). 
        # Ici on fait une approximation simple par score Z (loi normale) pour le mock
        current_val = rolling_sum.iloc[-1]
        mean_val = rolling_sum.mean()
        std_val = rolling_sum.std()

        if std_val == 0:
            return 0.0

        spi_val = (current_val - mean_val) / std_val
        return round(float(spi_val), 2)

    def markov_transition(self, threshold_mm: float = 1.0) -> dict:
        """
        Calcule les probabilités de transition de Markov (jour sec -> jour sec).

        :param threshold_mm: Seuil pour considérer un jour comme humide.
        :return: Dictionnaire avec les probabilités p00, p01, p10, p11.
        """
        if self.rainfall_historical.empty:
            raise InsufficientData("Aucune donnée historique pour le calcul des probabilités de transition de Markov.")

        # États : 0 = sec, 1 = humide
        states = (self.rainfall_historical >= threshold_mm).astype(int)

        # Transitions
        transitions = pd.DataFrame({'current': states.iloc[:-1].values, 'next': states.iloc[1:].values})

        counts = transitions.groupby(['current', 'next']).size().unstack(fill_value=0)

        # On s'assure d'avoir la matrice 2x2
        for i in [0, 1]:
            if i not in counts.index:
                counts.loc[i] = [0, 0]
            for j in [0, 1]:
                if j not in counts.columns:
                    counts[j] = 0

        # Probabilités
        p0 = counts.loc[0].sum()
        p1 = counts.loc[1].sum()

        p00 = counts.loc[0, 0] / p0 if p0 > 0 else 0
        p01 = counts.loc[0, 1] / p0 if p0 > 0 else 0
        p10 = counts.loc[1, 0] / p1 if p1 > 0 else 0
        p11 = counts.loc[1, 1] / p1 if p1 > 0 else 0

        return {
            'p_dry_dry': round(p00, 2),
            'p_dry_wet': round(p01, 2),
            'p_wet_dry': round(p10, 2),
            'p_wet_wet': round(p11, 2)
        }

    def hurst_exponent(self, window: int = 1095) -> float:
        """
        Calcule l'exposant de Hurst par la méthode de gamme rééchelonnée (R/S) multi-échelle.

        L'algorithme segmente la série en sous-fenêtres de tailles croissantes,
        calcule le rapport R/S moyen pour chaque taille, puis estime H par régression
        log-log. Un exposant H > 0.5 indique une persistance climatique (mémoire longue).

        :param window: Taille maximale de la fenêtre d'analyse (jours).
        :return: Exposant de Hurst H (compris entre 0.01 et 0.99).
        """
        if len(self.rainfall_historical) < 100:
            raise InsufficientData("Pas assez de données pour l'exposant de Hurst (minimum 100 jours requis).")

        data = self.rainfall_historical.values
        n_total = len(data)

        # Taille de la plus petite fenêtre (doit être assez grande pour un R/S stable)
        min_w = 10
        # Taille de la plus grande fenêtre (plafonnée à la moitié de la série)
        max_w = min(window, n_total // 2)

        window_sizes = []
        rs_means = []

        # Progression géométrique des tailles de fenêtre (facteur 1.5)
        w = min_w
        while w <= max_w:
            # Calcul du R/S moyen sur toutes les sous-séquences non chevauchantes
            num_segments = n_total // w
            rs_values_w = []

            for k in range(num_segments):
                segment = data[k * w: (k + 1) * w]
                seg_mean = np.mean(segment)
                centered = segment - seg_mean
                cum_dev = np.cumsum(centered)

                # Gamme (R) et écart-type (S) du segment
                r = np.max(cum_dev) - np.min(cum_dev)
                s = np.std(segment)

                if s > 0:
                    rs_values_w.append(r / s)

            if rs_values_w:
                window_sizes.append(w)
                rs_means.append(np.mean(rs_values_w))

            w = max(w + 1, int(w * 1.5))

        # Minimum de 3 points pour une régression fiable
        if len(window_sizes) < 3:
            return 0.5

        # Régression linéaire dans l'espace log-log : log(R/S) = H * log(N) + c
        log_n = np.log(np.array(window_sizes, dtype=float))
        log_rs = np.log(np.array(rs_means, dtype=float))
        coeffs = np.polyfit(log_n, log_rs, 1)

        # H est la pente de la droite de régression
        h = float(coeffs[0])
        return float(np.clip(h, 0.01, 0.99))

    def rain_probability(self, days_ahead: int = 1, min_rainfall_mm: float = 1.0) -> dict:
        """
        Prévoit la probabilité de pluie pour les prochains jours.

        Combine deux sources d'information pour plus de robustesse :
        1. Les prévisions API (Open-Meteo) pour le court terme.
        2. La probabilité de transition de Markov (calculée sur l'historique local)
           pour estimer la tendance climatique sous-jacente.
        La probabilité combinée pondère 70 % sur la prévision API et 30 % sur Markov.

        :param days_ahead: Nombre de jours d'avance (1 à 7).
        :param min_rainfall_mm: Seuil de précipitation pour considérer un jour comme humide.
        :return: Dictionnaire avec probabilités par jour et recommandations.
        """
        if self.forecast_data is None or self.forecast_data.empty:
            raise InsufficientData("Données de prévision indisponibles pour estimer la probabilité de pluie.")

        # Construction de la matrice de Markov depuis l'historique local
        try:
            markov = self.markov_transition(min_rainfall_mm)
            p_wet_if_wet = markov['p_wet_wet']
            p_wet_if_dry = markov['p_dry_wet']

            # État du dernier jour connu dans l'historique
            last_precip = self.rainfall_historical.iloc[-1] if not self.rainfall_historical.empty else 0.0
            current_p_wet = 1.0 if last_precip >= min_rainfall_mm else 0.0
        except Exception:
            # Si Markov échoue, on ne l'utilise pas
            markov = None
            p_wet_if_wet = 0.5
            p_wet_if_dry = 0.5
            current_p_wet = 0.0

        df = self.forecast_data.head(days_ahead)
        probs = {}
        max_prob = 0.0

        for i, (date, row) in enumerate(df.iterrows()):
            forecast_precip = row.get('precipitation', 0.0)

            # 1. Probabilité issue de la prévision API (heuristique sur la pluie prévue)
            prob_api = min(1.0, forecast_precip / (min_rainfall_mm * 5))
            if forecast_precip < min_rainfall_mm:
                prob_api *= 0.5

            # 2. Probabilité de Markov (probabilité conditionnelle d'un jour humide)
            prob_markov = current_p_wet * p_wet_if_wet + (1.0 - current_p_wet) * p_wet_if_dry

            # 3. Combinaison pondérée (API court terme = 70 %, Markov tendance = 30 %)
            prob_combined = 0.7 * prob_api + 0.3 * prob_markov

            # Mise à jour de l'état courant pour le jour suivant
            current_p_wet = prob_combined

            key = 'tomorrow' if i == 0 else f"{i + 1}_days"
            probs[key] = round(prob_combined, 2)
            if prob_combined > max_prob:
                max_prob = prob_combined

        # Recommandations agronomiques selon le risque maximal
        msg = f"{int(max_prob * 100)} % de chance de pluie dans les {days_ahead} prochains jours."
        if max_prob > 0.7:
            rec = "Risque de lessivage élevé. Repousser les traitements phytosanitaires."
        elif max_prob < 0.2:
            rec = "Conditions sèches attendues. Bon moment pour les traitements phyto."
        else:
            rec = "Vigilance recommandée pour les opérations au champ."

        return {
            **probs,
            'message': msg,
            'recommendation': rec,
        }

    def _get_severity_level(self, spi_value: float) -> str:
        """
        Interprète la valeur du SPI pour donner un niveau de sévérité.

        :param spi_value: Valeur de l'indice SPI.
        :return: Catégorie de sécheresse.
        """
        if spi_value > 1.0:
            return 'no_drought' # En réalité anormalement humide
        elif -1.0 <= spi_value <= 1.0:
            return 'mild'
        elif -1.5 <= spi_value < -1.0:
            return 'moderate'
        elif spi_value < -1.5:
            return 'severe'
        return 'unknown'

__init__ 

__init__(location: Location, rainfall_historical: Series, forecast_data: DataFrame)

Initialise les indicateurs de risque.

:param location: Instance de Location. :param rainfall_historical: Série historique de pluie. :param forecast_data: DataFrame de prévisions météorologiques.

Source code in kadi/weather/risk.py 
21
22
23
24
25
26
27
28
29
30
31
def __init__(self, location: Location, rainfall_historical: pd.Series, forecast_data: pd.DataFrame):
    """
    Initialise les indicateurs de risque.

    :param location: Instance de Location.
    :param rainfall_historical: Série historique de pluie.
    :param forecast_data: DataFrame de prévisions météorologiques.
    """
    self.location = location
    self.rainfall_historical = rainfall_historical
    self.forecast_data = forecast_data

drought_index 

drought_index(method: str = 'spi', window_months: int = 3) -> dict

Calcule l'indice de sécheresse avec la méthode spécifiée.

:param method: Méthode ('spi', 'markov', 'hurst', 'combined'). :param window_months: Fenêtre temporelle pour le calcul du SPI. :return: Dictionnaire avec les résultats de sécheresse.

Source code in kadi/weather/risk.py 
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
def drought_index(self, method: str = 'spi', window_months: int = 3) -> dict:
    """
    Calcule l'indice de sécheresse avec la méthode spécifiée.

    :param method: Méthode ('spi', 'markov', 'hurst', 'combined').
    :param window_months: Fenêtre temporelle pour le calcul du SPI.
    :return: Dictionnaire avec les résultats de sécheresse.
    """
    results = {}

    if method in ['spi', 'combined']:
        spi_val = self.spi(window_months)
        results[f'spi_{window_months}month'] = spi_val
        results['drought_severity'] = self._get_severity_level(spi_val)

    if method in ['markov', 'combined']:
        markov_res = self.markov_transition()
        results['markov_p_dry'] = markov_res.get('p_dry_dry', 0.0)

    if method in ['hurst', 'combined']:
        # L'exposant de Hurst nécessite une série assez longue (au moins qq années)
        hurst = self.hurst_exponent()
        results['hurst_exponent'] = round(hurst, 2)

    if method not in ['spi', 'markov', 'hurst', 'combined']:
        raise ValidationError(f"Méthode {method} non supportée pour l'indice de sécheresse.")

    return results

spi 

spi(window_months: int) -> float

Calcule le Standardized Precipitation Index (SPI) pour la période récente. Implémentation simplifiée pour le MVP.

:param window_months: Fenêtre d'accumulation en mois. :return: Valeur du SPI.

Source code in kadi/weather/risk.py 
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
def spi(self, window_months: int) -> float:
    """
    Calcule le Standardized Precipitation Index (SPI) pour la période récente.
    Implémentation simplifiée pour le MVP.

    :param window_months: Fenêtre d'accumulation en mois.
    :return: Valeur du SPI.
    """
    if self.rainfall_historical.empty:
        raise InsufficientData("Aucune donnée historique pour le calcul du SPI.")

    # Cumul de précipitation sur la fenêtre temporelle
    days = window_months * 30
    rolling_sum = self.rainfall_historical.rolling(window=days, min_periods=days//2).sum()
    rolling_sum = rolling_sum.dropna()

    if len(rolling_sum) < 30: # Pas assez de données pour faire une loi Gamma
        raise InsufficientData("Pas assez de jours de données pour ajuster le modèle SPI (minimum 30 requis).")

    # Ajustement d'une distribution Gamma simplifié
    # On évite les zéros pour la distribution Gamma
    valid_data = rolling_sum[rolling_sum > 0]
    if len(valid_data) == 0:
        return -3.0 # Extrême sécheresse

    # Dans un environnement de production, on utiliserait scipy.stats.gamma.fit
    # et une correction de probabilité mixte pour p(x=0). 
    # Ici on fait une approximation simple par score Z (loi normale) pour le mock
    current_val = rolling_sum.iloc[-1]
    mean_val = rolling_sum.mean()
    std_val = rolling_sum.std()

    if std_val == 0:
        return 0.0

    spi_val = (current_val - mean_val) / std_val
    return round(float(spi_val), 2)

markov_transition 

markov_transition(threshold_mm: float = 1.0) -> dict

Calcule les probabilités de transition de Markov (jour sec -> jour sec).

:param threshold_mm: Seuil pour considérer un jour comme humide. :return: Dictionnaire avec les probabilités p00, p01, p10, p11.

Source code in kadi/weather/risk.py 
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
def markov_transition(self, threshold_mm: float = 1.0) -> dict:
    """
    Calcule les probabilités de transition de Markov (jour sec -> jour sec).

    :param threshold_mm: Seuil pour considérer un jour comme humide.
    :return: Dictionnaire avec les probabilités p00, p01, p10, p11.
    """
    if self.rainfall_historical.empty:
        raise InsufficientData("Aucune donnée historique pour le calcul des probabilités de transition de Markov.")

    # États : 0 = sec, 1 = humide
    states = (self.rainfall_historical >= threshold_mm).astype(int)

    # Transitions
    transitions = pd.DataFrame({'current': states.iloc[:-1].values, 'next': states.iloc[1:].values})

    counts = transitions.groupby(['current', 'next']).size().unstack(fill_value=0)

    # On s'assure d'avoir la matrice 2x2
    for i in [0, 1]:
        if i not in counts.index:
            counts.loc[i] = [0, 0]
        for j in [0, 1]:
            if j not in counts.columns:
                counts[j] = 0

    # Probabilités
    p0 = counts.loc[0].sum()
    p1 = counts.loc[1].sum()

    p00 = counts.loc[0, 0] / p0 if p0 > 0 else 0
    p01 = counts.loc[0, 1] / p0 if p0 > 0 else 0
    p10 = counts.loc[1, 0] / p1 if p1 > 0 else 0
    p11 = counts.loc[1, 1] / p1 if p1 > 0 else 0

    return {
        'p_dry_dry': round(p00, 2),
        'p_dry_wet': round(p01, 2),
        'p_wet_dry': round(p10, 2),
        'p_wet_wet': round(p11, 2)
    }

hurst_exponent 

hurst_exponent(window: int = 1095) -> float

Calcule l'exposant de Hurst par la méthode de gamme rééchelonnée (R/S) multi-échelle.

L'algorithme segmente la série en sous-fenêtres de tailles croissantes, calcule le rapport R/S moyen pour chaque taille, puis estime H par régression log-log. Un exposant H > 0.5 indique une persistance climatique (mémoire longue).

:param window: Taille maximale de la fenêtre d'analyse (jours). :return: Exposant de Hurst H (compris entre 0.01 et 0.99).

Source code in kadi/weather/risk.py 
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
def hurst_exponent(self, window: int = 1095) -> float:
    """
    Calcule l'exposant de Hurst par la méthode de gamme rééchelonnée (R/S) multi-échelle.

    L'algorithme segmente la série en sous-fenêtres de tailles croissantes,
    calcule le rapport R/S moyen pour chaque taille, puis estime H par régression
    log-log. Un exposant H > 0.5 indique une persistance climatique (mémoire longue).

    :param window: Taille maximale de la fenêtre d'analyse (jours).
    :return: Exposant de Hurst H (compris entre 0.01 et 0.99).
    """
    if len(self.rainfall_historical) < 100:
        raise InsufficientData("Pas assez de données pour l'exposant de Hurst (minimum 100 jours requis).")

    data = self.rainfall_historical.values
    n_total = len(data)

    # Taille de la plus petite fenêtre (doit être assez grande pour un R/S stable)
    min_w = 10
    # Taille de la plus grande fenêtre (plafonnée à la moitié de la série)
    max_w = min(window, n_total // 2)

    window_sizes = []
    rs_means = []

    # Progression géométrique des tailles de fenêtre (facteur 1.5)
    w = min_w
    while w <= max_w:
        # Calcul du R/S moyen sur toutes les sous-séquences non chevauchantes
        num_segments = n_total // w
        rs_values_w = []

        for k in range(num_segments):
            segment = data[k * w: (k + 1) * w]
            seg_mean = np.mean(segment)
            centered = segment - seg_mean
            cum_dev = np.cumsum(centered)

            # Gamme (R) et écart-type (S) du segment
            r = np.max(cum_dev) - np.min(cum_dev)
            s = np.std(segment)

            if s > 0:
                rs_values_w.append(r / s)

        if rs_values_w:
            window_sizes.append(w)
            rs_means.append(np.mean(rs_values_w))

        w = max(w + 1, int(w * 1.5))

    # Minimum de 3 points pour une régression fiable
    if len(window_sizes) < 3:
        return 0.5

    # Régression linéaire dans l'espace log-log : log(R/S) = H * log(N) + c
    log_n = np.log(np.array(window_sizes, dtype=float))
    log_rs = np.log(np.array(rs_means, dtype=float))
    coeffs = np.polyfit(log_n, log_rs, 1)

    # H est la pente de la droite de régression
    h = float(coeffs[0])
    return float(np.clip(h, 0.01, 0.99))

rain_probability 

rain_probability(days_ahead: int = 1, min_rainfall_mm: float = 1.0) -> dict

Prévoit la probabilité de pluie pour les prochains jours.

Combine deux sources d'information pour plus de robustesse : 1. Les prévisions API (Open-Meteo) pour le court terme. 2. La probabilité de transition de Markov (calculée sur l'historique local) pour estimer la tendance climatique sous-jacente. La probabilité combinée pondère 70 % sur la prévision API et 30 % sur Markov.

:param days_ahead: Nombre de jours d'avance (1 à 7). :param min_rainfall_mm: Seuil de précipitation pour considérer un jour comme humide. :return: Dictionnaire avec probabilités par jour et recommandations.

Source code in kadi/weather/risk.py 
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
def rain_probability(self, days_ahead: int = 1, min_rainfall_mm: float = 1.0) -> dict:
    """
    Prévoit la probabilité de pluie pour les prochains jours.

    Combine deux sources d'information pour plus de robustesse :
    1. Les prévisions API (Open-Meteo) pour le court terme.
    2. La probabilité de transition de Markov (calculée sur l'historique local)
       pour estimer la tendance climatique sous-jacente.
    La probabilité combinée pondère 70 % sur la prévision API et 30 % sur Markov.

    :param days_ahead: Nombre de jours d'avance (1 à 7).
    :param min_rainfall_mm: Seuil de précipitation pour considérer un jour comme humide.
    :return: Dictionnaire avec probabilités par jour et recommandations.
    """
    if self.forecast_data is None or self.forecast_data.empty:
        raise InsufficientData("Données de prévision indisponibles pour estimer la probabilité de pluie.")

    # Construction de la matrice de Markov depuis l'historique local
    try:
        markov = self.markov_transition(min_rainfall_mm)
        p_wet_if_wet = markov['p_wet_wet']
        p_wet_if_dry = markov['p_dry_wet']

        # État du dernier jour connu dans l'historique
        last_precip = self.rainfall_historical.iloc[-1] if not self.rainfall_historical.empty else 0.0
        current_p_wet = 1.0 if last_precip >= min_rainfall_mm else 0.0
    except Exception:
        # Si Markov échoue, on ne l'utilise pas
        markov = None
        p_wet_if_wet = 0.5
        p_wet_if_dry = 0.5
        current_p_wet = 0.0

    df = self.forecast_data.head(days_ahead)
    probs = {}
    max_prob = 0.0

    for i, (date, row) in enumerate(df.iterrows()):
        forecast_precip = row.get('precipitation', 0.0)

        # 1. Probabilité issue de la prévision API (heuristique sur la pluie prévue)
        prob_api = min(1.0, forecast_precip / (min_rainfall_mm * 5))
        if forecast_precip < min_rainfall_mm:
            prob_api *= 0.5

        # 2. Probabilité de Markov (probabilité conditionnelle d'un jour humide)
        prob_markov = current_p_wet * p_wet_if_wet + (1.0 - current_p_wet) * p_wet_if_dry

        # 3. Combinaison pondérée (API court terme = 70 %, Markov tendance = 30 %)
        prob_combined = 0.7 * prob_api + 0.3 * prob_markov

        # Mise à jour de l'état courant pour le jour suivant
        current_p_wet = prob_combined

        key = 'tomorrow' if i == 0 else f"{i + 1}_days"
        probs[key] = round(prob_combined, 2)
        if prob_combined > max_prob:
            max_prob = prob_combined

    # Recommandations agronomiques selon le risque maximal
    msg = f"{int(max_prob * 100)} % de chance de pluie dans les {days_ahead} prochains jours."
    if max_prob > 0.7:
        rec = "Risque de lessivage élevé. Repousser les traitements phytosanitaires."
    elif max_prob < 0.2:
        rec = "Conditions sèches attendues. Bon moment pour les traitements phyto."
    else:
        rec = "Vigilance recommandée pour les opérations au champ."

    return {
        **probs,
        'message': msg,
        'recommendation': rec,
    }