X - Julie Carreau

THÉORIE DES VALEURS
EXTRÊMES
Julie Carreau IRD
HydroSciences Montpellier
www.pages-perso-julie-carreau.univ-montp2.fr
lundi 10 mars 2014
Événements météorologiques majeurs 2012-2013
Courrier international
lundi 10 mars 2014
OURAGAN SANDY : OCTOBRE 2012
Une des plus grosse tempête à frapper les États-Unis : évacuation, fermeture du
métro (4-5 jours), pénurie d’essence, ...
lundi 10 mars 2014
INONDATION DANS LE VAR : JUIN 2010
•
•
•
Facteurs météorologiques
mer chaude
relief montagneux
masses d’air chaud venant d’Afrique
lundi 10 mars 2014
•
•
•
400 mm de pluie en 24 h
équivalent de 5 mois de
précipitation
pas d’événement semblable
depuis 1827
INONDATIONS DU GARD EN SEPTEMBRE 2002
L’étude des pluies extrêmes est indispensable pour le bon dimensionnement des
ouvrages hydrauliques
Barrage écréteur de crues de la Rouvière :
• construit en 1978
• débit de contrôl de 1 m3/s
• niveau de retour de 50 ans
• hauteur de 18m
lundi 10 mars 2014
Septembre 2002 :
687 mm de pluie en 24h
Débit de 830 m3/s
TEMPÉRATURES EXTRÊMES
Montréal, Canada, 21 janvier 2013 : -27.3 0C
Froid extrême ?
lundi 10 mars 2014
TEMPÉRATURES EXTRÊMES
Marseille, record de froid le 12 février 1956 : -16.8 0C
Froid extrême ?
lundi 10 mars 2014
TEMPÉRATURES EXTRÊMES
Europe : Août 2003
•
•
•
plus de 20 000 morts en Europe
période la plus chaude depuis 500 ans
plusieurs pays européen : record de température
maximale
Impacts environnementaux
•
•
•
niveaux des cours d’eau et des lacs très bas
feux de forêts
fonte des glaciers : chute de pierres et de glace
lundi 10 mars 2014
TEMPÉRATURES EXTRÊMES
Europe : Août 2003
Vautard et al. 2007 Summertime European heat and drought waves induced by
wintertime Mediterranean rainfall deficit, GRL
58 ans d’observations météorologiques
Le risque de vague de chaleur extrême va vraisemblablement augmenter.
Des déficits de pluie en hiver dans le sud de l’Europe précède des étés chauds.
lundi 10 mars 2014
CHANGEMENT CLIMATIQUE
Rapport de la Banque Mondiale par le Postdam Institute for
Climate Impact Research and Climate Analytics, Novembre 2012
sans engagement, de nombreux scénarios prévoient
une augmentation de 4 ˚C d’ici la fin du siècle
augmentation du niveau de la mer de 0.5 à 1 m
Turn Down
the
Heat
Why a 4°C Warmer World
Must be Avoided
augmentation des températures élevées
extrêmes : impacts importants sur l’agriculture et
les écosystèmes
augmentation de l’intensité des cyclones
tropicaux
augmentation de l’aridité et des sécheresses
lundi 10 mars 2014
Théorie des Valeurs
Extrêmes
lundi 10 mars 2014
THÉORIE DES VALEURS EXTRÊMES
Applications en climat et environnement
précipitation : pluie, neige
débit des rivières
température
niveau de la mer
vitesse du vent
concentration des polluants
.....
Caractériser le comportement de très grandes valeurs
Caractériser le comportement de très petites valeurs
EVT fondée sur le comportement asymptotique des grandes valeurs :
maxima ou excédents
Les très grandes (petites) valeurs sont rares : modèles non-paramétriques
inadéquats
lundi 10 mars 2014
FONCTION DE RÉPARTITION EMPIRIQUE
Observations i.i.d. tirées d’une même variable aléatoire telle
que le niveau maximal annuel de la mer :
{X1 , . . . , Xn }
X∼F
Statistiques d’ordre : tri des observations en ordre croissant
X(1) ≤ · · · ≤ X(n)
{X(1) , . . . , X(n) }
Fonction de répartition empirique
n
�
1
p.s.
�
Fn (x) =
I{X(i) ≤x} −→ F (x)
n i=1
Loi forte des grands nombres
lundi 10 mars 2014
p.s.
�
Fn (x) −→ F (x)
lorsque
n −→ ∞
n
�
1
p.s.
�
Fn (x) =
I{X(i) ≤x} −→ F (x)
n i=1
i
�
Fn (X(i) ) = P (X ≤ X(i) ) =
n
i
P (X ≤ x) = ∀x ∈ [X(i) , X(i+1) )
n
P (X ≤ x) = 1 ∀x ≥ X(n)
1
P (X ≤ x) = 1 − p ∀p ≤
n
1
}
n
Queue légère
à décroissance exponentielle
lundi 10 mars 2014
Queue lourde
à décroissance polynômiale
Modèles pour les maxima de variables aléatoires
Maurice Fréchet
mathématicien français
1878 - 1973
Waloddi Weibull
ingénieur et
mathématicien suédois
1887 - 1979
Emil Julius Gumbel
mathématicien allemand
1891 - 1966
lundi 10 mars 2014
{X1 , . . . , Xn } variables aléatoires indépendantes avec fonction
Xi ∼ F
de répartition F
Par exemple : cumul journalier de précipitation
Objectif : modéliser les maxima
Mn = max{X1 , . . . , Xn }
Cumul maximal journalier de précipitation
Dérivons la fonction de répartition des maxima :
P (Mn ≤ x) = P (X1 ≤ x, . . . , Xn ≤ x)
= P (X1 ≤ x) × · · · × P (Xn ≤ x)
n
= {F (x)}
�n = {F� }n
Estimateur possible de la fonction de répartition des maxima : F
Une petite erreur dans l’estimation de F donnera une grande erreur dans Fn
➡ considérer des observations extrêmes et estimer Fn directement
lundi 10 mars 2014
Théorème central limite
(presque) toute somme de variables aléatoires indépendantes et identiquement
distribuées tend vers une variable aléatoire gaussienne
X 1 , X2 , . . . , X n
Soit
variables aléatoires indépendantes et identiquement
distribuées
Sn = X 1 + X 2 + · · · + X n
Si l’espérance et l’écart-type existent et
[Xi ] = µ < ∞ 0 < V ar[Xi ] = σ 2 < ∞
Alors
Sn d
√ −→ N (µ, σ 2 )
n
Les paramètres de position et de dispersion des Xi sont également ceux de la
distribution limite.
lundi 10 mars 2014
Convergence des cumuls de pluie
mensuel
journalier
trimestriel
Normale
semestriel
lundi 10 mars 2014
annuel
Maxima de données provenant d’une loi Log-Normale
• On retient le maximum un échantillon de taille 10 d’une loi Log-Normale
• On répète l’expérience 10 000 fois
échantillon de
taille 10
échantillon de
taille 100
échantillon de
taille 1000
• On recommence avec des échantillons de taille 100
• Et finalement avec des échantillons de taille 1000
Lorsque la taille de l’échantillon augmente, la distribution des maximas de
loi Log-Normale se rapproche de la loi de Gumbel.
lundi 10 mars 2014
Densité et fonction de répartition de la loi de Gumbel
Pour la Gumbel standard (centrée réduite) :
f (x) = exp{−x} exp{− exp{−x}}
F (x) = exp{− exp{−x}}
x−µ
Dans le cas général : x ←
σ
�
x−µ
f (x) = exp −
σ
lundi 10 mars 2014
�
�
�
x−µ
exp − exp −
σ
��
�
�
x−µ
F (x) = exp − exp −
σ
��
Loi de Gumbel
•
appropriée pour modéliser les maxima d’une variable suivant une loi LogNormale
par exemple : les maxima des cumuls journaliers de précipitation
•
mais aussi les maxima d’une variable suivant une loi Normale, Exponentielle,
Gamma, etc ....
Propriété des lois pour lesquelles la loi de Gumbel approxime bien les maxima :
la queue supérieure de la distribution décroît de façon exponentielle
1
1 − F (x) ≈ exp(−x) =
exp(x)
x↑∞
Pour de très grandes valeurs, la densité décroît très rapidement, à vitesse
exponentielle, vers zéro.
lundi 10 mars 2014
La loi de Gumbel est parfois inappropriée en présence de très
fortes valeurs
Maxima annuels de cumuls
journaliers de précipitation
lundi 10 mars 2014
Graphique quantile-quantile de
l’ajustement d’une loi de Gumbel
Maxima de données provenant d’une loi de Pareto généralisée
• On retient le maximum un échantillon de taille 10 d’une loi de Pareto généralisée
• On répète l’expérience 10 000 fois
• On recommence avec des échantillons de taille 100 et 1000
échantillon de
taille 10
échantillon de
taille 100
échantillon de
taille 1000
Lorsque la taille de l’échantillon augmente, la distribution des maximas de
loi de Pareto généralisée se rapproche de la loi de Fréchet.
lundi 10 mars 2014
Loi de Fréchet
modèle pour les données ayant de très fortes valeurs
✴ appropriée pour modéliser les maxima de variables suivant des lois de
Pareto, Student t, Cauchy
Propriété des lois pour lesquelles la loi de Fréchet approxime bien les maxima :
la queue supérieure de la distribution décroît de façon polynômiale
1 − F (x) ≈ x
−α
1
= α
x
x↑∞
On appelle α le paramètre de forme :
✦ plus α est petit, plus la décroissance est lente, plus des valeurs extrêmes sont
susceptibles de se produire
✦ plus α est grand, plus la décroissance est forte, moins des valeurs extrêmes vont
se présenter
lundi 10 mars 2014
Fonction de répartition de la loi de Fréchet
comportement en fonction du paramètre de forme α
Décroissance polynômiale
➡ α = 10
➡ α=5
➡ α=2
lundi 10 mars 2014
1 − F (x) ≈ x
−α
1
= α
x
x↑∞
Maxima de données provenant d’une loi Uniforme
échantillon de
taille 10
échantillon de
taille 100
échantillon de
taille 1000
Loi de Weibull : loi pour les données ayant une borne supérieure
• modèle réaliste lorsque les données ne peuvent pas dépasser une valeur maximale
• barrière physique qui borne les valeurs
• peut être le cas pour des données de température
✴ appropriée pour modéliser les maxima de variables suivant des lois
Uniforme, Bêta
lundi 10 mars 2014
Il y a trois et seulement trois distributions possibles pour les
maxima de variables aléatoires.
Décroissance exponentielle
Gumbel
Décroissance polynômiale
Fréchet
Densité bornée
Weibull
Ce sont les trois lois des valeurs extrêmes : elles correspondent aux trois
types de comportement de la queue de la distribution possibles
De façon analogue : la somme d’un grand nombre de variables aléatoires se
comportent comme une loi Normale (théorème central limite).
lundi 10 mars 2014
Théorème des types extrêmaux (Fisher-Tippett)
Soit Mn = max{X1 , X2 , . . . , Xn }
S’il existe des suites de constantes {an > 0} et {bn} telles que
P
�
Mn − bn
≤x
an
�
d
−→ G(x)
où G est une f.d.r. non dégénérée, alors G appartient à l’une des familles suivantes :
�
� �
���
I. Gumbel
x−b
G(x) = exp − exp −
a
II. Fréchet
� � � �
x −α
G(x) = exp −
a
x≥0
III. Weibull
� � x �α �
G(x) = 1 − exp −
a
−∞ < x < ∞
x≥0
pour des paramètres a > 0, b et α > 0
➡ Le comportement de la queue de distribution des Xi est lié au paramètre α de la
distribution limite.
lundi 10 mars 2014
MAX-STABILITÉ
Gumbel, Fréchet et Weibull sont les distributions de valeurs extrêmes.
Ce sont des distributions max-stables :
d
Mn = max{X1 , X2 , . . . , Xn } = an X + bn
�
�
Mn − bn
P
≤ x = F n (an x + bn ) = P (X ≤ x) = F (x)
an
Stabilité par l’addition :
d
X 1 + · · · + X n = a n X + bn
Quelle(s) loi(s) sont stables par l’addition ?
Relation entre les trois distributions max-stables :
X ∼ Fr´echet =⇒ Y = ln X ∼ Gumbel
X ∼ Weibull =⇒ Y = 1/X ∼ Fr´echet
X ∼ Weibull =⇒ Y = ln (1/X ) ∼ Gumbel
lundi 10 mars 2014
DOMAINE D’ATTRACTION
X est dans le domaine d’attraction maximum de la loi G si :
d
M
−b
n
n
P(
/an ≤ x) −→ G(x)
Queue légère
Gumbel distribution non bornée , tous les moments sont finis
Ex : Normale, Log-Normale, Exponentielle, Gamma
Caractéristique : queue supérieure décroît à vitesse exponentielle
Queue lourde
Fréchet distribution non bornée, certains moments ne sont
r
pas finis ; variance infinie si α > 2 et
[X ] < ∞ ↔ r < α
Ex : Cauchy, Student t, Pareto
−α
Caractéristique : queue supérieure décroît à vitesse polynômiale x
Queue finie
lundi 10 mars 2014
Weibull distribution bornée supérieurement
Ex : Uniforme, Bêta
Caractéristique : queue supérieure finie
DISTRIBUTION VALEUR-EXTRÊME GÉNÉRALISÉE (GEV)
G(x) =
�
� �
� x−µ ��−1/ξ �
exp − 1 + ξ σ
�
� x−µ ��
exp − exp − σ
position µ
échelle σ
forme
ξ = 1/α
Gumbel ξ → 0
lundi 10 mars 2014
si 1 + ξ(x − µ)/σ > 0 et ξ �= 0
avec x ∈
et ξ = 0
Condition sur les moments
[X r ] < ∞ ↔ ξ < 1/r
Fréchet ξ > 0
Weibull ξ < 0
EN PRATIQUE
Block-Maxima
On crée des blocs de longueur n à partir des observations avec n «grand»
Mn,1 , . . . , Mn,m
Choix fréquent : bloc = 1 an
Maxima annuel : cumul journalier, horaire ou mensuel
Maximum de vraisemblance
L(µ, σ, ξ; Mn,1 , . . . , Mn,m ) =
Où
∂G(x; µ, σ, ξ)
g(x; µ, σ, ξ) =
∂x
m
�
g(Mn,m ; µ, σ, ξ)
i=1
est la densité de la GEV
Minimisation de la log-vraisemblance négative
l(µ, σ, ξ; Mn,1 , . . . , Mn,m ) = − ln (L(µ, σ, ξ; Mn,1 , . . . , Mn,m ))
Sous contrainte :
1 + ξ(x − µ)/σ > 0
Cas à part :
➡ Bonnes propriétés de convergence en général et adaptabilité
lundi 10 mars 2014
ξ=0
FONCTION DE QUANTILES
xp est un quantile de niveau 1 − p ⇐⇒ G(xp ) = 1 − p
�
σ�
−ξ
xp = µ −
1 − {− log(1 − p)}
ξ
Niveau de retour xp de période de retour p :
niveau dépassé en moyenne une fois par T = 1/p années
T = 2 ←→ p = 0.5
p = 10%
Temps moyen entre les
dépassements : environ 10 ans
T = 5 ←→ p = 0.2
T = 10 ←→ p = 0.1
T = 20 ←→ p = 0.05
10 %
T = 50 ←→ p = 0.02
T = 100 ←→ p = 0.01
T = 1000 ←→ p = 0.001
lundi 10 mars 2014
90 %
Niveaux de retour
•
•
�
σ�
−ξ
xp = µ −
1 − {− log(1 − p)}
ξ
Pour des niveaux de quantiles modérés, le paramètre de position de la GEV a
le plus d’influence.
Pour des niveaux de quantiles élevés, c’est le paramètre de forme qui est
déterminant.
T = 10 T = 20 T = 50 T = 100 T = 1000
σ=1
Fréchet ξ = 0.2
Gumbel
ξ=0
µ = 100
µ = 10
ξ = −0.2
Weibull
➡ L’estimation du paramètre de forme dépend des observations extrêmes et
présente beaucoup de variance (sensibilité aux fortes observations).
lundi 10 mars 2014
Graphique Quantiles-Quantiles
Statistiques d’ordre : {X(1) , . . . , X(n) }
telles que
X(1) ≤ · · · ≤ X(n)
Fonction de répartition empirique
n
�
1
p.s.
F�n (x) =
I{X(i) ≤x} −→ F (x)
n i=1
Fréquences empiriques
� �n
�
�
i
1
n−1
=
,...,
,1
n i=1
n
n
Fréquences empiriques de Hazen
�
�n
�
�
i − 0.5
0.5
n − 0.5
=
,...,
n
n
n
i=1
➡ Correspondance entre les quantiles estimés par le modèle aux
fréquences empiriques et les statistiques d’ordre (quantiles empiriques)
lundi 10 mars 2014
Graphique Quantiles-Quantiles
��
�
�
��n
��
��n
�
�
i − 0.5
σ
�
−1
−ξ�
�
(i−0.5)
F
, X(i)
=
µ
�−
1 − {− log(
/n)}
, X(i)
n
ξ�
i=1
i=1
Si le modèle est juste, le graphique quantiles-quantiles
s’aligne sur la diagonale.
lundi 10 mars 2014
CUMUL DE PRÉCIPITATION JOURNALIER À MARSEILLE
1948
1948 - 1951
1948-2005 : 58 ans
Propriétés des précipitations
• intermittence
• variabilité saisonnière
• variabilité inter-annuelle
• valeurs extrêmes
lundi 10 mars 2014
MAX ANNUEL JOURNALIER À MARSEILLE
1948-2005 : 58 ans
Ajustement d’une loi GEV par maximum de vraisemblance
Paramètres estimés
location : 49.531 scale : 18.812 shape : 0.205
location : 49.461 scale : 18.182 shape : 0.133
Erreurs standard
location : 2.817 scale : 2.247 shape : 0.111
location : 2.756 scale : 2.129 shape : 0.111
lundi 10 mars 2014
sans le plus
grand max
DIAGNOSTIQUES
�
σ
��
−ξ�
x
�p = µ
�−
1 − {− log(1 − p)}
ξ�
Quelle est la période de retour de la plus forte pluie (environ 200 mm ) ?
➡ Incertitude plus grande pour les grands niveaux de retour (grands
quantiles) due à l’incertitude sur l’estimation du paramètre de forme
lundi 10 mars 2014
Comparaison avec un ajustement d’une loi de Gumbel ξ = 0
Paramètres estimés
Gumbel location : 51.79 (2.86) scale : 20.88 (2.25)
GEV
location : 49.53(2.82) scale : 18.81 (2.25) shape : 0.205 (0.111)
Log-vraisemblance nég. AIC (valeur élevée)
BIC (valeur élevée)
Gumbel 268.5355
2 x -268.5 - 2 x 2 = -541
2 x 268.5 - 2 x log(58) = -545
GEV
2 x -271.0 - 2 x 3 = -548
2 x 271,0 - 2 x log(58) = -554
270.9504
➡ Le paramètre de forme est-il significativement différent de zéro ?
➡ Quelle distribution est plus adéquate ?
lundi 10 mars 2014
NIVEAU DE LA MER À VENISE
Le acque alte
+1m
lundi 10 mars 2014
NIVEAU ANNUEL MAXIMAL
Ajustement d’une loi GEV par maximum de vraisemblance
Paramètres estimés
location : 111.092 (2.628) scale : 17.174 (1.803) shape : -0.077 (0.074)
lundi 10 mars 2014
DIAGNOSTIQUES
�
σ
��
−ξ�
x
�p = µ
�−
1 − {− log(1 − p)}
ξ�
Quel est le niveau de retour de 100 ans ?
lundi 10 mars 2014
TENDANCE ANNUELLE
Introduire une covariable pour le paramètre de position de la GEV : t année
µ(t) = µ0 + µ1 t
µ0 = 111.65 (2.25)
µ1 = 8.39 (2.07)
σ = 14.58 (1.58)
ξ = −0.027 (0.083)
Comment faire l’apprentissage des coefficients de la régression ?
Log-vraisemblance conditionnelle
l(µ0 , µ1 , σ, ξ|(t1 , x1 ), . . . , (tn , xn )) =
n
�
i=1
lundi 10 mars 2014
ln g(xi |µ = µ0 + µ + 1ti , σ, ξ)
Log-vraisemblance conditionnelle
Soit z, des variables indépendantes qui influencent la distribution de x
On fait l’hypothèse que :
µ(z) = h1 (z; β1 )
σ(z) = h2 (z; β2 )
ξ(z) = h3 (z; β3 )
où h1 (·; β1 ), h2 (·; β2 ), h3 (·; β3 ) sont des fonctions des variables z
avec les vecteurs de paramètres β1 , β1 , β3 à apprendre en maximisant :
l(β1 , β2 , β3 |(z 1 , x1 ), . . . , (z n , xn )) =
n
�
i=1
ln g(xi |µ = h1 (z i ; β1 ), σ = h2 (z i ; β2 ), ξ = h3 (z i ; β3 ))
➡ On fait donc l’hypothèse que les maxima sont distribués selon une loi GEV
conditionnellement à z
Niveaux de retour
�
�
σ(z)
−ξ(z)
1 − {− log(1 − p)}
xp (z) = µ(z) −
ξ(z)
Interprétation ?
lundi 10 mars 2014
Modèles pour les excédents de variables aléatoires
Approche «Peaks-over-Threshold» PoT développée par les hydrologues
Vilfredo Pareto, sociologue et économiste italien 1848 - 1923
Loi de Pareto : répartition des richesses
« 20 % de la population détient 80% des richesses »
lundi 10 mars 2014
Peaks-over-Threshold
Fixer un seuil et utiliser les excédents au-delà du seuil
Permet d’inclure plus d’observations dans l’estimation que l’approche des maxima
X∼F
V1
La distribution des excédents est donnée par :
Z1
Fu (y) = P {X ≤ u + y|X > u}
V3
Z5
V5
V4
u
Z9
Z6
V2
V7
Z12
V6
Z
Z3
11
Z
4
Z7
Z10
Z2
Z
8
∀0 < y < xF − u
Si F est connue alors :
F (u + y) − F (u)
Fu (y) =
1 − F (u)
Si on se sert d’une estimation de F pour estimer Fu (y) , on risque d’avoir des
problèmes d’extrapolation car il y a peu d’observations extrêmes.
Student Version of MATLAB
Distribution de Pareto généralisée GPD :
modélise les excédents au-delà d’un seuil
lundi 10 mars 2014
Distribution de Pareto généralisée
approximation pour la distribution des excédents
H(y) =
�
1−
1−
�
�
y −1/ξ
1 + ξ σ�
� y�
exp − σ�
si ξ > 0
si ξ = 0
ξ>0
ξ=0
ξ<0
Condition sur les moments
[X r ] < ∞ ↔ ξ < 1/r
Cas particuliers
�
lorsque ξ = 0 , H est la distribution exponentielle de paramètre σ
lorsque ξ < 0 , H est une distribution bornée supérieurement
lorsque ξ > 0 , H est une distribution de Pareto à queue lourde
➡ Le paramètre de forme définit trois classes de distributions tout comme les
domaines d’attraction maximale Gumbel, Weibull et Fréchet.
lundi 10 mars 2014
Théorème de Pickands-Balkema-de Haan
Soit X1, X2, ... Xn une suite de variables aléatoires indépendantes avec f.d.r. F et
considérons
M = max{X , X , . . . , X }
n
1
2
n
Soit X l’une des variables de la suite et supposons que, pour n grand, l’on ait :
P {Mn ≤ x} ≈ G(x)
où G est la f.d.r. de la GEV pour µ, σ et ξ. Alors, pour u suffisamment grand,
la distribution des excédents peut être approximée par
�
�
�
y −1/ξ
1 − 1 + ξ σ�
si ξ > 0
Fu (y) = P {X > u + y|X > u} ≈ H(y) =
� y�
si ξ = 0
1 − exp − σ�
définie sur {y : y > 0, (1 + ξy/σ�) > 0}
� = σ + ξ(u − µ)
où σ
➡ Si les maxima sont approximativement distribués selon la GEV, alors les
excédents sont approximativement distribués selon la GPD.
lundi 10 mars 2014
STABILITÉ PAR SEUILLAGE
H(y) =
�
1−
1−
�
�
y −1/ξ
1 + ξ σ�
� y�
exp − σ�
si ξ > 0
si ξ = 0
Pour une distribution donnée, si les maxima convergent vers la GEV, alors les
excédents convergent vers la GPD.
Les paramètres de la GPD sont déterminés par ceux de la GEV :
σ
� = σ + ξ(u − µ)
Le paramètre de forme est le même :
‣ ξ > 0 queue lourde (décroissance polynômiale)
‣ ξ = 0 queue légère (décroissance exponentielle)
‣ ξ < 0 queue finie (bornée)
La GPD possède la propriété de stabilité par l’opération seuil :
Y ∼ GP D =⇒ Y − u|Y > u ∼ GP D
lundi 10 mars 2014
∀u > 0
EN PRATIQUE
Peaks-over-Threshold
On définit les excédents à partir des observations avec un seuil u «grand»
Par exemple, on fixe k et on pose u = X(k)
Alors, les observations qui dépassent le seuil sont :
{X(1) , . . . , X(n) }
=⇒
{X(k+1) , . . . , X(n) }
Les excédents sont donnés par :
{Y1 = X(k+1) − u, . . . , Yn−k = X(n) − u}
Maximum de vraisemblance
nu
�
L(�
σ , ξ; Y1 , . . . , Ynu ) =
h(Yi − u; σ
�, ξ)
i=1
�, ξ) est la densité de la GPD
Où h(y; σ
Minimisation de la log-vraisemblance
négative
n
u
�
l(�
σ , ξ; Y1 , . . . , Ynu ) =
ln h(Yi − u; σ
�, ξ)
i=1
Sous contrainte : (1 + ξy/σ�) > 0
Cas à part :
➡ Bonnes propriétés de convergence en général et adaptabilité
lundi 10 mars 2014
ξ=0
CHOIX DU SEUIL
Dilemme biais-variance
‣ plus le seuil est élevé, plus l’approximation asymptotique de la
distribution des excédents par la GPD est juste et plus le biais diminue
‣ moins le seuil est élevé, plus le nombre d’excédents est grand pour
l’estimation des paramètres de la GPD et plus leur variance est petite
‣ dilemme similaire avec le choix de la taille du bloc dans l’approche des
maxima mais solution pragmatique bloc = année
➡ Adopter le seuil le plus bas tel que l’approximation asymptotique est à peu près
valide
lundi 10 mars 2014
Mean Residual Life Plot graphique de l’espérance excédentaire
Si Y ∼ GP D
et ξ < 1
alors
σ
[Y ] =
1−ξ
Supposons X − u0 |X > u0 ∼ GP D alors
σ u0
[X − u0 |X > u0 ] =
1−ξ
Aussi
∀u > u0
σu
σu0 + ξu
[X − u|X > u] =
=
1−ξ
1−ξ
Donc, pour u > u0 la moyenne est une fonction linéaire de u
➡ Tracer le graphique de l’espérance excédantaire et identifier la plus petite
valeur de u pour laquelle il devient linéaire :
��
lundi 10 mars 2014
1
u,
n−k
n
�
i=k+1
(X(i) − u)
��
Precipitation journalière positive à Marseille
graphique de l’espérance excédentaire
lundi 10 mars 2014
Stabilité des paramètres estimés par rapport au choix du seuil
��
lundi 10 mars 2014
u, ξ�u
��
{(u, σ
�u∗ )}
où
σu∗ = σu − ξu
NIVEAU DE RETOUR
Supposons que
�
P {X > x|X > u} = 1 + ξ
�
Il s’en suit :
P {X > x} = ζu 1 + ξ
Où
ζu = P {X > u}
�
�
x−u
σ
x−u
σ
��−1/ξ
��−1/ξ
∀x > u
∀x > u
Pour trouver le niveau qui sera dépassé une fois sur m, il faut résoudre :
�
�
��−1/ξ
xm − u
1
ζu 1 + ξ
=
σ
m
On obtient :
xm
σ
= u + [(mζu )ξ − 1]
ξ
xm > u
Pour obtenir le niveau dépassé en moyenne une fois par N ans avec ny le nombre
d’observations par année :
m=N ×n
y
k
Estimation de la probabilité d’excéder le seuil : ζ�u =
n
lundi 10 mars 2014
Precipitation journalière positive à Marseille
Fixons le seuil u = 10 mm : 1029 excédents soit moins de 5% des données
Estimateurs MLE des paramètres de la GPD
σ
� = 12.07(0.55) ξ = 0.13(0.033)
Estimateur MLE du paramètre de forme de la GEV (58 max annuels)
ξ� = 0.21 (0.11)
lundi 10 mars 2014
EXCÉDENTS CONDITIONNELLEMENT DÉPENDANTS
Log-vraisemblance conditionnelle
Soit z, des variables indépendantes qui influencent la distribution de x
On fait l’hypothèse que :
u(z) = h1 (z; β1 )
σ
�(z) = h2 (z; β2 )
ξ(z) = h3 (z; β3 )
où h1 (·; β1 ), h2 (·; β2 ), h3 (·; β3 ) sont des fonctions des variables z
avec les vecteurs de paramètres β1 , β1 , β3 à apprendre en maximisant :
l(β1 , β2 , β3 |(z 1 , x1 ), . . . , (z n , xn )) =
n
�
i=1
ln h(xi − u|u = h1 (z i ; β1 ), σ
� = h2 (z i ; β2 ), ξ = h3 (z i ; β3 ))
➡ On fait donc l’hypothèse que les maxima sont distribués selon une loi GPD
conditionnellement à z
lundi 10 mars 2014
Global Flood Risk under Climate Change
Nature Climate Change June 2013
‣ 11 modèles de climat de CMIP5
‣ Calcul du niveau de retour 100 ans pour les crues au 20ème siècle (1971-2000)
‣ Sous le scénario RCP85 au 21ème siècle (2071-2100), calcul de la période de
retour du niveau de retour 100 ans du 20ème siècle
‣ On considère la médiane des périodes de retour des 11 modèles.
lundi 10 mars 2014