Distributions Sinistres
Distributions Sinistres
Distributions Sinistres
de sinistres avec R
Vincent Goulet
Professeur titulaire | cole dactuariat | Universit Laval
Introduction
Ce document traite de modlisation statistique des montants de sinistres en
assurance de dommages. Il sagit de mes notes du cours Mathmatiques actuarielles
IARD I offert lcole dactuariat de lUniversit Laval auxquelles jai intgr du
code informatique ainsi que le recueil dexercices de Cossette et collab. (2009).
Les notes et les exercices de cet ouvrage ont atteint un niveau de maturit
enviable. Si cest le cas, cest parce que jai eu le privilge de btir sur du matriel
de grande qualit dvelopp au fil des ans par les professeurs Hlne Cossette,
Franois Dufresne, Jos Garrido, Michel Jacques et Jacques Rioux, ainsi que par le
charg de cours Mathieu Pigeon alors quil tudiait la matrise en actuariat au
milieu des annes 2000.
La matire est divise en sept chapitres. Le premier est un chapitre de rappels des notions de base en analyse, en probabilit et en statistique. Le chapitre 2
traite des fondements de la modlisation en assurance de dommages, en particulier le traitement mathmatique des franchises, limite suprieure et coassurance,
ainsi que de leffet de linflation sur la frquence et la svrit des sinistres. Les
aspects plus statistiques apparaissent au chapitre 3 avec la modlisation non paramtrique. Le chapitre 4 tudie les principales distributions utilises en assurance
de dommages et la cration de nouvelles distributions partir des lois usuelles.
Les chapitres 5 et 6 portent quant eux sur lestimation paramtrique et les tests
dadquation des modles. Enfin, le chapitre 7 propose une brve incursion dans la
modlisation des distributions de frquence des sinistres.
Les termes anglais ordinary deductible et franchise deductible utiliss dans les
textes de rfrence usuels en actuariat ont pos quelques soucis de traduction. Pour
le premier, jutilise lexpression franchise forfaitaire recommande par Bguin
(1990). Pour le second terme, beaucoup moins rpandu, jai opt pour lexpression
franchise atteinte suggre, entre autres, dans Charbonnier (2004).
La prsentation fait une large place aux considrations numriques auxquelles
tout actuaire praticien se retrouvera rapidement confront. Parce que jestime quil
sagit du meilleur outil aujourdhui disponible pour faire lanalyse, le traitement et la
modlisation de donnes de sinistres, les exemples et illustrations sont entirement
v
vi
Introduction
raliss dans lenvironnement statistique R (R Development Core Team, 2013). Aux
fonctionnalits de base du systme R pour nos fins sajoutent celles du package
actuar (Dutang et collab., 2008). Dailleurs, plusieurs fonctions du package ont t
dveloppes lorigine spcifiquement pour ce matriel.
Ltude du document implique quelques allers-retours entre le texte et les
sections de code informatique prsentes la fin de certains chapitre. Les sauts
vers ces sections sont clairement indiqus dans le texte par des mentions mises en
vidence par le symbole .
Chaque chapitre comporte des exercices. Les rponses de ceux-ci se trouvent
la fin de chacun des chapitres, alors que les solutions compltes sont regroupes
lannexe E. En consultation lectronique, le numro dun exercice est un hyperlien
vers sa solution, et vice versa.
On trouvera galement la fin de chaque chapitre (sauf le premier) une liste
non exhaustive dexercices proposs dans Klugman et collab. (2012a). Des solutions
de ces exercices sont offertes dans Klugman et collab. (2012b).
Lannexe A prsente la paramtrisation des lois de probabilit continues et
discrtes utilise dans les exercices. Linformation qui sy trouve est en plusieurs
points similaire celle des annexes A et B de Klugman et collab. (2008, 2012a),
mais la paramtrisation des lois est dans certains cas diffrente. Le lecteur est donc
fortement invit la consulter.
Lannexe B explique comment configurer R pour faciliter linstallation et ladministration de packages externes. Enfin, les annexes C et D offrent des tables de
quantiles des lois normale et khi carr.
Je remercie davance les lecteurs qui voudront bien me faire part de toute erreur
ou omission dans les notes, les exercices ou leurs solutions.
Vincent Goulet
Qubec, septembre 2013
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
1
2
6
9
16
17
19
29
33
35
40
42
43
45
49
.
.
.
.
.
.
.
.
.
.
55
56
56
60
62
63
68
73
74
76
80
viii
Inflation . . . . . . . . . . . . . . . . . . . .
Franchise, limite, coassurance et inflation
Code informatique . . . . . . . . . . . . .
Exercices . . . . . . . . . . . . . . . . . . .
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
81
87
87
90
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
97
98
100
111
125
128
133
141
.
.
.
.
.
.
.
149
149
153
155
156
167
167
173
.
.
.
.
.
.
.
.
181
182
184
185
191
203
212
216
222
.
.
.
.
.
233
233
238
247
248
250
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
ix
Modles de frquence
7.1 Fonction gnratrice des probabilits
7.2 Principales distributions discrtes . .
7.3 Estimation et tests . . . . . . . . . . . .
7.4 Familles (a, b, 0) et (a, b, 1) . . . . . . .
7.5 Exercices . . . . . . . . . . . . . . . . .
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
255
255
256
261
263
264
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
269
270
273
276
278
.
.
.
.
.
.
.
.
.
.
281
283
285
Bibliographie
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
.
287
287
300
313
333
345
364
371
379
Objectifs du chapitre
x Calculer la limite dune fonction.
x Calculer le dveloppement en srie de Taylor dune fonction et son intervalle de
convergence.
x Dfinir et calculer les principales quantits lies une variable alatoire (densit,
x
x
x
x
x
x
fonction de rpartition, esprance, densits conjointe et marginale, fonctions gnratrices des moments et des probabilits).
Dterminer la loi de la transformation dune variable alatoire.
Dfinir, calculer et dterminer la distribution des principales quantits lies un
chantillon alatoire (moyenne, variance, statistiques dordre).
Dfinir les principales caractristiques dun estimateur (biais, efficacit, convergence, exhaustivit).
Dfinir, calculer et tablir les caractristiques des estimateurs ponctuels des moments et du maximum de vraisemblance.
Dfinir et dterminer un intervalle de confiance pour une quantit.
tablir un test dhypothses pour une quantit donne. Dfinir et calculer les principales quantits lies au test (hypothses nulle et alternative, erreur de type I, erreur
de type II, niveau de confiance, valeur p).
1.1
Limites
1.1.1
Thorme du sandwich
xc
xc
alors
lim f (x) = L.
xc
xc
Alors, pour tout > 0, il existe > 0 tel que pour tout x, lingalit c < x < c +
implique
L < g (x) < L +
et
L < h(x) < L + .
Ces deux ingalits, combines g (x) f (x) h(x), donnent
L < g (x) f (x) h(x) < L +
L < f (x) < L +
< f (x) L < .
On peut donc en conclure que, pour tout x, lingalit c < x < c + implique
| f (x) L| < .
1.1. Limites
Comme on peut prendre aussi prs de 0 quon le dsire, on peut en conclure que
lim f (x) = L.
xc +
Les cas o la limite est approche par la gauche et par les deux cts se traitent de
la mme faon.
Exemple 1.1. Soit la relation
1
x2
x sin(x)
<
<1
6
2 2 cos(x)
x sin(x)
x0 2 2 cos(x)
x0
et donc
x2
=1
6
x sin(x)
= 1.
x0 2 2 cos(x)
lim
1.1.2
Rgle de lHpital
La rgle de lHpital permet de raliser des calculs de limites lorsque des formes
0
indtermines (
, 0 , 0 , etc.) sont rencontres.
Thorme 1.2. Soit f (x) et g (x) deux fonctions diffrentiables sur un intervalle
ouvert I , sauf ventuellement en un point x = c. On suppose que, sur I , on a
g (x) 6= 0 pour x 6= c. Si, dans ces conditions, limxc f (x) = limxc g (x) = 0, ou encore
l i m xc f (x) = et limxc g (x) = , alors la limite de
f (x)
g (x)
peut tre calcule laide de la rgle
lim
xc
f (x)
f 0 (x)
= lim 0
g (x) xc g (x)
0.85
0.65
0.75
f(x)
0.95
0.5
0.0
f(x)
0.5
1.0
1.1. Limites
0.0
0.5
1.0
1.5
2.0
2.5
3.0
x0+
x0+
x0
x0
= 0.
0.5
0.0
f(x)
1.0
0.0
0.5
1.0
1.5
2.0
1.2
Dveloppement en sries
a n (x c)n = a 0 + a 1 (x c) + a 2 (x c)2 + . . .
n=0
o x est une variable et o les a i et c sont des constantes. Le but de cette section est
dexprimer certaines fonctions sous la forme de sries de puissances. Ces dernires
sont utilises, entre autres, par les calculatrices pour valuer des fonctions trigonomtriques. Elles peuvent galement tre utilises pour calculer numriquement
certaines intgrales dfinies difficiles valuer autrement.
1.2.1
Srie de Taylor
Soit une fonction f (x) que lon souhaite dvelopper autour du point x = c. La
fonction f (x) doit tre diffrentiable indfiniment en x = c. On pose alors
f (x) = a 0 + a 1 (x c) + a 2 (x c)2 + a 3 (x c)3 + . . .
On calcule ensuite les drives :
f 0 (x) = a 1 + 2a 2 (x c) + 3a 3 (x c)2 + 4a 4 (x c)3 + . . .
f 00 (x) = 2a 2 + 6a 3 (x c) + 12a 4 (x c)2 + . . .
f 000 (x) = 6a 3 + 24a 4 (x c) + . . . ,
f 00 (c)
f 000 (c)
f (n) (c)
X
(x c)n .
n!
n=0
+...
2
6
24
X
(x 1)n
(1)n+1
=
.
n
n=1
ln(x) = (x 1)
1.2.2
Intervalle de convergence
Lorsque que lon attribue la variable x dune srie une valeur numrique, il
devient alors possible de dterminer si la srie diverge ou converge. On appelle
intervalle de convergence dune srie de puissances lensemble des valeurs de x
pour lesquelles la srie converge.
Pour dterminer cet intervalle, on value
u n+1
|u n+1 |
= |H |
lim
= lim
n |u n |
n u n
et on utilise comme critre
xn = 1 + x + x2 + x3 + x4 + . . .
n=0
On a
u n+1
|u n+1 |
lim
= lim
n |u n |
n u n
n+1
x
= lim n
n x
= lim |x|
n
= |x|.
Il faut vrifier pour x = 1 et x = 1. On a
11+11+...
qui diverge et
1+1+1+1+...
qui diverge galement. Lintervalle de convergence est donc (1, 1).
On conclut cette section avec trois thormes (prsents sans dmonstration)
concernant les dveloppements en sries de puissances.
Thorme 1.3. Une srie de puissances est une fonction continue lintrieur de son
intervalle de convergence.
Thorme 1.4. Une srie de puissance peut tre drive terme par terme lintrieur
de son intervalle de convergence.
Thorme 1.5. Un srie de puissances peut tre intgre terme par terme lintrieur
de son intervalle de convergence.
Exemple 1.6. On a vu que lintervalle de convergence de la srie
f (x) = 1 + x + x 2 + x 3 + x 4 + . . .
Exprience
Variable alatoire
Lancer de deux ds
Lancer de 25 pices de monnaie
crasement dun avion
1/2
0
x3 d x + . . .
x2 x3
+
+ . . . |1/2
0
2
3
1 1 1
+....
= + +
2 8 24
=x+
1.3
Variable alatoire
Soit X , une variable alatoire. Il sagit dune fonction qui, chaque lment
de lespace fondamental dune exprience alatoire, associe un nombre rel. Le
tableau 1.1 prsente quelques exemples de variables alatoires.
1.3.1
Fonction de rpartition
0.0
0.2
0.4
F(x)
0.6
0.8
1.0
10
1.3.2
Fonction de survie
11
0.0
0.2
0.4
S(x)
0.6
0.8
1.0
1.3.3
0.4
0.0
0.2
f(x)
0.6
12
1.3.4
Taux dincidence
Le taux dincidence (hazard rate en anglais) h(x) (ou (x)) dune variable alatoire X est dfini comme le ratio de la densit et de la fonction de survie, cest--dire
f (x)
S(x)
S 0 (x)
=
S(x)
d
=
ln S(x).
dx
h(x) =
Le taux dincidence porte aussi parfois le nom de force de mortalit et est alors
gnralement not (x).
La fonction de survie peut tre retrouve partir du taux dincidence :
S(x) = e
Rx
0
h(y) d y
13
= ln S(x).
1.3.5
Esprance et moments
g (x)Pr[X = x].
x=
Cette dernire expression peut tre vue comme une moyenne pondre des valeurs
de g (x). Si E [|g (X )|] = , on dit que E [g (X )] nexiste pas.
Le thorme suivant nonce quelques proprits importantes de loprateur
esprance.
Thorme 1.6. Soit X une variable alatoire, a, b et c des constantes et g 1 (x) et g 2 (x)
deux fonctions quelconques. On a alors les proprits suivantes :
i)
ii)
iii)
iv)
On nomme, lorsquil existe, moment dordre k de la variable alatoire X lesprance de cette dernire leve la puissance k. On a
0k = E [X k ]
Z
=
x k f (x) d x
14
x k Pr[X = x].
x=
ou
=
(x )k Pr[X = x].
E [(X )3 ]
Var[X ]3/2
3
3/2
2
15
2 (X ) =
E [(X )4 ] 4
= 2.
Var[X ]2
2
1.3.6
Pour une variable alatoire X , la fonction gnratrice des moments est donne
par
M X (t ) = E [e t X ]
lorsque lesprance existe (ce qui nest pas le cas pour de nombreuses distributions.
La fonction porte ce nom parce que sa k e drive value en 0 est gale au k e
moment de la variable alatoire. En effet,
M X(k) (t ) =
dk
E [e t X ]
d"t k
#
dk tX
e
=E
dtk
= E [X k e t X ],
do M X(k) (0) = E [X k ].
16
1.3.7
X
=
t x Pr[X = x]
x=0
X
x
=
e t x
x!
x=0
(t )x
X
= e
x=0 x!
= e e t
= e (t 1) .
1.4
Soient X et Y deux variables alatoires. Leur densit conjointe est note f X Y (x, y)
dans le cas continu et Pr[X = x, Y = y] dans le cas discret.
Dans le cas continu, les densits marginales sont
Z
f X (x) =
f (x, y) d y
17
et
f Y (y) =
f (x, y) d x.
f X Y (x, y)
f Y (y)
f Y |X (y|x) =
f X Y (x, y)
.
f X (x)
et
f Y |X (y|x) f X (x)
f Y (y)
Enfin, la loi des probabilits totales permet de calculer une densit marginale
partir dune distribution conditionnelle :
Z
f X (x) =
f X |Y (x|y) f Y (y) d y.
1.5
e (x)
0
(2)n/2 ||1/2
(x)/2
Concepts de convergence
E [g (X )]
.
r
18
E [g (X )] =
g (x) f X (x) d x.
= r Pr[g (X ) r ],
do
Pr[g (X ) r ]
1.5.1
E [g (X )]
.
r
Convergence en probabilit
Il faut donc qu partir dun certain point dans la suite, la probabilit que la distance
entre les deux variables alatoires soit aussi proche de 0 quon le dsire soit de 1.
Il est noter quil nest pas ncessaire que les variables alatoires X 1 , X 2 , . . .
soient indpendantes et/ou identiquement distribues.
1.5.2
Ce type de convergence est plus fort que la convergence en probabilit (cest-dire que la convergence presque certaine implique la convergence en probabilit,
mais que linverse nest pas ncessairement vrai). Une suite de variables alatoires
X 1 , X 2 , . . . converge presque certainement vers une variable alatoire X si, pour tout
> 0, on a
h
i
Pr lim |X n X | < = 1.
n
n
1X
Xi .
n i =1
1.5.3
Convergence en distribution
1.6
Cette section porte sur les transformations de variables alatoires, soit des
fonctions dune ou plusieurs variables alatoires. En termes mathmatiques, tant
donn la distribution conjointe des variables alatoires X 1 , . . . , X n , on cherche
dterminer la fonction de probabilit ou de densit de la variable alatoire Y =
u(X 1 , . . . , X n ).
Voici quelques exemples de transformations frquemment rencontres :
Y = X2
X E [X ]
Y = p
Var[X ]
X1 + + Xn
Y =
n
Y = F X (X ).
Il existe trois techniques principales pour dterminer la distribution de la transformation Y :
1. la technique de la fonction de rpartition ;
2. la technique du changement de variable ;
3. la technique de la fonction gnratrice des moments.
19
20
1.6.1
Cest la technique la plus simple, mais pas toujours la plus facile demploi. En
effet, la fonction de rpartition de certaines lois de probabilit est complique,
voire nexiste pas sous forme explicite (penser ici aux lois normale et gamma, par
exemple).
Lide consiste simplement calculer la fonction de rpartition de la transformation avec
F Y (y) = Pr[Y y]
= Pr[u(X 1 , . . . , X n ) y],
puis calculer la densit (ou la fonction de probabilit) par diffrenciation :
f Y (y) = F Y0 (y).
Il importe de noter que le domaine de dfinition de la transformation nest pas
ncessairement le mme que celui de ou des variables alatoires de dpart.
Exemple 1.8. Soit
f X (x) =
(
6x(1 x), 0 < x < 1
0,
ailleurs,
= 3y 2/3 2y,
0 < y < 1,
et donc
f Y (y) = F Y0 (y)
(
2y 1/3 2,
=
0,
0<y <1
ailleurs.
21
Exemple 1.9. Soit X une variable alatoire continue quelconque avec fonction de
rpartition F X (x) et Y = aX + b, o a et b sont des constantes relles. On a
F Y (y) = Pr[a X + b y]
y b
= FX
a
et, par consquent,
1
y b
f Y (y) = f X
.
a
a
f Y (y) =
f X (y) + f X (y),
y >0
0,
ailleurs.
y > 0.
22
Exemple 1.11. Soit Y = X 1 + X 2 , o X 1 et X 2 sont deux variables alatoires indpendantes chacune distribue uniformment sur lintervalle (0, 1). On a donc
f X 1 X 2 (x 1 , x 2 ) = f X 1 (x 1 ) = f X 2 (x 2 ) = 1,
F Y (y) =
=
d x1 d x2
1 2
y .
2
F Y (y) = (y 1)(1) +
yx 2
d x1 d x2
y1 0
1
= 1 (2 y)2 .
2
4. Enfin, si y 2, clairement F Y (y) = 1.
Par consquent, on a
F Y (y) =
0,
1 y 2,
y 0
0<y 1
1,
1
2
2 (2 y) ,
1<y <2
y >2
et
f Y (y) =
0,
y,
2 y,
0,
y 0
0<y 1
1<y <2
y > 2.
1.6.2
23
16 ,
4,
16
6
Pr[X = x] = 16
,
16 ,
1
16 ,
x =0
x =1
x =2
x =3
x = 4.
6
16
10
16
2
Pr[Y = 1] = Pr[X = 0] + Pr[X = 4] = .
16
Pr[Y = 1] = Pr[X = 1] + Pr[X = 3] =
24
Cas continu
Le cas continu est beaucoup plus dlicat. On considre toujours la transformation Y = u(X ) avec les hypothses suivantes :
y/2 1 y/2
f Y (y) = f X (e
) e
2
1
= e y/2
2
1/2 11 y/2
=
y e
,
y > 0,
(1)
w 0 (y) =
Par consquent,
f Y (y) = f X (F (y))
f X (x)
= f X (x)
f X (x)
1
= 1,
0 < y < 1,
p y 1/2
f Y (y) = f X ( y)
2
1
= p e y/2
2
( 12 )1/2 1/2 y/2
=
y e
, y > 0,
( 21 )
soit Y 2 (1).
25
26
1.6.3
Il sagit simplement ici de gnraliser les concepts tudis la section prcdente des transformations impliquant plusieurs variables alatoires.
Lide reste la mme sinon quil faut sassurer que la transformation compte
autant de nouvelles variables que danciennes. Ainsi, si lon part de la distribution
conjointe de deux variables alatoires X 1 et X 2 , il faudra trouver la distribution
conjointe de deux nouvelles variables alatoires Y1 = u 1 (X 1 , X 2 ) et Y2 = u 2 (X 1 , X 2 ).
Plus souvent quautrement, seule la distribution de la variable Y1 est dintrt.
Il suffit alors de dfinir Y2 comme une variable muette (dummy), par exemple
Y2 = X 2 .
Cas discret
Si la transformation Y1 = u 1 (X 1 , X 2 ) et Y2 = u 2 (X 1 , X 2 ) est bijective, alors simplement
Pr[Y1 = y 1 , Y2 = y 2 ] = Pr[u 1 (X 1 , X 2 ] = y 1 , u 2 (X 1 , X 2 ) = y 2 )
= Pr[X 1 = w 1 (y 1 , y 2 ], X 2 = w 2 (y 1 , y 2 )),
o
w 1 (y 1 , y 2 ) = u 11 (x 1 , x 2 )
w 2 (y 1 , y 2 ) = u 21 (x 1 , x 2 ).
Les fonctions de probabilit marginales sont alors obtenues en sommant :
Pr[Y1 = y 1 ] =
Pr[Y2 = y 2 ] =
Pr[Y1 = y 1 , Y2 = y 2 ]
y 2 =
Pr[Y1 = y 1 , Y2 = y 2 ].
y 1 =
X 1 = Y1 Y2
X 2 = Y2 .
27
11
22 e 2
(y 1 y 2 )!
y2!
et donc
y y 2
Pr[Y1 = y 1 ] =
y1 1
X
1
y 2 =0
22 e (1 +2 )
(y 1 y 2 )! y 2 !
(1 +2 )
y1!
y1
X
y1!
y y y
11 2 22
(y
y
)!
y
!
1
2
2
y 2 =0
!
y
1
X y 1 y 1 y 2 y 2
1
2
y 2 =0 y 2
e (1 +2 )
y1!
e (1 +2 )
(1 + 2 ) y 1 ,
y1!
soit Y Poisson(1 + 2 ).
Cas continu
On gnralise le thorme 1.10 afin de trouver la distribution conjointe (et
ventuellement les distributions marginales) de Y1 = u 1 (X 1 , X 2 ) et Y2 = u 2 (X 1 , X 2 ).
On suppose que
28
x 1
y 2
x 2
y
x 1
y
1
J =
x 2
y
1
1
1
x 1 x 2 e x1 x2 ,
()() 1
x 1 > 0, x 2 > 0.
De plus,
Y1 = X 1 + X 2
X1
Y2 =
X1 + X2
X 1 = Y1 Y2
X 2 = Y1 (1 Y2 )
et donc
x 1
= y1
y 2
x 2
= y 1 ,
y 2
x 1
= y2
y 1
x 2
= 1 y2
y 1
do le Jacobien de la transformation est
y2
J =
1 y
y 1
= y 1 .
y 1
Le domaine de Y1 est R+ alors que celui de Y2 est limit lintervalle (0, 1). Par
le thorme 1.11,
f Y1 Y2 (y 1 , y 2 ) = f X 1 X 2 (y 1 y 2 , y 1 (1 y 2 )) |y 1 |
1
=
y 1 (y 1 y 2 )1 (y 1 (1 y 2 ))1 e y 1 y 2 y 1 (1y 2 )
()()
1
+1 1
=
y1
y 2 (1 y 2 )1 e y 1
()()
= g (y 1 ) h(y 2 )
29
f Y1 (y 1 ) =
f Y1 Y2 (y 1 , y 2 ) d y 2
Z 1
1
+1 y 1
= y1
e
y 21 (1 y 2 )1 d y 2
0 ()()
1
+1 y 1
=
y
e ,
y1 > 0
( + ) 1
0
et
f Y2 (y 1 , y 2 ) =
f Y1 Y2 (y 1 , y 2 )
f Y1 (y 1 , y 2 )
( + ) 1
=
y
(1 y 2 )1 ,
()() 2
0 < y 2 < 1.
1.6.4
Cette technique savre tout spcialement puissante pour dterminer la distribution (marginale) dune combinaison linaire de variables alatoires indpendantes. La technique repose sur le thorme suivant.
Thorme 1.12. Soit X 1 , . . . , X n des variables alatoires indpendantes et Y = X 1 +
+ X n . Alors
n
Y
M Y (t ) =
M X i (t ).
i =1
1.7
chantillon alatoire
x statistique ;
x chantillon alatoire.
30
x On lance la pice 100 fois. Les valeurs x 1 , . . . , x 100 forment un chantillon alatoire
dune population Bernoulli.
X 1 + + X 100
100
x 1 + + x 100
.
100
1.7.1
Quelques dfinitions
Dfinition 1.2 (Statistique). Une statistique est une fonction dune ou plusieurs
variables alatoires qui ne dpend pas de paramtres inconnus.
Remarques.
31
1. Par exemple,
Y =
1
(X 1 + + X n )
n
n
1X
Xi
n i =1
n
1X
(X i X )2
n i =1
n
1 X
(X i X )2 .
n 1 i =1
32
n
1X
2,
(x i x)
n i =1
des ralisations de ces deux variables alatoires obtenues partir dune ralisation de lchantillon alatoire.
1.7.2
2
.
n
Z=
Pr[Z z] (z).
1.7.3
X N (, 2 /n)
ii)
nS 2 /2 2 (n 1)
iii)
X et S 2 sont indpendantes.
1.8
h0
Pour une distribution continue, cela se rsume au point x tel que F (x) = p.
33
34
1.8.1
Distribution du maximum
1.8.2
Distribution du minimum
1.8.3
x1 xn
n!
(F (x))k1
(k 1)!(r k 1)!(n r )!
(F (y) F (x))r k1 (1 F (y))nr f (x) f (y),
f X (k) (x) =
1.8.4
xy
n!
(F (x))k1 (1 F (y))nk f (x)
(k 1)!(n k)!
1.9
1.9.1
Biais
35
36
Enfin, lerreur quadratique moyenne (MSE) dun estimateur est dfinie comme
= E [
)2 ] = b()2 + Var[].
MSE()
Exemple 1.18. Soit X 1 , . . . , X n un chantillon alatoire issu dune loi exponentielle
translate avec fonction de densit de probabilit
f (x; ) = e (x) ,
x > .
On cherche un estimateur sans biais pour . Essayons dabord avec X . Pour faciliter
les calculs, dfinissons toutefois dabord la variable alatoire X = + Y , o Y
Exponentielle(1). Ainsi,
E [X ] = + E [Y ]
= +1
et par consquent
E [ X ] = E [X ] = + 1.
La moyenne de lchantillon est donc un estimateur biais de . Par contre, X 1
est un estimateur sans biais.
Intuitivement, X (1) = min(X 1 , . . . , X n ) devrait aussi constituer un estimateur
raisonnable de . Or, par le rsultat de la sous-section 1.8.2, on a f Y(1) (y) = ne n y et
donc Y(1) Exponentielle(n) E [Y(1) ] = n 1 . Par consquent,
E [X (1) ] = + E [Y(1) ]
1
= + .
n
Le minimum est donc un estimateur biais, mais asymptotiquement sans biais du
paramtre dune loi exponentielle translate.
nS 2
= n 1
2
et donc
E [S 2 ] =
n 1 2
,
n
S2 =
garantie que g () sera sans biais pour g (). Par exemple, S2 est sans biais pour
2 , mais S est biais pour .
2. Comme on a pu le voir en partie lexemple 1.18, il existe une multitude destimateurs sans biais dun paramtre . Quel estimateur choisir alors ?
3. Il peut exister un meilleur estimateur biais, cest--dire un estimateur qui, bien
que biais, a une plus petite variance.
37
38
1.9.2
Efficacit
1
2 i .
h
ln f (X )
nE
h
2 i
ln f (X )
1.9.3
39
Convergence
Var[]
,
2
= 0 limn Pr[|
| ] = 0.
do limn Var[]
1.9.4
Exhaustivit
On peut omettre cette section lors dune premire lecture. Son contenu est plac ici titre de
rfrence.
40
( + 1) 1
x
()
= x 1 ,
On dmontre que
Qn
i =1 X i
f (x 1 , . . . , x n ; ) =
n
Y
x 1
i =1
n
= (x 1 x 2 x n )1
= n (x 1 x 2 x n )
1
x1 x2 xn
= u(y; ) v(x 1 , . . . , x n ),
o y = i =1 x i , u(y; ) = n y et v(x 1 , . . . , x n ) = (x 1 x 2 x n )1 . Par le thorme de
Q
factorisation, Y = ni=1 X i est donc une statistique exhaustive pour .
Qn
1.10
Estimation : mthodes
On se tourne maintenant vers ltude de mthodes pour dvelopper des estimateurs du ou des paramtres dune distribution.
1.10.1
k = 1, . . . , r.
1.10.2
41
x > 0.
n
Y
e xi
i =1
= n e
Pn
i =1 x i
n
X
i =1
do
n
n X
d
l () =
xi .
d
i =1
xi ,
42
i =1 x i
1
.
x
1
.
X
On remarque que :
P
= X Gamma(n, n), do
1. Puisque ni=1 X i Gamma(n, ), alors
n
n
=
n 1 n 1
n 2 2
=
Var[]
.
(n 1)2 (n 2)
=
E []
Par consquent,
est asymptotiquement sans biais ;
x
= 0, donc
est convergent.
x limn Var[]
En fait, pour tous les cas intressants lEMV est un estimateur convergent.
2. Il est facile de vrifier que la borne de RaoCramr est 2 /n. Par consquent
=
Efficacit()
n3
(n 1)2 (n 2)
1.11
1.11.1
Invariance
1.11.2
Exhaustivit
1.11.3
Efficacit
1.12
43
44
X n
Pr z /2 < p
< z /2 = 1
n(1 )
45
Pr
X
z /2
n
(1 )
X
< < + z /2
n
n
(1 )
= 1 .
n
Afin de pouvoir obtenir simplement des valeurs numriques pour les bornes
infrieure et suprieure de lintervalle de confiance pour , on remplace ce paramtre dans les racines carres par son estimateur du maximum de vraisemblance,
= X /n. Un intervalle de confiance approximatif de niveau 1 est donc
z /2
(1
+ z /2
,
n
(1
,.
n
Dans cet exemple o n = 400 est grand, on peut utiliser lapproximation cidessus. On a donc 1 = 0,95 z 0,025 = 1,96 et x = 136, do
s
0,34 1,96
(0,34)(0,66)
< < 0,34 + 1,96
400
(0,34)(0,66)
400
m
0,294 < < 0,386
avec approximativement une probabilit de 95 %. Exprim autrement, la probabilit de dcs est de 0,34 0,046, 19 fois sur 20.
1.13
Tests dhypothses
Reprenons le contexte de lexemple 1.22, savoir que 136 des 400 personnes
atteintes dune maladie grave dcdent dans lanne suivant le diagnostic de la
maladie, soit une proportion de = 0,34. Un nouveau mdicament que lon espre
efficace contre la maladie est dcouvert. Afin dvaluer si le mdicament est efficace ou non, celui-ci est administr aux nouveaux cas de la maladie. On voudra
par la suite comparer la nouvelle proportion de dcs avec celle observe avant
lutilisation du mdicament.
En termes plus mathmatiques, on voudra tester si < 0,34 chez les gens ayant
pris le nouveau mdicament. On formule donc deux hypothses statistiques :
H0 : = 0,34 (le mdicament est inefficace)
H1 : < 0,34 (le mdicament est efficace).
46
(1 )
N ,
.
n
Ainsi, pour n = 400,
0,3 0,34
Pr Z p
= 0,046,
(0,34)(0.66)/400
o Z N (0, 1).
Inversement, si > 0,30 on ne rejettera pas H0 . Le mdicament sera jug inefficace, alors quil peut ltre. Le fait de ne pas rejeter lhypothse H0 alors quelle est
fausse est appel une erreur de type II. Cette probabilit est note :
= Pr[faire une erreur de type II]
= Pr[ne pas rejeter H0 alors quelle est fausse].
47
X
> 0,30; , < 0,34 .
n
0,28 0,34
0,28; = 0,34] Pr Z p
Pr[
(0,34)(0,66)/400
= 0,0057.
Cette valeur est appele la valeur p du test, et elle reprsente le seuil de signification
maximal auquel on peut rejeter H0 tant donn la valeur de la statistique obtenue
avec un chantillon. Autrement dit, tant donn = 0,28, on serait prt rejeter H0
avec un niveau de confiance dau moins 99,43 %.
La valeur p doit-elle tre grande ou petite ? Cela dpend si lon souhaite ou non
rejeter lhypothse nulle ! Par exemple, si
H0 : ajustement dun modle des donnes est bon
H1 : ajustement nest pas bon,
alors on cherchera avoir une valeur p trs grande, de faon ne pas rejeter H0 .
Par contre, dans une situation o
H0 : droite de rgression nexplique pas les donnes
H1 : droite de rgression explique les donnes,
on souhaitera habituellement rejeter H0 , alors la valeur p devrait tre petite pour
une statistique donne.
48
1.14. Exercices
1.14
49
Exercices
1.1 On a lingalit
1 x 2 1 cos(x) 1
<
<
2 24
x2
2
vraie pour toutes valeurs de x prs de 0. Calculer
lim
x0
1 cos(x)
x2
x0 ln(x + 1)
i = (i 2 )(i 2 )
= (1)(1)
=1
5
i =i
i 6 = 1
50
x2 x3 x4
+
+
+...,
2!
3!
4!
1
1 + e x
< x < .
,
(x + )+1
1.14. Exercices
51
1 1
,
1 + x2
< x < .
0 < x < ,
0 < x < 3.
x > 0.
52
1 < x < 1.
2x
,
c2
0 < x < c.
moyenne de .
1.24 Soit X 1 , . . . , X n , un chantillon alatoire dune population avec moyenne et
variance 2 .
P
a) Dmontrer que lestimateur T (X ) = ni=1 a i X i est un estimateur sans biais
P
de si ni=1 a i = 1.
b) On nomme les estimateurs de la forme en a) des estimateurs sans biais
linaires. Parmi ceux-ci, trouver celui avec la plus petite variance.
1.25 Soit X 1 , . . . , X n un chantillon alatoire dune distribution avec moyenne et
P
variance 2 . Dmontrer que n 1 ni=1 (X i )2 est un estimateur sans biais de
2 .
1.26 Soit X , une observation dune population dont la densit est
|x|
f (x; ) =
(1 )1|x| , x = 1, 0, 1; 0 1.
2
Soit lestimateur
T (X ) =
(
2,
0,
x =1
ailleurs.
1.14. Exercices
53
X
X
1
n
n
n
1 (1)/
x
,
Rponses
1.1
1
2
1.2 1
1.3 e
54
y >3
1
27 , 0 < y < 27
Gamma( 12 , 12 2 )
1.15 f Y (y) =
1.16
1.18 2
1.19 9/5
1.20 2(ct )2 (ct 2t c e t c + 1)
1.21 332
1.22 a) Gamma(2 500, 50) b) 0,682722 c) 0,6826
1.23 Biais : 500 ; MSE : 250 750
1.24 a) X
1.28 1
1.29 = 0,6122, = 0,5102.
1.30 a) Bta(1/, 1)
x Connatre les dfinitions et leffet sur les donnes de sinistres des modifications
x
x
x
x
x
56
2.1
Processus de modlisation
Un modle est une description plus ou moins simplifie dune ralit. Le modles actuariels peuvent tre de plusieurs types :
2.2
Donnes
2.2. Donnes
57
Intervalle
Nombre dobservations
(0, 25]
(25, 50]
(50, 100]
(100, 150]
(150, 250]
(250, 500]
30
31
57
42
65
84
x une franchise qui limine les petites rclamations et permet dviter une frquence trop leve de sinistres entranant beaucoup de frais administratifs ;
58
2.2. Donnes
59
0.0000
0.0010
Density
0.0020
0.0030
Histogram of x
200
400
600
800
1000
500
1000
1500
y
2000
2500
Density
0.0000
0.0010
0.0012
0.0008
Density
0.0000
0.0004
0.0008
0.0004
0.0000
Density
Histogram of y
0.0020
Histogram of y
0.0012
Histogram of y
500
1000
y
1500
2000
200
400
600
800
1000
60
Y =
X,
X 700
700,
X > 700.
0 y < 700
f X (y),
f Y (y) = 1 F X (700), y = 700
0,
y > 700.
,
0 y < 700
e
=
1 e 700 ,
0,
y = 700
y > 700.
1
.
637,51
2.3
Les donnes sont rarement entres directement au clavier dans R. Pour les
importer dune source externe, on utilisera prfrablement la fonction scan. Dans
61
Density
100
200
300
400
500
600
700
son utilisation la plus simple, celle-ci lit un fichier texte dans lequel se trouvent les
donnes, une par ligne ou spares par des blancs, en omettant les lignes dbutant
par le caractre spcifi dans largument comment.char. Le rsultat est plac dans
un vecteur.
Par exemple, pour importer les donnes du fichier donnees.txt reprsent la
figure 2.4, on fera simplement 1 :
> x <- scan("donnees.txt", comment.char = "#")
> x
[1] 141
16
46
40 351
62
noter que ces fonctions retournent un data frame. Lindicage et le traitement des
data frames tant plus lent que celui des vecteurs simples, il est recommand de
transformer ses donnes dans ce dernier format, si possible.
Pour plus dinformations sur limportation de donnes dans R, consulter le
manuel R Data Import/Export, distribu avec R et disponible dans le site CRAN.
2.4
Consulter le code informatique de la section 2.13 pour des illustrations de la cration, de limportation et du traitement des donnes
individuelles et groupes avec R.
x Y S : la variable alatoire du montant pay par sinistre par lassureur. Sil ny pas
de paiement pour un sinistre, la variable alatoire prend une valeur de 0. Elle est
donc dfinie mme si aucun paiement nest effectu.
Exemple 2.2. Soit un portefeuille de polices dassurance habitation des franchises
de 1 000 $. Ces polices gnrent les montants de sinistres suivants au cours dune
2.5. Franchises
63
Variable alatoire
X
YP
YS
Valeurs
750
900
1 200
200
200
5 000
4 000
4 000
10 000
9 000
9 000
priode : 1 200, 5 000, 750, 10 000, 900. Le tableau 2.1 prsente les valeurs des
variables alatoires X , YP et YS .
2.5
Franchises
2.5.1
Franchise forfaitaire
La variable du montant pay par paiement est une troncature par le bas et une
translation de la variable alatoire du montant du sinistre, conditionnellement ce
que ce sinistre dpasse la franchise. On a donc
n
Y P = X d, X > d
ou, pour faire ressortir plus clairement la nature conditionnelle de la variable
alatoire Y P :
Y P = X d |X > d .
Par consquent,
F Y P (x) = Pr[Y P x]
= Pr[X d x|X > d ]
0.00
0.20
0.10
0.00
0.10
fX(x)
fYP(x)
0.20
0.30
0.30
64
10
10
Pr[d < X x + d ]
Pr[X > d ]
0,
x <0
= F X (x + d ) F X (d )
, x 0
1 F X (d )
=
et, en drivant,
0,
x <0
f Y P (x) = f X (x + d )
, x 0.
1 F X (d )
2.5. Franchises
65
donc
S
Y =
(
0,
X d,
X d
X >d
(
F X (d ),
x =0
f X (x + d ), x > 0.
La figure 2.6 montre leffet dune franchise forfaitaire sur la densit de Y S . Cette
fois, suite lintroduction de la franchise : 1) laire sous la densit dorigine entre 0
et la franchise est maintenant entirement concentre en 0 ; 2) la densit au-dessus
de la franchise est translate vers la gauche ; 3) les valeurs de la densit ne sont pas
autrement modifies puisque laire totale sous la densit est toujours gale 1.
Remarque. On a la relation suivante entre les variables alatoires du montant pay
par sinistre et du montant pay par paiement :
Y P = Y S |Y S > 0.
2.5.2
Franchise atteinte
0.20
0.10
0.00
0.00
0.10
fX(x)
fYS(x)
0.20
0.30
0.30
66
10
10
0,
x <d
F Y P (x) = F X (x) F X (d )
, x d
1 F X (d )
et
f Y P (x) =
0,
x <d
f X (x)
, x d.
1 F X (d )
67
0.00
0.20
0.10
0.00
0.10
fX(x)
fYP(x)
0.20
0.30
0.30
2.5. Franchises
10
10
et
f Y S (x) =
(
F X (d ), x = 0
f X (x),
x > d.
On a donc la relation
Y P = Y S |Y S > d .
La figure 2.8 montre leffet dune franchise forfaitaire sur la densit de Y S .
La figure 2.9 runit pour fins de comparaison les densits de Y P et Y S pour les
deux types de franchise.
2.5.3
Fonction coverage
La fonction coverage du package actuar permet dobtenir facilement les fonctions de rpartition et fonctions de densit sous diverses combinaisons de modifications de couverture. On donne en argument la fonction
x f X () (si ncessaire) ;
x F X () (si ncessaire) ;
x la franchise d ;
0.20
0.10
0.00
0.00
0.10
fX(x)
fYS(x)
0.20
0.30
0.30
68
10
10
2.6
Esprance limite
69
fYS(x)
0.00
0.00
0.10
0.10
fYP(x)
0.20
0.20
0.30
0.30
10
10
fYS(x)
0.00
0.00
0.10
0.10
fYP(x)
0.20
0.20
0.30
0.30
10
10
70
E [X ; u] =
=
Z0 u
0
x f (x) d x + u
f (x) d x
u
De plus,
uZ x
E [X ; u] =
Z0 u Z0 u
E [X ] =
(1 F (x)) d x,
F (x) = 1
, x >0
+x
et fonction de densit
f (x) =
,
( + x)+1
x >0
71
est, si 6= 1,
u
E [X ; u] =
(1 F (x)) d x
Z u
dx
=
+x
0
0
1
=
1 x +
u
=
1
.
1
u +
0
ln x
,
=
1
ln x
n 1 ln x 2 o
1 1
=p
exp
.
2
2 x
f X (x) =
On a donc
x
E [X ; x] =
=
Z0 x
0
ln y
d y + x(1 F (x)).
Z ln x
n 1 z 2 o
1
dz
=p
e z exp
2
2
Z ln x
n z 2 2z + 2 22 z o
1
=p
exp
dz
22
2
72
Z ln x
n 22 + 4 o
n 1 z 2 2 o
1
=p
exp
exp
dz
22
2
Z ln x
2
z 2
dz
= e + /2
ln x 2
+2 /2
.
=e
Z
0
ln x
x f (y) d y = p
2
Ainsi, on obtient
E [X ; x] = e
+2 /2
ln x 2
ln x
+x 1
,
x > 0.
Exemple 2.5. Soit Y S le montant pay par sinistre pour une franchise forfaitaire
de d . Puisque
(
F X (d ),
y =0
f Y S (y) =
f X (y + d ), y > 0,
on a
Z
y f X (y + d ) d y
E [Y S ] = 0 F X (d ) +
0
Z
=
(x d ) f X (x) d x
Zd
=
x f X (y) d x d (1 F X (d ))
Z
x f X (x) d x
x f X (x) d x d (1 F X (d ))
= E [X ] E [X ; d ].
On peut aussi astucieusement procder ainsi : on a
S
Y =
(
0,
X d,
X d
X >d
73
Y =
X X,
X d
X d,
X >d
= X (X d ),
do E [Y S ] = E [X ] E [X ; d ].
Il est laiss en exercice de vrifier que
E [Y p ] =
2.7
E [X ] E [X ; d ]
.
1 F X (d )
Esprance rsiduelle
do
e(x) =
E [X ] E [X ; x]
.
1 F X (x)
74
2.8
Le rapport dlimination de perte (loss elimination ratio, LER) est une mesure de
lefficacit pour lassureur de lintroduction ou de laugmentation dune franchise
(et, par extension, de toute autre modification apporte un contrat).
Le LER est le rapport entre la diminution des cots esprs suite une modification apporte au contrat et les cots esprs sans celle-ci :
LER =
LER =
Exemple 2.7. Les sinistres suivants sont survenus pour un contrat dassurance :
{5, 10, 15, 20}.
Si lassureur introduit une franchise de 10 dans le contrat, il conomisera les montants suivants sur chacun des sinistres ci-dessus :
{5, 10, 10, 10}.
Son rapport dlimination de perte sera donc
5 + 10 + 10 + 10
5 + 10 + 15 + 20
= 0,70.
LER =
75
LER = 1
E [X ; d ]
.
E [X ]
E [X ; d ] =
(1 F (x)) d x
e x d x
1 e u
=
.
1 000(1 e 0,5 )
= 0,3935.
1 000
76
2.9
Limite suprieure
(
Pr[X x], 0 x < u
1,
x u,
(
F X (x), 0 x < u
1,
x u,
et
f Y (x) =
f X (x),
0x <u
1 F X (u), x = u
0,
x > u.
77
0.20
fY(x)
0.6
0.10
0.4
0.2
FY(x)
0.8
1.0
0.30
0.0
0.00
10
10
F IG . 2.10: Fonction de rpartition et fonction de densit de probabilit de la distribution du montant pay par lassureur avec une limite suprieure
Y P = min(X , u) d |X > d
(
X d, d < X u
=
u d , X > u.
On dit que la variable alatoire Y P est gale la variable alatoire X tronque vers
le bas d , translate gauche et censure par le haut u. La fonction de rpartition
78
F
(x
+
d ) F X (d )
X
, 0 x < u d
1 F X (d )
=
1 F X (d )
,
x u d
1 F X (d )
F X (x + d ) F X (d ) , 0 x < u d
1 F X (d )
=
1,
x u d.
=
0,
f X (x + d )
1 F X (d )
f Y P (x) =
1 F X (u)
1 F X (d )
0,
x <0
0 x < u d
x = u d
x > u d.
La figure 2.11 montre leffet combin dune franchise forfaitaire et dune limite
sur la fonction de rpartition de Y P .
La variable alatoire du montant pay par sinistre, Y S , est
Y S = max(min(X , u) d , 0)
X d
0,
=
X d, d < X u
u d , X > u.
On a donc
F Y S (x) =
F X (d ),
x =0
F X (x + d ), 0 < x < u d
1,
x u d,
0.8
0.6
0.0
0.2
0.4
FYP(x)
0.6
0.4
0.0
0.2
FX(x)
0.8
1.0
79
1.0
10
10
et
F X (d ),
f (x + d ),
X
f Y S (x) =
1 F X (u),
0,
x =0
0 < x < u d
x = u d
x > u d.
La figure 2.12 montre leffet combin dune franchise forfaitaire et dune limite
sur la fonction de rpartition de Y S .
Exemple 2.9. Soient X la variable alatoire reprsentant le montant des sinistres,
u = 100, une limite suprieure ajout au contrat et Y , la variable alatoire reprsentant le montant pay par lassureur. On a donc
Y = min(X , 100).
Calculer le rapport dlimination de perte si les sinistres suivants sont survenus :
{50, 100, 150, 200}.
0.8
0.6
0.0
0.2
0.4
FYS(x)
0.6
0.4
0.0
0.2
FX(x)
0.8
1.0
1.0
80
10
10
On a
0 + 0 + 50 + 100
50 + 100 + 150 + 200
50 + 100 + 150 + 200 50 + 100 + 100 + 100
=
LER =
50+100+150+200
4
50+100+150+200
4
50+100+100+100
4
50+100+150+200
4
E [X ; 100]
E [X ]
= 0,30.
= 1
Lajout dune limite suprieure de 100 au contrat permet donc lassureur dliminer 30 % de la charge financire lie ce contrat.
2.10
Coassurance
0 1,
0.2
0.0
0.1
fY(x)
0.2
0.0
0.1
fX(x)
0.3
81
0.3
2.11. Inflation
10
10
do
f Y (x) =
x
1
fX
.
2.11
Inflation
Les modles pour les montants de sinistres sont construits partir dobservations survenues un certain nombre dannes dans le pass. En gnral, il est
ncessaire dajuster les modles obtenus afin de les amener au niveau actuel de
lexprience de pertes. De plus, on peut tre intress faire une projection du
modle afin de reflter les sinistres anticips dans une priode future.
fZ(x)
0.002
0.000
0.001
fX(x)
200
400
600
800
0.003
82
200
400
600
800
2.11.1
Effet de linflation
x
1
fX
.
1+r
1+r
2.11. Inflation
83
trois ans, trouver la probabilit que, dans trois ans, le montant dun sinistre soit
infrieur ou gal 500 $.
On a X Pareto(2, 100) et Z = (1,05)3 X = 1,157625X . Il est laiss en exercice de
vrifier que Z Pareto(2, 115,7625). Par consquent,
Pr[Z 500] = 1
115,7625
500 + 115,7625
= 0,9647.
On pourrait aussi procder sans trouver la distribution de Z :
500
Pr[Z 500] = Pr X
1,157625
2
100
= 1
100 + 431,9188
= 0,9647.
On peut vrifier le rsultat ci-dessus avec R et actuar :
> F <- coverage(cdf = ppareto, inflation = 1.05^3 - 1)
> F(500, 2, 100)
[1] 0.9646565
Exemple 2.11. Un assureur modlise le montant dun sinistre par une Pareto de
paramtres = 2 et = k. Le contrat comporte une franchise de 2k. Calculer le
rapport dlimination de perte avant et aprs une inflation de 100 %.
Avant linflation, on a
LER =
E [X ; 2k]
E [X ]
21
k
k
1
21
2k+k
2k
2k + k
2
= .
3
k
21
84
LER =
2.11.2
Linflation augmente les cots des sinistres. Si les contrats comportent une
franchise et que celle-ci demeure inchange, alors leffet de linflation est amplifi :
V P = Z d |Z > d ,
o Z = (1 + r )X . On a donc
V P = (1 + r )X d |X >
d
1+r
2.11. Inflation
85
et
E [Y p ] =
E [X ] E [X ; d ]
.
1 F X (d )
E [V ] =
d
])
(1 + r )(E [X ] E [X ; 1+r
d
1 F X ( 1+r
)
Tournons-nous maintenant vers le volet frquence des sinistres. Soit p, lesprance de la frquence des sinistres. On a alors que p(1 F X (d )) est lesprance
de la frquence des paiements (la frquence des sinistres pour lesquels lassureur
effectue un paiement). Aprs inflation, lesprance de la frquence des paiements
devient
d
p(1 F Z (d )) = p(1 F X ( 1+r
)).
Par consquent, si la prime pure avant inflation est
= p(1 F X (d ))E [Y P ]
E [X ] E [X ; d ]
= p(1 F X (d ))
1 F X (d )
= p(E [X ] E [X ; d ]),
elle devient, aprs linflation,
= p(1 F Z (d ))E [V P ]
d
= p(1 F X ( 1+r
))
d
(1 + r )(E [X ] E [X ; 1+r
])
d
1 F X ( 1+r
)
d
= p(1 + r )(E [X ] E [X ; 1+r
]).
86
on a que
(1 + r )p(E [X ] E [X ; d ])
et donc
(1 + r ).
En rsum, si un contrat dassurance comporte une franchise et que celle-ci
demeure inchange alors que les sinistres subissent une certaine inflation, la prime
pure subit une augmentation suprieure au taux dinflation. Cest leffet de levier
de linflation.
On remarquera en terminant que si le montant de la franchise d est ajust pour
devenir (1 + r )d , la prime pure subit une augmentation gale linflation.
2.11.3
Limposition dune limite suprieure dans un contrat dassurance protge lassureur contre les effets de linflation. En effet, si la limite demeure constante alors
que le cot des sinistres augmente, le rsultat net quivaut abaisser la limite.
Intuitivement, leffet de levier de linflation pour une limite suprieure sera lavantage de lassureur. Il ne reste qu tablir que la baisse des cots sera suprieure
linflation.
Clairement, la frquence des sinistres nest pas affecte par linflation dans un
contrat ne comportant quune limite suprieure u. Par contre, le montant pay par
paiement par lassureur passe de
Y = min(X , d )
V = min(Z , u)
= min((1 + r )X , u)
u
= (1 + r ) min(X , 1+r
).
Par consquent,
u
E [V ] = (1 + r )E [X ; 1+r
]
(1 + r )E [X ; u],
87
do
E [V ] (1 + r )E [Y ].
Une limite suprieure a donc un effet modrateur sur la hausse de la prime
pure lorsque les sinistres augmentent avec linflation. Par contre, si le montant de
la limite u est ajust pour devenir (1 + r )u, la prime pure subit une augmentation
gale linflation.
2.12
X < 1+r
0,
=
((1 + r )X d ),
(u d ),
d
1+r
X <
u
1+r
u
1+r
2.13
Code informatique
###
### DONNES INDIVIDUELLES ET GROUPES
###
## Le vecteur simple (ou atomique) cr avec la fonction 'c' est encore
88
89
# changement de frquences
# changement de bornes
###
### FONCTION 'coverage'
###
## La fonction provient du package actuar. Il faut le charger en
## mmoire si ce n'tait pas fait prcdemment.
library(actuar)
## La fonction 'coverage' retourne une fonction pour calculer la f.r.
## ou la f.d.p. de la variable alatoire du montant pay par paiement
## ou pay par sinistre modifie par une franchise (ordinaire ou
## atteindre), une limite, de la coassurance ou l'inflation.
##
## On examine les quatre combinaisons possibles de variables
## alatoires et de franchise. On suppose que la variable alatoire X
## a une distribution gamma. Le montant de la franchise est d = 1.
##
## 1. Montant par paiement et franchise ordinaire (cas par dfaut)
f <- coverage(dgamma, pgamma, deductible = 1) # cration de l'objet
mode(f)
# c'est une fonction
f
# code de la fonction
f(0, 3)
# calcul en x = 0
f(5, 3)
# calcul en x = 5
dgamma(5 + 1, 3)/pgamma(1, 3, lower = FALSE) # idem
curve(dgamma(x, 3), from = 0, to = 10, ylim = c(0, 0.3)) # originale
curve(f(x, 3), from = 0.1, col = "blue", add = TRUE)
# modifie
## 2. Montant par sinistre et franchise ordinaire
f <- coverage(dgamma, pgamma, deductible = 1, per.loss =
f(0, 3)
pgamma(1, 3)
curve(dgamma(x, 3), from = 0, to = 10, ylim = c(0, 0.3))
curve(f(x, 3), from = 0.1, col = "blue", add = TRUE)
points(0, f(0, 3), pch = 16, col = "blue")
TRUE)
# masse 0
# idem
# originale
# modifie
# masse 0
90
2.14
Exercices
2.1 Les montants suivants reprsentent les cots associs aux rparations automobiles de 12 contrats :
{579, 110, 842, 213, 98, 445, 1 332, 162, 131, 276, 312, 482}.
Les contrats prsentent une franchise forfaitaire de 250 $. Calculer le rapport
dlimination de perte (LER) de lassureur.
2.2 Les montants suivants reprsentent les cots associs des accidents automobiles pour huit contrats :
{86 000, 123 000, 423 000, 43 000, 213 000, 28 000, 52 000, 178 000}.
Les contrats prsentent une limite suprieure de 100 000 $. Calculer le rapport
dlimination de perte de lassureur.
2.3 Pour un portefeuille dont le montant dun sinistre obit une loi exponentielle
de paramtre 0,02, trouver le rapport dlimination de perte dcoulant de
lintroduction des limites de couvertures suivantes.
a) Une franchise atteinte de 10.
b) Une franchise forfaitaire de 10.
2.14. Exercices
91
2.4 On suppose que le montant dun sinistre obit une distribution gamma de
paramtres = 4 et = 0,1. Un assureur a sign un trait avec un rassureur o
ce dernier sengage payer lexcdent de 100 sur chacun des sinistres. Trouver
le rapport dlimination de perte de lassureur.
2.5 Dans un groupe dassurs, les sinistres suivants sont survenus :
{20, 50, 80, 80, 80, 85, 90, 110, 150, 240, 360, 400}.
Trouver le rapport dlimination de perte de lassureur si celui-ci a instaur
une franchise forfaitaire de 70 et sil limite ses paiements 200.
2.6 Soit X , la variable alatoire reprsentant le montant dun sinistre. On sait que
E [X ] = 2 000, que E [X ; 30 000] = 1 640,79 et que le rapport dlimination de
perte de lassureur pour un contrat avec une franchise forfaitaire de 100 est de
0,0465. Trouver le rapport dlimination de perte de lassureur pour un contrat
avec une franchise forfaitaire de 100 et une limite suprieure de 30 000.
2.7 Soit X , une variable alatoire reprsentant le montant dun sinistre tel que
f X (x) = e 2x +
e x
,
2
x > 0.
a) Trouver E [X ; d ].
b) Soit N , une variable alatoire reprsentant la frquence des sinistres. Calculer la prime pure (frquence moyenne multiplie par la svrit moyenne)
pour une franchise de d = 0,25 et une frquence moyenne de un sinistre
tous les 10 ans.
c) Si on observe un taux dinflation de 5 %, que devient la prime pure ?
2.8 On suppose que le montant dun sinistre obit une loi Pareto de paramtres
= 1,5 et = 2 500.
a) Calculer le montant moyen des sinistres pay par un assureur pour un
contrat de rassurance avec une rtention de 50 000.
b) Trouver le rapport dlimination de perte pour le rassureur si la rtention
est de 100 000.
2.9 Soit Y P la variable alatoire du montant pay par paiement pour un contrat
dassurance avec une franchise forfaitaire de d et X est la variable alatoire du
montant dun sinistre. Dmontrer que
E [Y P ] =
E [X ] E [X ; d ]
,
1 F X (d )
92
2.14. Exercices
93
x le montant dun sinistre pour lanne 1990 obit une loi Pareto de paramtres = 1,5 et = 1 500 ;
Nombre
Montant moyen
0 2 500
2 500 7 500
7 500 12 500
12 500 17 500
17 500 22 500
22 500 32 500
32 500 47 500
47 500 67 500
67 500 87 500
87 500 125 000
125 000 225 000
225 000 300 000
300 000
41
48
24
18
15
14
16
12
6
11
5
4
3
1 389
4 661
9 991
15 482
20 232
26 616
40 278
56 414
74 985
106 851
184 735
264 025
300 000
Pour modliser les donnes, on utilise une distribution log-normale de paramtres et 2 . laide dune technique destimation quelconque, on trouve
= 1,596.
que = 9,356 et
a) Estimer le montant pay espr.
b) Estimer le pourcentage de changement dans le montant pay par paiement espr si lon observe une inflation de 10 % des sinistres.
94
x > 0.
X <d
0,
S
Y = (X d ), d X < u
(u d ), X u.
2.14. Exercices
a) Dmontrer que E [Y S ] = (E [X ; u] E [X ; d ]).
b) Trouver Var[Y S ].
c) Trouver lexpression gnrale de lesprance du montant pay par sinistre
la suite dune inflation de 100r %.
2.18 Soient Y S , la variable alatoire du montant pay par sinistre, X , la variable
alatoire du montant dun sinistre, d une franchise forfaitaire et u, une limite
suprieure. Dmontrer la relation E [Y S ] = E [X ; u] E [X ; d ] laide dintgrales, et non par une dfinition astucieuse de la variable alatoire Y S .
2.19 Le ratio de perte (loss ratio) R est dfini comme tant le montant total des
sinistres pays pendant lanne, S, divis par le montant total des primes
reues pendant lanne, . Une compagnie dassurance souhaite bien entendu conserver ce ratio sous un certain niveau pour ne pas tre en difficult
financire. Pour ce faire, elle offre un bonus B ses agents la fin de lanne
si le ratio de perte pour lanne est infrieur 75 %. Le montant du bonus est
calcul comme suit :
0,75 R
.
B = max 0,
3
Calculer le montant espr du bonus si = 600 000 et que la distribution de
la variable alatoire S est une Pareto avec paramtres = 3 et = 700 000.
2.20 Soit X , une variable alatoire reprsentant le montant dun sinistre. Un assureur souhaite connatre les paiements sa charge pour un contrat dassurance
incluant une franchise dcroissante (disappearing deductible). Dans ce type
de contrat, lassur assume en entier tout sinistre infrieur d et lassureur
assume en entier tout sinistre suprieur d . Entre d et d , le paiement
effectu par lassureur est une fonction linaire du montant dun sinistre.
a) Dfinir la variable alatoire Y P reprsentant le montant pay par paiement
pour un contrat avec une franchise dcroissante.
b) Trouver lexpression gnrale en termes de E [X ], E [X ; x] et F X (x) du montant pay par paiement espr.
95
96
Rponses
2.1 0,4946
2.2 0,4686
2.3 a) 0,0175 b) 0,1813
2.4 0,0034
2.5 0,567
2.6 0,226
2.7 a) (3 e 2d 2e d )/4 b) 0,0541 c) 0,0576
2.8 a) 1 091,09 b) 0,8438
2.11 a) 0,1069 b) 0,4107 c) 1 255,23
2.12 a) 33 962 b) +8,04 % c) 2,87 %
2.14 0,22
2.15 3 857
2.16 8,0925
2.17 a) 2 (E [X 2 ; u 2 ] E [X 2 ; d 2 ] 2d E [X ; u] + 2d E [X ; d ]) 2 (E [X ; u] E [X ; d ])2
b) (1 + r )(E [X ; u/(1 + r )] E [X ; d /(1 + r )])
2.19 76 559,55
2.20 a) (E [X ] + d /(d d )E [X ; d ] d /(d d )E [X ; d ])/(1 F X (d ))
x
x
x
x
x
x
x
x
Dans ce chapitre, on dveloppe directement partir des donnes des estimateurs de la fonction de rpartition, de la fonction de densit et de certaines quantits lies dune variable alatoire sans faire dhypothse quant la distribution
complte de cette variable alatoire. Cest lapproche non paramtrique.
Cette approche a comme avantages dtre flexible et de bien prendre en compte
la disparit des donnes. De plus, elle peut tre trs prcise lorsque le nombre de
donnes est grand. Par contre, elle est souvent moins efficace quune approche
paramtrique et linfrence statistique qui en rsulte est plus complique.
97
98
3.1
Statistiques descriptives
3.1.1
Moments
Lorsquil existe, on nomme moment dordre k de la variable alatoire X lesprance de cette dernire leve la puissance k :
0k = E [X k ]
Z
x k d F (x).
=
0
3.1.2
Moments centraux
x variance :
2 = Var[X ]
= 2
= 02 2 ;
x coefficient de variation :
p
Var[X ]
CV(X ) =
E [X ]
= ;
99
E [(X )3 ]
Var[X ]3/2
3
= 3;
2 (X ) =
3.1.3
Lesprance rsiduelle est utile pour mesurer lpaisseur relative des queues des
distributions. Nous avons tabli son lien avec lesprance limite la section 2.7,
soit :
e(x) = E [X x|X > x]
Z
1
=
(y x)d F (y)
1 F (x) x
E [X ] E [X ; x]
=
.
1 F (x)
Lesprance limite est la moyenne des valeurs censures x dune variable
alatoire.
Plus intressante pour dcrire des donnes, lesprance rsiduelle reprsente
quant elle lexcdent moyen des valeurs dune variable alatoire au-dessus de
x. Ainsi, plus cette esprance est grande (pour un mme x), plus la queue de la
distribution est lourde. Voir la figure 3.1, inspire de Hogg et Klugman (1984), pour
une comparaison de la fonction desprance rsiduelle pour quelques lois usuelles.
Nous verrons plus loin comment estimer les esprances limite et rsiduelle
partir des donnes. On pourra ensuite utiliser ces estimations pour identifier
des modles potentiels pour les donnes. Par exemple, si lesprance rsiduelle
empirique crot plus ou moins linairement, on favorisera la loi de Pareto comme
modle.
3.1.4
Quantiles
h0
Pa
r
et
o
1500
100
1000
1),
ll ( <
Weibu ormale
logn
e(x)
Gamma
(
> 1)
500
Exponentielle
1)
Weibull (
>
500
1000
1500
Ainsi, tout nombre pour lequel la fonction de rpartition vaut au moins p est le
100p e quantile de la distribution. Pour une distribution continue, cela se rsume
au point p tel que F (p ) = p. Par contre, le quantile nest pas unique pour les fonctions de rpartition non strictement croissantes, comme les fonctions en escalier
des variables alatoires discrtes.
Par exemple, pour la fonction de rpartition reprsente la figure 3.2, on a
3 0,6 < 4.
3.2
1.0
0.6
0.4
0.2
F(x)
0.8
0.0
3.2.1
Donnes individuelles
#X j x
n
n
1 X
I {X x}
n j =1 j
n
1 X
I {X =x} .
n j =1 j
#X j = x
n
101
102
Histogram of dental
0.004
1.0
ecdf(dental)
0.003
0.6
Density
0.4
0.001
Fn(x)
0.2
0.0
0.000
0.002
0.8
500
1000
1500
500
1000
1500
2000
dental
Tel que mentionn ci-dessus, ces fonctions sont des estimateurs de la fonction de
rpartition et de la fonction de densit de probabilit, dans lordre, cest--dire que
nous considrerons que
F (x) = F n (x)
f(x) = f n (x).
On peut dmontrer quelques proprits de lestimateur F n de la fonction de
rpartition. Soit
(
1, X j x
B j = I {X j x} =
0, sinon
et
Y =
n
X
Bj
j =1
= nF n (x).
On a Pr[B j = 1] = F (x) et B j Bernoulli(F (x)), do
Y = nF n (x) Binomiale(n, F (x)).
E [F n (x)] =
et
Var[Y ]
n2
nF (x)(1 F (x))
=
n2
F (x)(1 F (x))
=
.
n
Var[F n (x)] =
Par consquent,
0,
x <0
5/10, 0 x < 1
F 10 (x) =
9/10, 1 x < 3
1,
x 3.
5/10, x = 0
4/10, x = 1
f 10 (x) =
1/10, x = 3
0,
ailleurs.
103
104
1.0
ecdf(x)
fn(x)
0.6
0.4
0.0
0.2
Fn(x)
0.8
0.0
0.5
1.0
1.5
2.0
2.5
3.0
3.2.2
Donnes groupes
Classe
Frquence
(c 0 , c 1 ]
(c 1 , c 2 ]
..
.
(c j 1 , c j ]
..
.
(c r 1 , c r ]
n1
n2
..
.
nj
..
.
nr
On a donc
F n (c 0 ) = 0
F n (c j ) =
j
1X
ni ,
n i =1
j = 1, . . . , r.
cj x
x c j 1
F n (c j 1 ) +
F n (c j )
Fn (x) =
c j c j 1
c j c j 1
F n (c j ) F n (c j 1 )
= F n (c j 1 ) +
(x c j 1 )
c j c j 1
pour c j 1 x c j .
Par consquent, logive Fn (x) est dfinie ainsi :
Fn (x) =
0,
(c j x)F n (c j 1 ) + (x c j 1 )F n (c j )
1,
c j c j 1
x c0
, c j 1 < x c j
x > cr .
105
106
1.0
ogive(gdental)
0.6
0.8
F(x)
0.4
0.2
0.0
1000
2000
3000
4000
f(x) =
0,
F n (c j ) F n (c j 1 )
0,
c j c j 1
x c0
=
nj
n(c j c j 1 )
, c j 1 < x c j
x > cr .
107
1.0
ecdf(dental)
0.8
0.6
Fn(x)
0.4
0.2
0.0
500
1000
1500
Soit
Y =
n
X
I {X i c j 1 }
i =1
= nF n (c j 1 )
Binomiale(n, F X (c j 1 ))
et
Z=
n
X
I {c j 1 <X i c j }
i =1
= n(F n (c j 1 ) F n (c j ))
Binomiale(n, F X (c j 1 ) F n (c j )).
108
Y (c j c j 1 ) + Z (x c j 1 )
n(c j c j 1 )
do
E [Fn (x)] =
F (c j 1 )(c j c j 1 ) + (F (c j ) F (c j 1 ))(x c j 1 )
(c j c j 1 )
6= F (x).
Ainsi, Fn (x) est un estimateur biais de F (x).
Pour lhistogramme, on a
fn (x) =
Z
,
n(c j c j 1 )
do
E [ fn (x)] =
F (c j ) F (c j 1 )
c j c j 1
6= f (x)
et f(x) est un estimateur biais de f (x).
De plus, on a
x c j 1 2 Var[Z ]
Var[Y ]
+
n2
c j c j 1
n2
x c j 1 Cov(Y , Z )
,
+2
c j c j 1
n2
Var[Fn (x)] =
o
Var[Y ] = nF (c j 1 )(1 F (c j 1 )),
Var[Z ] = n(F (c j ) F (c j 1 ))(1 F (c j ) + F (c j 1 ))
et
Cov(Y , Z ) = nF (c j 1 )(F (c j ) F (c j 1 )).
Pour lhistogramme, on a
Var[ f(x)] =
(F (c j ) F (c j 1 ))(1 F (c j ) + F (c j 1 ))
n(c j c j 1 )2
0,
3/7,
F 7 (x) = 5/7,
6/7,
1,
3/7,
2/7,
f 7 (x) = 1/7,
1/7,
0,
x < 11
11 x < 13
13 x < 17
17 x < 19
x 19
x = 11
x = 11
x = 17
x = 19
ailleurs.
Frquence
(10, 12]
(12, 14]
(14, 18]
(18, 23]
3
2
1
1
0,
14 (x 10),
1 (x 9),
F7 (x) = 71
(x + 6),
28
35 (x + 12),
1,
x 10
10 < x 12
12 < x 14
14 < x 18
18 < x 23
x > 23
109
110
Histogram of xg
1.0
1.0
ecdf(x)
0.8
0.8
0.6
0.2
0.0
0.4
Density
0.6
0.4
0.2
0.0
Fn(x)
10
12
14
16
18
20
10
12
14
16
18
20
22
xg
0, x 10
14 , 10 < x 12
1 , 12 < x 14
f7 (x) = 71
28 , 14 < x 18
35 , 18 < x 23
0, x > 23.
Complter lexemple en excutant le code informatique de la section 3.6 correspondant ce segment de matire.
3.3
3.3.1
n
X
I {xi =y j } ,
i =1
P
le nombre de fois que la donne y j apparat dans lchantillon. On a donc kj=1 s j =
n. La variable s j reprsente le nombre de dcs la dure y j . Enfin, on dfinit r j
comme le nombre de sujets dans le groupe tmoin juste avant les s j dcs la
111
112
yj
sj
rj
1
2
3
4
5
6
0,8
2,9
3,1
4,0
4,1
4,8
1
2
1
2
1
1
8
7
5
4
2
1
n
X
I {xi y j }
i =1
k
X
si .
i=j
n
X
I {xi y j } +
i =1
k
X
n
X
I {ui y j }
i =1
si +
i=j
n
X
I {ui y j } .
i =1
0,
x < y1
rj
F n (x) = 1 , y j 1 x < y j , j = 2, . . . , k
1,
x > yk .
3.3.2
Estimateur de Nelson-alen
h(x) =
= ln S(x).
La fonction de taux dincidence cumul est donc troitement lie la fonction de
survie S(x). Lestimateur de Nelson-alen exploite ce lien pour estimer la fonction
de survie.
Lestimateur de Nelsonalen du taux dincidence cumul est
0,
x < y1
jX
1
si
, y j 1 x < y j , j = 2, . . . , k
H (x) = i =1 r i
k s
X
, x yk .
i =1 r i
Par la suite, un estimateur de la fonction de survie est
= e H (x) .
S(x)
113
114
F (x) = 1 e H (x) .
H (x) =
0,
j 1
x < y1
f n (y i )
, y j 1 x < y j , j = 2, . . . , k
1
F n (y i 1 )
i =1
k
X
f n (y i )
, x yk ,
i =1 1 F n (y i 1 )
0,
8 = 0,125,
0,125 + 7 = 0,411,
H (x) = 0,411 + 15 = 0,611,
0,411 + 24 = 1,111,
1,111 + 12 = 1,611,
1,611 + 11 = 2,611,
0 x < 0,8
0,8 x < 2,9
2,9 x < 3,1
3,1 x < 4,0
4,0 x < 4,1
4,1 x < 4,8
x 4,8
F IG . 3.8: Fonction R qui retourne une fonction pour calculer lestimateur de Nelson
alen H (x) en tout x
0,
0,8825,
0,6632,
= 0,5430,
0,3293,
0,1997,
0,0735,
0 x < 0,8
0,8 x < 2,9
2,9 x < 3,1
3,1 x < 4,0
4,0 x < 4,1
4,1 x < 4,8
x 4,8.
115
0.6
0.2
0.4
^
S(x)
1.5
1.0
0.0
0.0
0.5
^
H(x)
2.0
0.8
2.5
1.0
116
Klugman et collab. (2012a, p. 233) justifient succinctement que, sous des hypothses raisonnables, H (x) est un estimateur sans biais de H (x) et que
d H (y j )]
Var[
3.3.3
j
X
si
2
i =1 r i
Estimateur de Kaplan-Meier
Lestimateur de KaplanMeier aussi communment appel estimateur produitlimite est principalement utilis comme estimateur de la fonction de survie en
prsence de donnes tronques et censures.
Lide, assez simple au demeurant, derrire lestimateur de KaplanMeier est la
suivante. Sil y a r 1 sujets dans le groupe tmoin juste avant la dure y 1 et que s 1
117
r 1 s1 r 2 s2
,
r1
r2
et ainsi de suite jusqu y k . On obtient alors un estimateur de la fonction de survie.
la notation introduite la sous-section 3.3.1, on ajoute d j , le point de troncature de la donne j , ou lge dentre dans le groupe tmoin du sujet j . Ainsi, pour
des donnes tronques et censures, le groupe risque y j est
r j = # sujets sortis du groupe ou aprs y j
# sujets non encore entrs dans le groupe y j
n
n
n
X
X
X
=
I {xi y j } +
I {ui y j }
I {di y j }
i =1
i =1
i =1
i =1
i =1
S n (t ) =
1,
Y j 1 r i s i
i =1
i =1
0 t < y1
,
y j 1 t < y j , j = 2, . . . , k
ri
r i si
ou 0, t y k .
ri
k r s
Y
i
i
r
i
i =1
118
Y
k r s t /w
i
i
,
r
i
i =1
o w = max(x 1 , . . . , x n , u 1 , . . . , u n ), pour t y k .
Remarque. Si les donnes de lchantillon ne sont ni tronqus, ni censures, alors
r 1 = n et r j = r j 1 s j 1 . Lestimateur de KaplanMeier est alors identique la
fonction de survie empirique 1 F n (t ).
Klugman et collab. (2012a, p. 229231) dmontrent que
E [S n (y j )] =
S(y j )
S(y 1 )
j
S(y j )2 Y
S(y j 1 ) S(y j )
Var[S n (y j )] =
1+
1 .
S(y 1 )2 i =1
r j S(y j )
Or, cette formule est complexe et difficile utiliser. On peut lui prfrer lapproximation
S(y j )2 X
S(y i 1 ) S(y i )
Var S n (y j )
.
S(y 1 )2 i =1
r i S(y i )
Ceci dit, cette formule nest pas dun plus grand secours en pratique puisque
la fonction de survie est inconnue. On utilisera donc plutt lapproximation de
Greenwood
j
X
si
d n (y j )] = S n (y j )2
Var[S
.
i =1 r i (r i s i )
Exemple 3.5. On rutilise les donnes de lexemple 3.3, mais en ajoutant comme
information que les donnes j = 1, 3, 8 sont censures. Autrement dit, on a les
donnes
{0,8+ , 2,9, 2,9+ , 3,1+ , 4,0, 4,0, 4,1, 4,8+ },
o un signe + en indice indique une donne censure. Aucune donne nest tronque, donc d j = 0 pour tout j . Lchantillon contient un total de huit donnes, soit
quatre observations x 2 = 2,9, x 5 = x 6 = 4,0, x 7 = 4,1 et quatre donnes censures
u 1 = 0,8, u 3 = 2,9, u 4 = 3,1, u 8 = 4,8. Le tableau 3.3 prsente les calculs des groupes
risque, soit les valeurs des quantits y j , s j et r j . On a
119
yj
sj
rj
1
2
3
2,9
4,0
4,1
1
2
1
4+40 = 7
3+10 = 4
1+10 = 2
TAB. 3.3: Calculs des groupes risque pour les donnes de lexemple 3.5
1,
71 = 6 ,
7
S 8 (t ) = 67 42
3
7 4 = 7 ,
3 21
3
= 14
,
7
2
0 t < 2,9
2,9 t < 4,0
4,0 t < 4,1
t 4,1.
Exemple 3.7. On a les donnes de survie suivantes, prsentes sous forme dintervalles (d j , x j ] ou (d j , u j ], les donnes censures tant identifies par un signe + en
indice :
(0, 0,1]
(0, 0,5]
(0, 0,8]
(0, 0,8+ ]
(0, 1,8]
(0, 1,8]
(0, 2,1]
(0, 2,5]
(0, 2,8]
(0, 2,9+ ]
(0, 2,9+ ]
(0, 3,9]
(0, 4,0+ ]
(0, 4,0]
(0, 4,1]
(0, 4,8+ ]
(0, 4,8]
(0, 4,8]
(0, 5,1]
(0, 5,1]
(0, 5,1]
(0, 5,1]
(0, 5,1]
(0, 5,1]
(0, 5,1]
(0, 5,1]
(0, 5,1]
(0, 5,1]
(0, 5,1]
(0, 5,1]
(0,3, 5,1]
(0,7, 5,1]
(1,0, 4,1+ ]
(1,8, 3,1+ ]
(2,1, 3,9]
(2,9, 5,1]
(2,9, 4,8]
(3,2, 4,0+ ]
(3,4, 5,1]
(3,9, 5,1]
0.0
0.2
0.4
S8(t)
0.6
0.8
1.0
120
donc :
1,
29
30 = 0,967,
0,967( 29
30 ) = 0,934,
29
0,934( 30
) = 0,903,
27
0,903( 29 ) = 0,841,
0,841( 27 ) = 0,811,
28
27
S 40 (t ) = 0,811( 28
) = 0,782,
26
0,782( 27 ) = 0,753,
25
0,753( 27
) = 0,697,
25
0,697( 26 ) = 0,670,
22
0,670( 23
) = 0,641,
20
0,641( 21 ) = 0,550,
0,550( 0 ) = 0,
17
0 t < 0,1
0,1 t < 0,5
0,5 t < 0,8
0,8 t < 1,8
1,8 t < 2,1
2,1 t < 2,5
2,5 t < 2,8
2,8 t < 3,9
3,9 t < 4,0
4,0 t < 4,1
4,1 t < 4,8
4,8 t < 5,1
t 5,1.
121
yj
sj
rj
1
2
3
4
5
6
7
8
9
10
11
12
0,1
0,5
0,8
1,8
2,1
2,5
2,8
3,9
4,0
4,1
4,8
5,1
1
1
1
2
1
1
1
2
1
1
3
17
30
31 1 = 30
32 2 = 30
33 4 = 29
34 6 = 28
35 7 = 28
35 8 = 27
39 12 = 27
40 14 = 26
40 17 = 23
40 19 = 21
40 23 = 17
TAB. 3.4: Calculs des groupes risque pour les donnes de lexemple 3.7
S(3) S(5)
.
S(3)
S 40 (3) S 40 (5)
S 40 (3)
0,753 0,550
=
0,753
= 0,270.
=
Var[S 40 (5)]
.
(0,753)2
On utilise lapproximation de Greenwood pour Var[S 40 (5)]. Seulement, on considre maintenant la survie jusquau temps 3 comme accomplie. De plus, 5 nest pas
0.0
0.2
0.4
S40(t)
0.6
0.8
1.0
122
0,550 2
2
1
1
3
=
+
+
+
0,753
27(25) 26(25) 23(22) 22(19)
d 2 q3 |S 40 (3) = 0,753] =
Var[
= 0,00728.
(Cet exemple est lexercice 12.16 de Klugman et collab. (2012a) avec une lgre
modification des donnes pour clarifier le contexte.)
3.3.4
123
j =0
i =0
j =0
di
P j 1
i =0
(u i + x i ),
j = 1, . . . , k.
Lorsque les points de troncature et de censure se trouvent lintrieur des intervalles plutt quaux bornes, on utilise parfois plutt
rj =
(
(d 0 u 0 )/2,
(d j u j )/2 +
j =0
P j 1
i =0
(d i u i x i ),
j = 1, . . . , k.
j =0
1,
j ) = Y j 1 r j x j
S(c
, j = 1, . . . , k.
i =0
rj
Exemple 3.8. On a observ 19 sinistres, dont six avaient une franchise de 250, six
une franchise de 500 et sept une franchise de 1 000. De plus, trois sinistres, de
montants 1 000, 2 750 et 5 500, ont atteint la limite de la police dassurance. Pour les
124
(c j , c j +1 )
dj
xj
uj
rj
0
1
2
3
4
5
(250, 500)
(500, 1 000)
(1 000, 2 750)
(2 750, 5 500)
(5 500, 6 000)
(6 000, 10 000)
6
6
7
0
0
0
1
2
4
7
1
1
0
1
1
1
0
0
6
6 + 6 1 = 11
11 + 7 3 = 15
15 + 0 5 = 10
10 + 0 8 = 2
2+01 = 1
TAB. 3.5: Calculs des groupes-risque pour les donnes de lexemple 3.8
16 autres sinistres, un se trouve dans lintervalle (250, 500), deux dans (500, 1 000),
quatre dans (1 000, 2 750), sept dans (2 750, 5 500), un dans (5 500, 6 000) et un dans
(6 000, 10 000). On veut estimer la probabilit de payer plus de 5 500 pour un contrat
avec une franchise de 500.
Puisquun paiement de 5 500 correspond un sinistre de 6 000, on cherche
Pr[X > 6 000|X > 500] =
S(6 000)
.
S(500)
S(250)
=1
5
S(500)
= = 0,833
6
9
= 0,682
S(1 000) = (0,833)
11
11
750) = (0,682)
S(2
= 0,500
15
3
500) = (0,500)
S(5
= 0,150
10
1
000) = (0,150)
S(6
= 0,075
2
0
S(10
000) = (0,075)
= 0.
1
Ainsi, une estimation de la probabilit recherche est
000) 0,075
S(6
=
= 0,090.
0,833
S(500)
125
3.4
Les techniques destimation de la fonction de rpartition tudies jusqu maintenant rsultent en des fonctions discrtes (exception faite de logive, bien sr).
Cela peut tre vu comme un problme lorsque la variable alatoire sous-jacente
est continue et que le nombre de donnes est faible. La technique destimation
par noyaux permet destimer la fonction de rpartition (ou la fonction de densit)
par une fonction lisse et continue. De plus, elle est flexible, dune criture mathmatique relativement simple et conviviale et elle permet une meilleure prise en
compte de la contribution des observations situes dans le voisinage dun point x.
Un estimateur liss par noyaux est obtenu en remplaant chaque donne de
lchantillon par une variable alatoire continue et en attribuant un poids de 1/n
chacune de ces variables alatoires. Celles-ci doivent tre identiques un changement dchelle ou de position prs (fonction de la donne associe). La distribution
continue utilise pour chaque variable alatoire est appele le noyau. Ultimement,
cest comme si une donne tire de cet estimateur tait dabord tire de la densit
empirique, puis ensuite dune densit centre la valeur obtenue.
Formellement, lide est dattribuer chacune des donnes une distribution
dont la moyenne sera la donne (ce nest pas obligatoire, mais cela permet la distribution lisse et la distribution empirique davoir la mme moyenne) plutt que
son poids 1/n. Comme dans les sections prcdentes, y 1 < < y k , reprsentent
les k (k n) donnes uniques et tries de lchantillon x 1 , . . . , x n . Soit k y (x) une
fonction de densit de probabilit de moyenne y (le noyau) et K y (x) la fonction de
rpartition correspondante. Les estimateurs par noyaux sont de la forme
f(x) =
s
X
f n (y j )k y j (x)
j =1
et
F (x) =
s
X
f n (y j )K y j (x).
j =1
Ainsi, f(x) est un mlange discret des noyaux k y 1 (x), . . . , k y k (x). La figure 3.12 fournit
un exemple.
Les noyaux les plus courants sont les suivants.
0.10
0.00
0.05
Density
0.15
126
10
12
F IG . 3.12: Six noyaux gaussiens (traits briss) et leur somme (trait plein)
0,
x < y b
1
k y (x) =
, y b x y +b
2b
0,
x > y + b,
0,
x < y b
x y +b
K y (x) =
, y b x y +b
2b
1,
x > y + b.
127
0,
x < y b
y
+
b
, y b x y
2
k y (x) = y bx + b
, y x y +b
b2
0,
x > y + b,
0,
x < y b
(x y + b) , y b x y
K y (x) =
2b 2
(yx+b)2
, y x y +b
2b 2
1,
x > y + b.
3. Noyau gamma Le noyau la donne y est une distribution gamma de paramtre
de forme et de moyenne y :
k y (x) =
(/y) x 1 e x/y
,
()
x > 0.
x < y b
0,
x y 2
3
k y (x) =
1
, y b x y +b
4b
b
0,
x > y + b.
5. Noyau normal, ou gaussien Le noyau la donne y est une distribution normale
de moyenne y et de variance b 2 :
x y
k y (x) =
, < x <
x b
y
, < x < ,
K y (x) =
b
o () et () sont, dans lordre, la fonction de densit de probabilit et la
fonction de rpartition de la distribution N (0, 1).
Dans tous les cas sauf le noyau gamma, la valeur b reprsente la largeur de
bande de lestimateur par noyaux. Silverman (1986) suggre dutiliser, du moins
dans un premier temps,
R/1,34)
0,9 min(,
,
b=
1/5
n
128
3.5
Voir le code informatique de la section 3.6 pour des exemples dutilisation et la comparaison de divers noyaux avec lhistogramme pour
un jeu de donnes.
3.5.1
Moments
n
X
i =1
x kj f n (x j )
n
1X
xk .
n i =1 j
En particulier,
= x
2 = 02 2
n
1 X
=
x 2 x 2
n j =1 j
Pn
2
j =1 (x j x)
=
= s2.
n
129
n j (c 2j c 2j 1 )
2n(c j c j 1 )
n j (c j + c j 1 )
2n
do
=
n n c +c
X
j
j 1
j
j =1
n
1X
k
(x j )
=
n i =1
et
0k =
3.5.2
r Z
X
cj
j =1 c j 1
k fn (x) d x
(x j )
Esprance limite
1 X
=
x j + u(1 F n (u))
n x j <u
=
n
1 X
min(x j , u).
n j =1
130
E [X ; u] =
=
.=
0
r Z ci
X
i =1 c i 1
jX
1 Z c i
i =1 c i 1
Z cj
+u
min(x, u) fn (x) d x
x fn (x) d x +
fn (x) d x + u
u
c j 1
x fn (x) d x
r Z
X
ci
i = j +1 c i 1
fn (x) d x
2
2
n i c i + c i 1 n j u c j 1
+
2
n 2(c j c j 1 )
i =1 n
r
X
n j u(c j u)
ni
+
+u
n c j c j 1
i = j +1 n
jX
1
n i c i + c i 1 n j c j 1 + u u c j 1
=
+
2
n
2
c j c j 1
i =1 n
r
X ni
cj u
nj
+u
.
+u
n c j c j 1
i = j +1 n
jX
1
Exemple 3.10. Les classes [0, 100), [100, 200) et [200, 500) contiennent respectivement 5, 10, et 2 donnes. Pour estimer E [X ; 120], on spare la classe [100, 200) en
deux classes : [100, 120) et [120, 200). Les donnes sont rparties proportionnellement entre les deux classes. On obtient alors
10 120 100
5
(50) +
(100)
17
17 200 100
10 200 120
2
+
(120) + (120)
17 200 100
17
1 670
=
.
17
E [X ; 120] =
3.5.3
Esprance rsiduelle
=
e(x)
3.5.4
E [X ; x]
.
1 F n (x)
Quantiles
131
132
3.5.5
Soit X 1 , . . . , X n un chantillon alatoire dune variable X , X (1) < < X (n) les
statistiques dordre cet chantillon et
N=
n
X
I {X j p }
j =1
133
1
1
Pr[N = i ] Pr i W < i +
,
2
2
1
1
.
Pr[Yi p < Y j ] Pr i W < j 1 +
2
2
3.6
Code informatique
###
### FONCTIONS DE RPARTITION ET DE PROBABILIT EMPIRIQUES
###
## DONNES INDIVIDUELLES
## On va travailler avec les donnes "dental" fournies dans la
## premire dition de Loss Models. Il s'agit de 10 montants de
## rclamation en assurance dentaire avec une franchise de 50. Ces
## donnes sont distribues avec actuar.
data(dental, package = "actuar") # charger les donnes
dental
# les donnes
##
##
##
##
##
134
#
#
#
#
#
135
136
137
##
## On dfinit un objet contenant les donnes D2 de Loss Models.
d <- c(rep(0, 30), 0.3, 0.7, 1.0, 1.8, 2.1, 2.9, 2.9, 3.2, 3.4, 3.9)
x <- c(0.1, 0.5, 0.8, 0.8, 1.8, 1.8, 2.1, 2.5, 2.8, 2.9, 2.9, 3.9,
4.0, 4.0, 4.1, 4.8, 4.8, 4.8, rep(5.0, 12), 5.0, 5.0, 4.1,
3.1, 3.9, 5.0, 4.8, 4.0, 5.0, 5.0)
e <- rep(0, 40); e[c(4, 10, 11, 13, 16, 33, 34, 38)] <- 1
Surv(d, x, e)
# objet contenant les donnes
## On peut ensuite passer ces donnes
## Par dfaut, la fonction ajuste les
## Kaplan-Meier.
(fit <- survfit(Surv(d, x, e) ~ 1)) #
(sfit <- summary(fit))
#
plot(fit)
#
138
col
col
col
col
=
=
=
=
"red")
"red")
"red")
"red")
#
#
#
#
uniforme
triangulaire
Epanechnikov
gaussien
139
avec le noyau
= "e", adjust
= "e", adjust
= "e", adjust
= "e", adjust
d'Epanechnikov.
= 0.25)) # 25 % du dfaut
= 0.50)) # 50 % du dfaut
= 1.25)) # 125 % du dfaut
= 1.50)) # 150 % du dfaut
140
## individuelles ou groupes.
(mu <- emm(dental))
emm(dental, 2) - emm(dental)^2
emm(dental - mu, 3)/sd(dental)^3
emm(gdental, 2) - emm(gdental)^2
#
#
#
#
moyenne
variance (biaise)
asymtrie
variance donnes groupes
## ESPRANCE LIMITE
## actuar fournit la fonction 'elev' pour calculer la fonction
## d'esprance limite empirique en n'importe quelle limite (donnes
## individuelles et groupes). Le fonctionnement est similaire
## 'ecdf' et 'ogive'.
(lev <- elev(dental))
# dfinition
lev(knots(lev))
# calcul en chacun des points
lev(1200)
# calcul en un point quelconque
mean(pmin(dental, 1200))
# vrification
plot(lev)
# graphique
## Pour illustrer avec des donnes groupes, on vrifie le rsultat de
## l'exemple 3.10.
x <- grouped.data(cj = c(0, 100, 200, 500), nj = c(5, 10, 2))
elev(x)(120)
# = 1670/17
## QUANTILES
## La fonction 'quantile' de R calcule les quantiles empiriques pour
## des donnes individuelles selon l'un de sept (!) algorithmes
## diffrents. Sans autre argument, la fonction retourne les quartiles
## de l'chantillon selon l'algorithme 7 prsent dans les remarques
## de la section 3.5.4. Notre algorithme est le numro 6.
quantile(dental)
# utilisation simple
quantile(dental, c(0.45, 0.80)) # autres quantiles
y <- sort(dental)
# statistiques d'ordre
0.3 * y[3] + 0.7 * y[4]
# algorithme de R
quantile(dental, 0.30)
# vrification
0.7 * y[3] + 0.3 * y[4]
# notre algorithme
quantile(dental, 0.3, type = 6) # numro 6
## Pour des donnes groupes, la mthode de 'quantile' dfinie par
## actuar retourne simplement l'inverse de l'ogive.
quantile(gdental)
# utilisation simple
quantile(gdental, c(0.45, 0.80)) # autres quantiles
ogive(gdental)(quantile(gdental)) # inverse de l'ogive
3.7. Exercices
3.7
141
Exercices
3.1 Un assureur prsente les cots (en millions de $) crs par les crasements de
mtorites :
{3, 5, 5, 6, 8, 8, 8, 8, 9, 10, 10, 11, 11, 11, 16, 21, 23, 26, 29, 36}.
a) Faire des graphiques de la fonction de rpartition empirique et de la fonction de masse de probabilit empirique du cot des crasements.
b) partir des bornes c 0 = 2, c 1 = 7, c 2 = 12, c 3 = 22 et c 4 = 38, crire lquation
de logive.
c) En utilisant les mmes bornes quen b), crire lquation de lhistogramme.
3.2 Le tableau ci-dessous prsente les sinistres enregistrs par un assureur.
Classe
Nombre de sinistres
(0, 50]
(50, 150]
(150, 250]
(250, 500]
(500, 1 000]
(1 000, )
36
x
y
84
80
0
Total
Nombre de sinistres
(0, 500]
(500, 1 000]
(1 000, 2 000]
(2 000, 5 000]
(5 000, 10 000]
(10 000, 25 000]
(25 000, )
200
110
x
y
Soit Fn () logive correspondant ces donnes. Sachant que F500 (1 500) = 0,689
et F500 (3 500) = 0,839, calculer la valeur de y.
142
Nombre de sinistres
0 1 000
1 000 3 000
3 000 5 000
5 000 10 000
10 000 25 000
25 000 50 000
50 000 100 000
100 000 et plus
16
22
25
18
10
5
3
1
3.7. Exercices
143
Nombre de sinistres
Groupe-risque
3
10
11
2
52
40
19
6
Nombre de sinistres
Groupe-risque
1
1
1
1
1
1
20
19
18
17
16
15
144
en dfinissant k j sur lintervalle [1, 1]. Par analogie, pour le taux dincidence, on va utiliser
x yj
s
1X
h n (y j )k j
h(x)
=
,
b j =1
b
Frquence
100
200
300
400
500
1
4
10
4
1
3.7. Exercices
145
Nombre de donnes
(0, 50]
(50, 100]
(100, 200]
(200, 500]
20
34
22
24
1 y
y
e
,
()
a) Dterminer la distribution de X .
y > 0.
146
X
X 1
(j)
distribution de .
e) Calculer les 45e et 70e quantiles lisss des donnes de la partie c).
3.7. Exercices
147
Rponses
3.1 a)
F20 (x) =
0,
(x 2)/25,
(x 5)/10,
x 2
2<x 7
7 < x 12
(x + 58)/100, 12 < x 22
(x + 42)/80, 22 < x 38
1,
x > 38
b)
f20 (x) =
0,
1/25,
1/10,
x 2
2<x 7
7 < x 12
1/100, 12 < x 22
1/80, 22 < x 38
0,
x > 38.
3.2 120
3.3 81
3.4 0,396
3.5 0,17
3.6 a) 225
3.7 a) 0 b) 0,05 c) 0,125 d) 0,1333 e) 0 f) 0,02 g) 0,095
3.8 1,0264
3.9 a) 9,225 et 10,2369
3.10 (0,0964, 0,7036)
3.11 0,5880
3.12 a) 0,05, 0,1026, 0,1582, 0,2170, 0,2795, 0,3462 b) 0,00001449
3.13 0,559
3.14 1 = 0, 2 = 3,125
3.15 a) 38,6 b) 42,5
148
4 Modles paramtriques
potentiels
Objectifs du chapitre
x Dfinir et connatre les principales caractristiques (densit, fonction de rpartition,
x
x
x
moments, esprance limite) des lois gamma et bta la base des familles de
distributions du mme nom.
Dterminer la distribution dune fonction dune ou plusieurs variables alatoires et
en calculer les principales caractristiques.
Dmontrer leffet dun changement de paramtre sur une distribution de svrit.
Appliquer les techniques suivantes pour crer de nouvelles distributions : multiplication par une constante, lvation une puissance, exponentielle et logarithme,
mlanges continus et discrets, raccordement.
Caractriser les modles de montants de sinistres en fonction de leurs valeurs
extrmes.
Dans ce chapitre, nous allons dvelopper des modles qui pourront tre utiliss
pour reprsenter les montants des sinistres. Les modles paramtriques potentiels
devront permettre de reprsenter une variable alatoire qui est toujours positive,
qui peut prendre des valeurs sur lintervalle (0, ) et qui est continue. De plus, il
faudra trouver des faons de caractriser les valeurs extrmes.
4.1
Famille gamma
Les distributions de la famille gamma sont parmi les plus utilises en actuariat.
On rappelle ici quelques rsultats importants leur propos.
La distribution gamma est continue et dfinie sur lintervalle (0, ), ce qui
en fait un choix intressant pour la modlisation des montants de sinistres en
149
150
f X (x) =
x > 0,
o
() =
t 1 e t d t
F X (x) =
1 y
y
e
d y.
()
o
1
(; x) =
()
t 1 e t d t ,
> 0, x > 0
1
()
x 1 e x/ ,
x > 0.
151
0.2
0.4
fX(x)
1.0
0.0
0.0
0.5
fX(x)
1.5
0.6
2.0
(a) 1
(b) > 1
t
= t (e )|0 +
t 1 e t d t
0
= 0 + ()
= ().
1
X
j =0
x j e x
.
j!
152
F X (x) = 1
1
X (x) j e x
j =0
j!
25 125
5
= 1e
1+5+
+
2
6
F X (2,5) = 1 e 5
= 0,7349741.
On peut vrifier ce rsultat avec R :
> pgamma(2.5, 4, 2)
[1] 0.7349741
E [X ] =
x +k1 e x d x
() 0
Z
( + k) +k
=
x +k1 e x d x
() +k
0 ( + k)
( + k)
=
()k
k
( + 1) ( + k 1)
,
k 1
k =0
E [X k ] = 1,
|k|
, k 1, |k| < .
( 1)( 2) ( |k|)
153
= ( + 1; x) + x(1 (; x)).
4.2
Famille bta
La distribution bta est continue et dfinie sur un intervalle fini. Elle est donc
rarement utilise pour la modlisation des montants de sinistres. Par contre, elle
joue un rle important dans la cration de nouvelles distributions, comme nous
le verrons plus loin. Lorsque son domaine est restreint lintervalle (0, 1), elle est
souvent utilise comme modle pour une probabilit.
On dbute la prsentation avec une version de la distribution bta trois
paramtres, o le troisime paramtre est la borne suprieure du support de la
distribution. Ainsi, soit X Bta(, , ) une variable alatoire avec fonction de
densit de probabilit
( + ) 1 x 1
x 1
1
()()
1
1 x
1
x 1
, 0 < x < ,
=
1
(, )
f (x) =
o
(a, b) =
=
Z
0
t a1 (1 t )b1 d t
(a)(b)
(a + b)
1
(a, b)
Z
0
t a1 (1 t )b1 d t ,
1.5
0.0
0.95
0.5
1.0
fX(x)
1.15
1.05
fX(x)
1.25
154
0.0
0.2
0.4
0.6
0.8
1.0
0.0
0.2
0.4
0.6
0.8
1.0
(a) = < 1
x k x 1
x 1
1
dx
Z
0
( + )
()()
u +k1 (1 u)1 d u
k ( + )( + k)
=
()( + + k)
pour k > . Si k Z+ (entier positif ), ce rsultat se rduit
E [X k ] =
k ( + 1) ( + k 1)
.
( + )( + + 1) ( + + k 1)
=
( + 1, ; x/) + x(1 (, ; x/)).
+
E [X ; x] =
4.3
155
156
= 1
+ x/c
c
= 1
,
c + x
do Y Pareto(, c). Le paramtre constitue donc un paramtre dchelle de
la distribution de Pareto.
4.4
4.4.1
1
f X (y/c).
c
4.4.2
X 1/
,
> 0,
= 1 e (y) ,
et
f Y (y) = y 1 e (y) ,
do Y Weibull(, ). Attention : la Weibull est parfois paramtre de telle sorte
que
f (x) = y 1 e y ,
ce qui correspond plutt la transformation Y = (X /)1/ .
Si < 0, on a
F Y (y) = 1 F X ((y) )
= e (y)
= e (y)
157
158
Si > 0, alors
F Y (y) = (; (y) )
et
f Y (y) =
y 1 e (y) ,
()
F Y (y) = 1 (; (y) )
et
f Y (y) =
e (y)
,
y +1 ()
f(x)
0.2
0.4
0.6
0.8
0.6
0.4
0.0
0.2
0.0
f(x)
159
0.8
1.0
(a) = 1
(b) > 1
Gamma transforme
Gamma
,
, , H
HH = 1
= 1
HH
H
j
H
Weibull
,
H
H
HH
=1
=1
H
H
j
H
Exponentielle
160
X 1/
Y =
.
1 X
Par la technique du changement de variable, on obtient
f Y (y) =
( + )
(y/)
,
()() y(1 + (y/) )+
y > 0,
do
( + )
F Y (y) =
()()
Z
0
(x/)
d x.
y(1 + (x/) )+
Or, avec le changement de variable u = v/(1+ v), v = (x/) (de sorte que 0 < u < 1),
on obtient
Z
( + ) u 1
F Y (y) =
u (1 u)1 d u
()() 0
= (, ; u).
La distribution de la variable alatoire Y est une bta transforme de paramtres
, , et . La famille bta transforme compte plusieurs membres dont, entre
autres :
1. Burr(, , ) lorsque = 1 ;
2. Pareto gnralise(, , ) lorsque = 1 ;
3. Pareto(, ) lorsque = = 1.
Ces relations sont illustres la figure 4.5.
De lexemple 1.17, si X 1 Gamma(, 1), X 2 Gamma(, 1) et X 1 et X 2 sont
indpendantes, alors Z = X 1 /(X 1 + X 2 ) Bta(, ). Puisque
Z
X1
=
,
1 Z
X2
on peut aussi dfinir la bta transforme comme la distribution de
1/
X1
Y =
.
X2
Pareto gnralise
H
=
HH
1
H
?
HH
H
Log-logistique
, ,
Burr inverse
HH
H= 1
HH
, ,
Burr
, ,
H Paralogistique inverse
HH =
,
=1
H1
HH
6
=
H
?
j
H
, , ,
,
= 6
1
=
Paralogistique
Bta transforme
=1
H
j
H
?
,
Pareto
161
HH
j
Pareto inverse
,
4.4.3
1
f X (ln y),
y
y > 0.
(ln y)1
y +1 ()
y >1
F (x) = (; ln x).
On dit que Y Log-gamma(, ). Pour modliser des montants de sinistre, on
utilise normalement Y 1.
162
4.4.4
Mlanges de distributions
f X (x) =
f (x|)u() d
1 e e x d
() 0
Z
=
+11 e (x+) d
() 0
Z
( + 1)
(x + )+1 +11 (x+)
=
e
d
+1
() (x + )
( + 1)
0
=
,
(x + )+1
f X (x) =
do X Pareto(, ).
Remarques.
1. Lorsque est une variable alatoire continue, on parle de mlange continu.
2. Un autre mlange (continu) bien connu est le suivant : si
X | = Poisson()
Gamma(, )
alors X Binomiale ngative(, /( + 1)).
3. On a les relations
E [X ] = E [E [X |]]
Var[X ] = E [Var[X |]] + Var[E [X |]].
1
X
f (x|)p (1 p)1
=0
= p f 1 (x) + (1 p) f 2 (x).
Ce mlange est appel un mlange discret.
De manire gnrale, la fonction de densit dun mlange discret de n densits
est
n
X
f X (x) =
p i f i (x),
i =1
o p 1 , . . . , p n [0, 1] et
Pn
i =1 p i
= 1. On remarque que
F X (x) =
n
X
p i F i (x)
i =1
et
E [X k ] =
n
X
i =1
p i E [X ik ].
Pn
163
164
0.010
0.000
f(x)
0.020
50
100
150
200
250
0.010
0.000
f(x)
50
100
150
200
250
200
250
0.006
0.000
f(x)
0.012
Mlange
50
100
150
x
165
0.10
0.00
f(x)
0.20
Poisson(2)
10
15
20
0.15
0.00
f(x)
0.30
10
15
20
15
20
0.10
0.00
f(x)
0.20
Mlange
10
x
166
4.4.5
Raccordement de distributions
ci
c i 1
f i (x) d x = 1
..
0,
ailleurs.
x >0
16
g 2 (x) =
0,4
x 15 e 0,4x ,
(16)
x >0
et
g 1 (x)
,
0 < x < 20
0,6
G 1 (20)
f X (x) =
g 2 (x)
, 30 x 50.
0,4
G 2 (50) G 2 (30)
4.5
4.6
Valeurs extrmes
167
168
0.02
0.00
g1(x)
Gamma(2, 0.1)
10
20
30
40
50
40
50
40
50
0.02
0.00
g2(x)
0.04
Gamma(16, 0.4)
10
20
30
x
0.02
0.00
f(x)
Raccordement
10
20
30
x
169
Famille
Loi
Racine (alias)
Bta transforme
Bta transforme
Burr
Log-logistique
Paralogistique
Pareto gnralise
Pareto
Burr inverse
Pareto inverse
Paralogistique inverse
trbeta (pearson6)
Gamma transforme
Gamma transforme inverse
Gamma*
Gamma inverse
Weibull*
Weibull inverse
Exponentielle*
Exponentielle inverse
trgamma
Log-normale*
Log-gamma
Pareto translate
Bta gnralise
lnorm
Gamma transforme
Autre
burr
llogis
paralogis
genpareto
pareto (pareto2)
invburr
invpareto
invparalogis
invtrgamma
gamma
invgamma
weibull
invweibull (lgompertz)
exp
invexp
lgamma
pareto1
genbeta
Cette section prsente diverses manires de classer les lois de probabilit selon
le poids de leur queue.
4.6.1
Moments
On a
k
E [X ] =
x k f (x) d x.
Si E [X k ] < pour toutes les valeurs de k, alors cest que f (x) dcrot linfini plus
rapidement que x k ne crot. Par contre, si E [X k ] = pour k k 0 , alors f (x) dcrot
170
4.6.2
Lide est que si S(x) est grand pour x grand, on est en prsence dune queue
paisse. Puisque, par dfinition,
lim S X (x) = 0,
il sagit dun concept relatif. Il faut donc comparer les taux de dcroissance des
fonctions de survie de deux distributions. On value
lim
S 1 (x)
x S 2 (x)
= lim
f 1 (x)
x f 2 (x)
Si le rsultat est une constante c, les deux distributions ont des queues comparables.
Si le rsultat est 0, la queue de f 1 est plus lgre. Si le rsultat est , la queue de f 2
est plus lgre.
Exemple 4.13. Soit f 1 la fonction de densit de probabilit dune loi Gamma(, )
et f 2 celle dune loi Pareto(, ). En liminant les termes qui ne dpendent pas de x,
on a
lim
x 1 e x
f 1 (x)
= lim
f 2 (x) x (x + )1
= lim e (+1) ln(x+)+(1) ln xx .
x
f 1 (x)
x f 2 (x)
= 0.
4.6.3
171
Taux dincidence
Il est galement possible dutiliser le taux dincidence pour comparer les queues
des distributions. On a
f (x)
h(x) =
S(x)
F 0 (x)
=
S(x)
F (x + t ) F (x)
= lim
t 0
t S(x)
Pr[x < X x + t |X > x]
= lim
.
t 0
t
Ainsi, on peut interprter h(x)t comme la probabilit davoir un incident (sinistre,
dcs) x sachant que lon a survcu jusqu x. Donc, si le taux dincidence est
petit, la probabilit davoir un incident est faible ou, inversement, la probabilit de
survivre est grande. On a donc une queue lourde.
De manire gnrale, si h(x) est dcroissant (DFR), on est en prsence dune
distribution avec queue lourde. Si h(x) est croissant (IFR), on est en prsence dune
distribution avec queue lgre.
De plus, puisque
h(x) =
d
ln S(x)
dx
et
S(x) = e
Rx
0
h(y) d y
x 0,
172
1 x
x
e
,
()
x 0
y 1 y
= 1+
e
.
x
Cette dernire fonction est croissante en x si < 1 et dcroissante en x si > 1.
Ainsi, la distribution gamma est IFR (queue lgre) si > 1 et DFR (queue lourde)
si < 1. Si = 1 (distribution exponentielle), le taux dincidence est constant.
4.6.4
Esprance rsiduelle
E [X ] E [X ; x]
.
1 F X (x)
Or, puisque
x
E [X ; x] =
=
Z0 x
(1 F (y)) d y
S(y) d y,
4.7. Exercices
173
Si la fonction desprance rsiduelle, e(x), est croissante, on dit que f est IMRL
et quelle possde une queue paisse. Dans le cas contraire, f est dite DMRL et elle
possde une queue lgre.
En rsum, f est la fonction de densit de probabilit dune loi avec queue
paisse si
e(x) =
=
e y d y
e x
1
,
une constante.
Exemple 4.16. Soit X Pareto(, ), cest--dire que
S(x) =
+x
On a
R
e(x) =
( + y) d y
( + x)
( + x)+1
=
( 1)( + x)
+x
=
,
1
4.7
Exercices
4.1 Soit X , une variable alatoire avec densit Pareto(, ) reprsentant le montant
dun sinistre et c > 0, une constante. Dmontrer que la distribution de Y = c X
est une distribution Pareto(, c).
174
1 |x/|
e
,
2
< x < .
X
X +
4.7. Exercices
175
4.9 Un assureur modlise des donnes laide de la variable alatoire X qui a une
distribution de Pareto de paramtres et . On pose
Y = ln(1 + X /).
Dterminer la distribution de Y .
4.10 Un assureur automobile a dans sa base de donnes les montants des sinistres
de 2004. Il estime que les sinistres obissaient alors une loi Burr( = 0,5, =
2, = 3). Pour sen servir le premier janvier 2007, il se doit de les mettre jour
selon les considrations suivantes :
x 2005 : inflation de 4 % ;
x 2006 : inflation de 4,5 % ; et
x nouvelles taxes de 16 %.
Quelle est la probabilit davoir un sinistre suprieur 4 en 2007 ?
4.11 Soit X , la variable alatoire reprsentant le montant dun sinistre (en millions)
pour lanne 2006. Sa fonction de densit de probabilit est
f (x) = 3x 4 ,
x 1.
x = 1, 2, . . .
( + )( + 1)( + x 1)
.
()()( + + x)
176
f (x| = ) = x 1 e x , x > 0.
Aussi, on suppose que Gamma(, ). Dmontrer que la distribution marginale de X est une Burr(, , 1/ ).
4.16 On suppose que le montant dun sinistre pour un groupe dassurs a une
distribution Burr(5, 1, ). Si est une ralisation de la variable alatoire
pour ce groupe dassurs et que lon suppose que Gamma(10, 2), trouver
lesprance et la variance du montant dun sinistre pour un assur pris au
hasard.
4.17 Soit le taux dchec suivant pour le montant dun sinistre pour une valeur
donne de ,
3
,
(x|) =
x +
o x est la ralisation de la variable alatoire X reprsentant le montant dun
sinistre et est la ralisation de la variable alatoire o Gamma(10, 0,01).
Trouver lesprance et la variance du montant dun sinistre pris au hasard.
4.18 Comparer les queues des lois Gamma(, ) et Log-normale(, 2 ).
4.19 Soit X , une variable alatoire reprsentant le montant dun sinistre et lesprance de vie rsiduelle suivante
e(x) = 2 000 + 2x.
Pour un contrat dassurance comportant une limite suprieure de 10 000,
trouver le ratio dlimination de perte (LER) de lassureur.
4.20 Le tableau ci-dessous prsente lesprance de vie rsiduelle pour certaines
valeurs de x.
x
e(x)
0
4
9
14
4
7
10,75
14,5
4.7. Exercices
177
4.21 On construit une distribution raccorde sur les sous-intervalles (0, 2), (2, 8)
et (8, 16) avec les poids respectifs 0,5, 0,20 et 0,30. Dans chacun des sousintervalles, on utilise une distribution gamma, de moyenne gale au point
milieu du sous-intervalle et de variance gale 1. crire la densit de probabilit obtenue sur (0, 16). La rponse sera en fonction de la gamma incomplte.
4.22 On construit un modle raccord avec une distribution uniforme sur lintervalle (0, 10) et une loi de Pareto de paramtres = 3 et = 100 sur le reste des
valeurs positives. Quels poids doivent tre accords aux distributions pour
que la densit obtenue soit continue ?
4.23 a) Comparer les queues dune distribution Weibull(, ) et dune distribution
Weibull inverse(, ) en utilisant les critres suivants :
i) lexistence des moments ; et
ii) la comparaison des fonctions de survie.
b) En utilisant une distribution Weibull et une distribution Weibull inverse
dont les moyennes et variances sont gales, comparer graphiquement les
queues des distributions.
4.24 Soit Y , une variable alatoire telle que
f Y (y) =
S X (y)
E [X ]
pour une variable alatoire X quelconque. On dit quune telle distribution est
quilibre. Dmontrer que
M Y (t ) =
M X (t ) 1
t E [X ]
converge.
4.25 Un assureur modlise ses sinistres par une variable alatoire X avec densit
f (x) = (1 + 2x 2 )e 2x ,
a) Calculer la fonction de survie S X (x).
b) Calculer le taux dincidence h(x).
x 0.
178
Rponses
4.2 F Y (y) = 12 e ln(y)/ I {0<y<1} + (1 12 e ln(y)/ )I {y1}
4.4 Bta(, )
4.5 Burr inverse(, 4, 5)
4.6 a) Log-gamma(, ) b) E [Y ] = (/( 1)) , Var[Y ] = (/( 2)) (( 1))2
c) Non
4.7 a) Pareto(, (1 + i )) b) Burr(, , (1 + i )) c) f Y (y) = (1 + i ) (ln(y) ln(1 +
i ))1 y 1 /()
4.8 Burr(, , 1/ )
4.9 Exponentielle()
4.10 0,6870
4.11 a) F (x) = 1 1,331x 3 , x 1,1 b) 0,125
4.12 4/11
4.13 X Pareto gnralise(, , )
4.16 5/4 et 145/48
4.17 500, 850 000
4.18 La queues de la distribution log-normale est plus lourde que celle de la distribution gamma.
4.19 0,30
4.20 a) Pareto(7/3, 16/3) b) 3,0215
4.7. Exercices
179
4.21
f X (x) =
0,5e x
(1; 2)
0<x 2
0,2
525 x 251 e 5x
,
2<x 8
(25)
12 x
e
0,3
, 8 < x 16
(144; 192) (144; 96)
(144)
4.22 3/14
4.25 a) (1 + x + x 2 )e 2x
b) 2 (1 + 2x)/(1 + x + x 2 )
c) (1 + x + 0,5x 2 )/(1 + x + x 2 )
5 Modlisation paramtrique
Objectifs du chapitre
x Estimer les paramtres dune distribution du montant des sinistres partir de
donnes individuelles ou groupes, modifies ou non, en utilisant lune ou lautre
des mthodes suivantes : la mthode des moments, la mthode des quantiles, la
mthode du maximum de vraisemblance, la mthode de la distance minimale.
x Estimer la variance de lestimateur du maximum de vraisemblance et sen servir
pour dterminer et calculer un intervalle de confiance pour des paramtres dune
distribution.
x Pour les cas requrant une optimisation numrique, calculer les divers estimateurs
ci-dessus avec R.
Lapproche paramtrique consiste choisir un modle connu pour le phnomne sous tude. Ce modle comportera habituellement des paramtres quil
faudra dterminer dune manire ou une autre. En gnral, on optera pour une
technique ayant un certain degr dobjectivit et se basant sur des observations du
phnomne.
En termes plus statistiques, on cherche estimer le vecteur de paramtres
= (1 , . . . , p )0 dune fonction de densit de probabilit ou fonction de masse de
probabilit f (x; ) partir dun chantillon alatoire X 1 , . . . , X n de la population.
Un estimateur ponctuel de est une fonction = T (X 1 , . . . , X n ) de lchantillon.
Toute statistique est un estimateur ponctuel.
De plus, on se rappelle quun estimateur, T (X 1 , . . . , X n ), est une fonction de
lchantillon (une variable alatoire), alors quune estimation, T (x 1 , . . . , x n ), est une
ralisation dun estimateur (un nombre) obtenue en considrant un chantillon en
particulier.
On classe les techniques destimation tudies dans ce chapitre en deux grandes
classes :
181
182
Modlisation paramtrique
1. techniques bases sur la reproduction de certaines caractristiques des donnes
par les mmes caractristiques du modle (moments, quantiles, etc.) ;
2. techniques bases sur loptimisation dun critre objectif (vraisemblance, distance).
Une autre technique destimation qui repose sur des principes forts diffrents
est lestimation bayesienne. Celle-ci nest toutefois pas aborde dans cet ouvrage
puisquil est davantage pertinent den discuter lors de ltude de la thorie de la
crdibilit. Pour les intresss, Klugman et collab. (2012a, chapitre 15) traitent
brivement du sujet.
5.1
n
1X
X k,
n i =1 i
k = 1, . . . , p,
n
1X
(X i X )2 .
=
n i =1
183
1
= 01 = X ,
do = X 1 .
b) Gamma(, ). On pose
= X
Var[X ] = 2
1
2,
=
=
E [X ] =
do
=
=
X 2
2
X
2
c) Pareto(, ). On a
= X
1
2
Var[X ] =
( 1)2 ( 2)
2
2.
=
=
1
2
E [X ] =
184
Modlisation paramtrique
En rsolvant pour et , on obtient
=
=
5.2
2
2
2 X 2
2 + X 2 )
X (
2 X 2
i = 1, . . . , p.
5.3
Lestimation laide des mthodes des moments et des quantiles est souvent
facile raliser, mais les estimateurs ainsi obtenus nont pas souvent de bonnes
proprits, la principale raison tant que ces techniques nutilisent quune petite
partie de linformation contenue dans lchantillon pour raliser lestimation. La
mthode du maximum de vraisemblance permet dutiliser une plus grande partie
de linformation contenue dans lchantillon. Ceci est particulirement important
lorsquune population a une queue paisse.
5.3.1
Donnes individuelles
i =1
n
Y
i =1
(F (x i ; ) F (x i ; )).
Si F est continue, Pr[X i = x i ] = 0 et f (x i ; ) nest pas une probabilit. Nanmoins, on peut trouver un intervalle (a i , b i ] tel que a i < x i b i , do
L() =
n
Y
i =1
(F (b i ; ) F (a i ; )).
n
Y
f (x i ; )i
i =1
n
Y
i =1
f (x i ; ).
185
186
Modlisation paramtrique
Lestimateur du maximum de vraisemblance (EMV) de est la valeur qui
maximise la fonction de vraisemblance L() sur son domaine de dfinition.
De plus, maximise L() si, et seulement si, maximise la fonction de logvraisemblance l () = ln L().
On dfinit les fonctions de score
S j () =
ln L(),
j
j = 1, . . . , p.
j = 1, . . . , p.
En gnral, ces quations ne sont pas linaires. La rsolution du systme dquations ncessite donc des mthodes numriques de type NewtonRaphson. Nous
verrons comment cela peut se faire avec R la section 5.4.
Exemple 5.3. Soit X 1 , . . . , X n un chantillon alatoire tir dune distribution exponentielle de paramtre . La fonction de vraisemblance du paramtre est
donc
L() =
n
Y
e xi
i =1
P
n ni=1 x i
= e
do
l () = n ln
n
X
xi
i =1
et donc
S() =
n
n X
xi .
i =1
i =1 X i
1
.
X
187
Exemple 5.4. Soit X 1 , . . . , X n un chantillon alatoire tir dune distribution Lognormale(, 2 ), cest--dire que
n 1 ln x 2 o
1
f (x) = p
exp
,
2
22 x
x > 0.
!1
n 1X
n
n ln x 2 o
Y
i
2
2 n/2
L(, ) = (2 )
xi
exp
2
i =1
i =1
ou, de manire quivalente,
l (,) =
n
n
X
n
1 X
ln x i 2 (ln x i )2 .
ln(22 )
2
2 1
i =1
Ainsi,
n
1 X
(ln x i )
2 1
n
n
1 X
S 2 (,2 ) = 2 + 4 (ln x i ).
2
2 1
S 1 (,2 ) =
De S 1 (, 2 ) = S 2 (, 2 ) = 0, on trouve
n
1X
ln x i
n 1
n
1X
2 =
2.
(ln x i )
n 1
f (x) =
n (x )2 o
exp
,
2x 3
22 x
x > 0.
L(,) =
n/2 n
Y
1
n
n (x )2 o
X
i
x i3/2 exp 2
2
xi
i =1
188
Modlisation paramtrique
ou encore
l (,) = ln L(, )
n
n (x )2
3X
X
n
ln x i 2
.
= ln
2
2
2 i =1
2 i =1
xi
On a donc
S 1 (,) =
=
l (,)
n (x )2
X
i
3 i =1
xi
n
X
1
x
n
i
2
1
et
l (,)
n (x )2
n
1 X
i
=
2
.
2 2 i =1
xi
S 2 (,) =
Or, puisque
n (x )2
n
X
X
i
= n x 2 n + 2
x i1 ,
x
i
i =1
i =1
189
et
l () = 4 ln() 130,
do
l 0 () =
4
130.
Y =
0,
X <d
X d
u d
d X <u
X u
= max(min(X , u) d , 0).
On a donc
f Y (x) =
F X (d ),
x =0
f X (x + d ), 0 x < u
1 F (u), x = u d .
X
n
Y
f Y (x i )
i =1
= (F X (d ))
s
Y
i =1
f X (x r +i + d ) (1 F X (u))t .
190
Modlisation paramtrique
5.3.2
Donnes groupes
r
Y
[F (c j ) F (c j 1 )]n j ,
j =1
r
X
n j ln[F (c j ) F (c j 1 )].
j =1
5.3.3
g
() = g ().
p
2.
Par exemple, lEMV de est
Soit X 1 , . . . , X n un chantillon alatoire dune distribution avec fonction de
densit de probabilit f (x; ) (on suppose p = 1 pour le moment afin de simplifier
la prsentation). Soit aussi
n = T (X 1 , . . . , X n ),
un estimateur sans biais quelconque de . Sous des conditions gnralement
satisfaites mais souvent difficiles tablir 1 , on a que
Var[n ]
1
,
I ()
I () = nE
ln f (X ; )
= nE
ln f (X ; )
2
= information de Fisher.
1. Voir le thorme 13.5 de Klugman et collab. (2012a) pour les conditions exactes.
1
,
I ()
l (; X 1 , . . . , X n )
l (; X 1 , . . . , X n ) .
= E
2
I () = E
do lEMV est asymptotiquement sans biais et efficace. De plus, on a que, asymptotiquement toujours,
n N (, I 1 ()).
Lorsque p > 1, les rsultats demeurent les mmes, sauf que la distribution est
une normale multivarie :
n N (, I 1 ()),
o I () est la matrice dinformation p p dont llment en position (r, s), r, s =
1, . . . , p, est
I ()r s = nE
5.4
2
2
ln f (X ; ) = E
ln l (; X 1 , . . . , X n ) .
r s
r s
191
192
Modlisation paramtrique
Exemple 5.8. Soit X 1 , . . . , X n un chantillon alatoire dune loi gamma avec fonction de densit de probabilit
f (x; , ) =
1 x
x
e
,
()
n
X
i =1
n
X
ln f (x i ; , )
( ln ln () + ( 1) ln x i x i )
i =1
= n ln n ln () + ( 1)
n
X
ln x i
i =1
n
X
xi .
i =1
l (, ) = n ln n
+
ln x i = 0
() i =1
n
n X
l (, ) =
x i = 0.
S 2 (, ) =
i =1
S 1 (, ) =
d
0 (x)
ln (x) =
dx
(x)
est appele fonction digamma. Sa valeur numrique est obtenue dans R avec
digamma(). Dans le mme ordre dide, la fonction
(x) =
1 (x) =
d
d2
(x) =
ln (x)
dx
d x2
5.4.1
La mthode de NewtonRaphson multivarie est similaire la mthode univarie, sauf quelle fait appel des notions simples de calcul matriciel.
On considre le systme n quations non linaires et n inconnues
f 1 (x 1 , . . . , x n ) = 0
f 2 (x 1 , . . . , x n ) = 0
..
.
f n (x 1 , . . . , x n ) = 0.
Soit x = (x 1 , . . . , x n ). On peut alors crire, de manire plus compacte,
f 1 (x) = 0
f 2 (x) = 0
..
.
f n (x) = 0
ou encore f (x) = 0, avec
f 1 (x)
.
f (x) = ..
f n (x)
= ( f 1 (x), . . . , f n (x))0 .
De plus, soit
f 1 (x) f 1 (x)
x 1
x 2
f 2 (x) f 2 (x)
x 2
J (x) = x 1
.
..
..
.
f n (x) f n (x)
x 1
x 2
f i (x)
=
.
x j nn
...
...
...
f 1 (x)
x n
f 2 (x)
x n
..
.
f n (x)
x n
193
194
Modlisation paramtrique
o x 0 = (x 10 , . . . , x n0 ) est un vecteur de valeurs de dpart.
On peut utiliser comme critre darrt de la procdure itrative
kx k x k1 k = kJ (x k1 )1 f (x k1 )k < ,
o
kxk = max |x i |.
1i n
J () =
S i ()
j
pp
par la matrice dinformation
I () = E [J ()]
2
= E
l () .
i j
5.4.2
Virtuellement tous les logiciels mathmatiques ou statistiques modernes proposent des outils doptimisation de fonctions linaires ou non et de sommes
de carrs. Ceci rend la programmation de lalgorithme de NewtonRaphson inutile
dans la grande majorit des cas. Le Solveur dExcel est un tel outil. Le chapitre 7 de
Goulet (2012) prsente les fonctions doptimisation disponibles dans R. Plusieurs
des exemples de ce chapitre portent dailleurs sur lvaluation numrique destimateurs du maximum de vraisemblance. La principale fonction doptimisation de R
est optim.
Exemple 5.9. On simule un chantillon de taille 100 dune loi gamma de paramtres = 5 et = 2, puis on estime ensuite les paramtres de la distribution par
le maximum de vraisemblance laide de la fonction optim de R. La fonction
minimiser est
n
X
l (, ) = ln f (x i ; , ).
i =1
NA
$convergence
[1] 0
$message
NULL
195
196
Modlisation paramtrique
de fournir la fonction de densit (ou de masse) de probabilit en argument
fitdistr().
Exemple 5.10. On rpte lestimation des paramtres de la loi gamma pour lchantillon alatoire simul lexemple 5.9. La distribution gamma tant reconnue par
fitdistr(), il suffit den spcifier le nom :
> fitdistr(x, "gamma")
shape
5.1044237
rate
2.0414616
(0.6995568) (0.2940095)
Les valeurs entre parenthses sont les carts types estims des estimateurs.
i =1
n
X
n
( + 1) (x i + )1 = 0.
i =1
(5.1)
(5.2)
Ce systme dquations na pas de solution explicite. On aura donc recours loptimisation numrique avec la fonction fitdistr.
Les estimateurs des moments de et seront utiliss comme valeurs de dpart
dans fitdistr(). De lannexe A de Klugman et collab. (2008, 2012a), on sait que
2(m 2 m 12 )
m 2 2m 12
m2 m1
=
,
m 2 2m 12
o m k = n 1
Pn
k
i =1 X i . Pour lchantillon simul ci-dessus, on obtient
scale = lambda.mme))
shape
7.170763
scale
6256.725402
21.931933) (21840.466242)
Lorsque lon tente destimer des paramtres strictement positifs comme ceux
dune Pareto, il nest pas rare dtre confront des soucis de convergence lorsque
la fonction doptimisation sgare dans les valeurs ngatives. On peut pallier ce
problme avec le truc suivant : plutt que destimer les paramtres eux-mmes, on
estime leurs logarithmes. Ceux-ci demeurent alors valides sur tout laxe des rels.
Cest lAstuce Ripley explique la section 11.6 de Goulet (2013).
Bien que pas ncessaire pour lchantillon simul ci-dessus, on illustre cette
approche en dfinissant dabord une nouvelle fonction de densit paramtre en
fonction de ln et ln :
> f <- function(x, lshape, lscale)
+
197
198
Modlisation paramtrique
lshape
1.8165661
lscale
8.5646370
(0.7320815) (0.8431280)
Les valeurs obtenues ci-dessus sont les logarithmes des paramtres. Les estimations
de et elles-mmes sont :
> exp(unname(fit$estimate))
[1]
6.150701 5242.936212
(La fonction unname enlve les noms dun objet qui, sils taient encore prsents,
porteraient ici confusion.)
Il peut arriver que lestimation des paramtres ne se passe pas aussi bien que cidessus. Dans un tel cas, on peut essayer dutiliser linformation dont nous disposons
sur les drives partielles (5.1) et (5.2) pour trouver les estimations du maximum de
vraisemblance. Dune part, on peut passer ces drives (le gradient de la fonction
minimiser) directement optim() par le biais de largument gr ; consulter la
rubrique daide pour de plus amples informations. Dautre part, si lon isole dans
(5.2), on obtient ici
P
1
ni=1 (x i + )
=
P
1
n ni=1 (x i + )
et lquation (5.1) devient alors
n
X
n
= 0,
+ n ln
ln(x i + )
i =1
n <- length(x)
x <- x + b
u <- b * sum(1/x)
n * (n - u)/u + n * log(b) - sum(log(x))
+
+ }
x 1
,
(1 + (x/) )+1
x >0
et
F (x) = 1
1
,
(1 + (x/) )
x > 0.
Pour simplifier la notation, on suppose que les donnes x 1 , . . . , x nm sont celles qui
ne sont pas censures.
On cherche donc estimer les paramtres de la loi de la variable alatoire
X Burr(, , ) partir dobservations de la variable alatoire Y = min(X , u) dont
la fonction de densit de probabilit est
(
g (x) =
f (x),
y <u
1 F (u), x u.
L(, , ) =
nm
Y
i =1
x i
(1 + (x i /) )+1
1
(1 + (u/) )
199
200
Modlisation paramtrique
et la fonction de log-vraisemblance est
l (, , ) = (n m)(ln ln )
nm
X
ln x i
+ ( 1)
( + 1)
i =1
nm
X
ln(1 + (x i /) )
i =1
m ln(1 + (u/) ).
Ainsi,
n m
m ln(1 + (x i /) )
m
X
ln(1 + (u/) ),
S 1 (, , ) =
i =1
nm
X
n m
ln x i
S 2 (, , ) =
(n m) ln +
i =1
nm
X (x i /) ln(x i /)
( + 1)
1 + (x i /)
i =1
(u/) ln(u/)
1 + (u/)
et
S 3 (, , ) =
+
X (x i /)
(n m) ( + 1) nm
+
i =1 1 + (x i /)
m (u/)
.
1 + (u/)
Pour trouver les estimations des paramtres , et , il faut rsoudre le systme trois quations et trois inconnues avec la mthode de NewtonRaphson
multivarie. On peut galement minimiser numriquement la fonction de logvraisemblance ngative avec les fonctions optim ou fitdistr. Dans ce cas, la
fonction coverage permet dobtenir facilement la densit g (x) ci-dessus :
> g <- coverage(dburr, pburr, limit = u)
460,26
45,56
48,60
588,26
240,00
350,30
253,15
725
93,39
52,36
60,41
118,10
204,77
27,90
389,96
41,68
452,59
2,39
725
111,30
725,00
493,15
710,67
197,68
725
260,80
14,80
160,67
42,20
40,30
340,34
239,10
8,75
281,80
100,45
369,10
327,20
86,55
138,50
70,50
166,08
241,55
225,37
370,60
564,15
134,54
119,42
55,43
191,94
81,01
201
485,60
218,60
123,80
371,34
29,40
312,08
68,25
369,39
266,05
725
scale
3125.70221
15.65693) (4526.98265)
202
Modlisation paramtrique
Classe
Nombre de donnes
(0, 1]
(1, 2]
(2, 4]
(4, 8]
(8, 12]
(12, 16]
(16, 20]
(20, ]
4
3
3
4
2
0
1
3
r
X
n j ln(F (c j ; ) F (c j 1 ; )).
j =1
nj = c(4, 3, 3, 4, 2, 0, 1, 3))
203
$par
[1] 0.58211645 0.06674133
$value
[1] 39.26542
$counts
function gradient
79
NA
$convergence
[1] 0
$message
NULL
5.5
5.5.1
N (, Var[]),
204
Modlisation paramtrique
avec
=
Var[]
1
.
I ()
Prz /2 < q
< z /2 = 1 ,
Var[]
o z est le 100(1 )e centile dune N (0, 1). On peut rcrire cette expression sous
la forme
q
q
q
q
z /2 Var[], + z /2 Var[]
Var[]
I ()
En principe, un intervalle de confiance pour est alors
t /2,n1
d ].
Var[
d ].
Var[
205
I()
avec
2
n
X
ln f (x i ; )
=
i =1
2
n
X
ln
f
(x
;
)
=
i
2
=
i =1
=
I()
l (; x i , . . . , x n )
= 2 l (; x i , . . . , x n )
=
I()
f (x) = x 1 e (x/)
ln f (x) = ln ln + ( 1) ln x
et
x
ln f (x) = +
2
( + 1) x
ln f (x) = 2
.
2
Par consquent,
l ()
n
X
=
ln f (x i )
i =1
Pn
i =1 x i
= n
.
S() =
206
Modlisation paramtrique
En posant S() = 0, on obtient
=
1/
i =1 x i
Pn
et
=
d ]
Var[
2
.
n2
Pn
z /2 p .
n
=
= 541,36
10
p
q
= 2 930 683 = 85,5962,
d ]
Var[
20
do (373,59, 709,13). On obtient les mmes rsultats avec fitdistr :
> fitdistr(dental, "weibull", shape = 2)
scale
541.31341
( 85.57799)
En pratique, il nest pas rare que lon souhaite estimer non pas pour , mais pour
une fonction h() de . Par la proprit dinvariance de lestimateur du maximum
de vraisemblance, on sait que
= h(),
h()
mais quen est-il dun intervalle de confiance pour cette estimation ? En gnral, il
peut tre trs complique. Pour
sagit dun calcul difficile car la distribution de h()
207
simplifier le travail, il existe une mthode, la mthode delta, qui est valide pour n
grand.
En ngligeant les termes de degr deux et plus dans le dveloppement de Taylor
autour de , on a
de h()
h() + h 0 ()( ).
h()
Ainsi, puisque est asymptotiquement sans biais,
= h()
E [h()]
et
= [h 0 ()]2 Var[].
Var[h()]
On a donc, lorsque n ,
h()
Exemple 5.16. Trouver un intervalle de confiance approximatif de niveau 1
pour la probabilit quune observation dune distribution exponentielle dpasse
700 partir des donnes dental.
On cherche un estimateur par intervalle pour
h() = Pr[X > 700] = e 700 .
On sait que
1
=
X
=
Var[]
2
.
n
=
Var[h()]
208
Modlisation paramtrique
et
=
d )]
Var[h(
p
Un intervalle de confiance 95 % pour Pr[X > 700] est donc 0,124 1,96 0,006707
ou 0,124 0,161.
5.5.2
Les ides restent essentiellement les mmes lorsque plusieurs paramtres sont
estims simultanment par la mthode du maximum de vraisemblance. On sait
que, toujours asymptotiquement,
Normale multivarie(, I 1 ()).
Si lon considre une fonction h(1 , . . . , p ) h(), la mthode delta dit que
Normale multivarie(h(), h T I 1 ()h),
h()
o
h = h()
h()
1
..
=
.
.
h()
p
Un intervalle de confiance de niveau 1 pour h() est donc
z /2
h()
p
h T I 1 ()h.
209
/2
n 1
2 ) = e +
h(,
/2
De plus,
2
1
ln f (X ; , 2 ) = 2
2
2
1
(ln X )2
2
ln
f
(X
;
,
)
=
(2 )2
24
6
2
ln X
ln f (X ; , 2 ) =
,
2
4
do, en ralisant que ln X N (, 2 ),
n
2
I (, 2 ) =
0
0
n
24
n
I 1 (, 2 ) =
et
Enfin,
"
h=
2
h(, )
h(, 2 )
2
=e
0
.
24
n
+2 /2
1
1
2
210
Modlisation paramtrique
On a donc que, lorsque n ,
2 )] = e +
E [h(,
/2
2 )] = h T I 1 (, 2 )h
Var[h(,
2
4
2+2
=e
,
+
n
2n
do un estimateur par intervalle avec un niveau de confiance approximatif de 95 %
de la moyenne du log-normale est
s
e +
/2
1,96e +
/2
2
4
+
.
n
2n
Exemple 5.18. Soit X 1 , . . . , X n un chantillon alatoire tir dune distribution Pareto inverse de paramtres et . On cherche un intervalle de confiance pour
lestimateur du maximum de vraisemblance du mode de cette distribution. Tout
dabord, on a
x 1
f (x; t , ) =
(x + )+1
et le mode est
m=
( 1)
.
2
n
X
ln f (x i ; , )
i =1
= n ln + n ln + ( 1)
n
X
i =1
( + 1)
n
X
i =1
ln(x i + )
ln x i
211
et
n
n
X
n X
ln(x i + )
ln x i
+
i =1
i =1
n
X
n
1
S 2 (, ) = ( + 1)
i =1 x i +
S 1 (, ) =
n
S 1 (, ) = 2
n
X
1
n
S 2 (, ) = 2 + ( + 1)
2
(x
+
i )
i =1
n
X
2
1
l (, ) =
S 1 (, ) =
.
x
i =1 i +
( + k)x
=
x k
dx
+k 0
(x + )+k+1
E [X k ],
=
+k
o X Pareto inverse( + k, ). Ainsi, en supposant que > 2, on obtient
1 ()(2)
1
E [(X + ) ] =
+1
( + 1)
1
=
( + 1)
et
2 ()(3)
E [(X + ) ] =
+2
( + 2)
2
= 2
.
( + 1)( + 2)
2
212
Modlisation paramtrique
Par consquent,
2
n
E
l (, ; X 1 , . . . , X n ) = 2
2
2
n
l (, ; X 1 , . . . , X n ) = 2
E
2
( + 2)
2
n
E
l (, ; X 1 , . . . , X n ) =
( + 1)
n
2
I (, ) =
n
( + 1)
n
( + 1)
.
n
2 ( + 2)
h=
h(, ) =
h(, )
1
.
2 1
Un intervalle de confiance approximatif de niveau 1 pour le mode de la distribution Pareto inverse est donc
q
1
h,
)
z /2 h T I 1 (,
2
o
5.6
h=
.
2 1
Lestimation par distance minimale est une autre technique destimation paramtrique base sur un critre optimiser. Lide, ici, consiste choisir les paramtres dune distribution de manire minimiser une certaine distance entre la
fonction de rpartition thorique F (x; ) et la fonction de rpartition empirique
5.6.1
Soit y 1 < y 2 < < y k des points arbitraires et w 1 , . . . , w k des nombres positifs
arbitraires qui serviront de poids. La forme gnrale de la statistique (ou distance)
de Cramrvon Mises est
d () =
k
X
w j [F (y j , ) F n (y j )]2 .
j =1
Les poids permettent de donner plus dimportance certains points par rapport
dautres. En gnral, on choisit des poids gaux, soit w 1 = = w k = 1. Si lon
souhaite donner plus de poids aux donnes dans les queues de la distribution, un
choix populaire est
n
wj =
.
F (y j )(1 F (y j ))
La figure 5.1 reprsente cette fonction pour une distribution gamma.
Pour des donnes individuelles, on prend y j = x j , j = 1, . . . , n. Lestimateur de
Cramr-von Mises de est donc la valeur qui minimise
d () =
n
X
w j [F (x j , ) F n (x j )]2 .
j =1
Pour des donnes groupes, on minimise la distance aux bornes des classes,
cest--dire que y j = c j , j = 1, . . . , r . On a alors
d () =
r
X
w j [F (c j , ) F n (c j )]2 .
j =1
213
500000
1000000 1500000
Modlisation paramtrique
wj
214
0.0
0.5
1.0
1.5
2.0
2.5
3.0
yj
F IG . 5.1: Fonction de poids donnant plus dimportance aux donnes dans les
queues
et
0,
1/8,
2/8,
3/8,
F 8 (x) = 4/8,
5/8,
6/8,
7/8,
1,
x <1
x =1
x =3
x =7
x = 11 ,
x = 22
x = 23
x = 44
x 51
7
d () =
e
8
6
+ e 3
8
+ + (e 51 )2 .
5.6.2
215
Cette distance est base sur une statistique utilise dans les tests dadquation
des modles. Elle ne semploie que pour les donnes groupes (ou aprs avoir
group des donnes individuelles). Lestimateur est obtenu en minimisant
d () =
r
X
w j [n(F (c j ; ) F (c j 1 ; )) n j ]2 ,
j =1
Pr
o n =
Soit
j =1 n j
r (E () n )2
X
j
j
j =1
nj
On voit ainsi plus clairement que les paramtres sont choisis de manire minimiser lcart entre le nombre observ de donnes dans chacune des classes de
lchantillon et le nombre de donnes selon le modle thorique.
Dans la mesure o le terme n j se retrouve au dnominateur, il nest pas possible
davoir des classes vides. Si le cas se prsente, il faut regrouper de manire arbitraire
certaines classes de faon avoir au moins une donne dans chacune des classes.
5.6.3
r
X
j =1
o
LAS(c j 1 , c j ; ) = E [min(X , c j ); ] E [min(X , c j 1 ); ]
LASn (c j 1 , c j ) = E [X ; c j ] E [X ; c j 1 ]
216
Modlisation paramtrique
et w 1 , . . . , w r sont des poids arbitraires. Il sagit donc dune comparaison entre les
sinistres moyens par classe entre le modle thorique et le modle empirique.
5.7
Code informatique
###
### ESTIMATION PAR LA MTHODE DU MAXIMUM DE VRAISEMBLANCE
###
## EXEMPLE 5.9
##
## Estimation du maximum de vraisemblance pour les paramtre d'une
## distribution gamma avec 'optim'. Aux fins de l'valuation
## numrique, on simule un chantillon alatoire d'une loi gamma.
x <- rgamma(100, 5, 2)
## On minimise la fonction de log-vraisemblance ngative directement
## avec la fonction 'optim'.
##
## Les principaux arguments de 'optim' sont:
##
## par: un vecteur contenant les valeurs initiales des paramtres;
##
fn: fonction minimiser. Le premier argument de fn doit tre
##
le vecteur des paramtres.
##
## On peut passer des valeurs 'fn' (les donnes, par exemple) par le
## biais de l'argument '...' de 'optim'.
##
## Note: l'argument 'log = TRUE' de 'dgamma' permet de calculer plus
## prcisment le logarithme de la densit.
f <- function(p, x) -sum(dgamma(x, p[1], p[2], log = TRUE))
optim(c(1, 1), f, x = x)
##
##
##
##
##
EXEMPLE 5.10
Estimation du maximum de vraisemblance pour les paramtre d'une
distribution gamma avec 'fitdistr'. La fonction se trouve dans le
package MASS.
library(MASS)
## La distribution gamma tant reconnue par 'fitdistr', il suffit d'en
## spcifier le nom en argument. La fonction prparera la fonction de
## log-vraisemblance et choisira des valeurs de dpart appropries.
fitdistr(x, "gamma")
## EXEMPLE 5.11
##
## On estime les paramtres d'une loi Pareto par la mthode du maximum
## de vraisemblance pour l'chantillon simul ci-dessous avec la
## fonction 'rapreto' du package actuar.
library(actuar)
# charger le package
set.seed(2)
# toujours le mme chantillon
x <- rpareto(100, 5, 4000)
## Les estimateurs des moments de $\alpha$ et $\lambda$ seront
## utiliss comme valeurs de dpart dans 'fitdistr'.
m <- emm(x, c(1, 2))
# deux premiers moments
(alpha.mme <- 2 * (m[2] - m[1]^2)/(m[2] - 2 * m[1]^2))
(lambda.mme <- (prod(m))/(m[2] - 2 * m[1]^2))
## Estimation avec 'fitdistr'. Pas de problme de convergence ici.
fitdistr(x, dpareto,
start = list(shape = alpha.mme,
scale = lambda.mme))
## Lorsque la fonction d'optimisation s'gare dans les valeurs
## ngatives des paramtres, on peut utiliser le truc suivant: plutt
## que d'estimer les paramtres eux-mmes, on estime leurs
## logarithmes. Ceux-ci demeurent alors valides sur tout l'axe des
## rels.
f <- function(x, lshape, lscale)
dpareto(x, exp(lshape), exp(lscale))
(fit <- fitdistr(x, f,
start = list(lshape = log(alpha.mme),
lscale = log(lambda.mme))))
## Les valeurs obtenues ci-dessus sont les logarithmes des paramtres.
exp(unname(fit$estimate))
## Dans le cas prsent, on peut rduire le problme de maximisation en
## deux dimensions un problme de recherche de racine en une seule
## dimension. On dfinit d'abord la fonction dont on cherche la racine.
g <- function(b, x)
217
218
Modlisation paramtrique
{
n
x
u
n
<- length(x)
<- x + b
<- b * sum(1/x)
* (n - u)/u + n * log(b) - sum(log(x))
}
## On trouve ensuite la racine avec la fonction 'uniroot' sur
## l'intervalle (0, max(x)).
(lambda <- uniroot(g, lower = 1E-10, upper = max(x), x = x)$root)
(alpha <- 1/(-1 + length(x)/(lambda * sum(1/(x + lambda)))))
## EXEMPLE 5.13
##
## On a des donnes de contrats d'assurance comportant une franchise
## de 75 et une limite de 800.
x <- c(280.48, 460.26, 60.41, 725.00, 340.34, 166.08, 485.60,
199.15, 45.56, 118.11, 493.15, 239.14, 241.55, 218.62,
725.00, 48.61, 204.77, 710.67, 8.75, 225.37, 123.82,
25.31, 588.26, 27.90, 197.68, 281.83, 370.61, 371.34, 229.94,
240.12, 389.96, 725.00, 100.45, 564.15, 29.40, 40.30,
142.98, 350.30, 41.68, 260.81, 369.12, 134.54, 312.08,
725.00, 253.15, 70.50, 452.59, 14.82, 327.22, 119.42,
81.01, 68.25, 90.84, 725.00, 2.39, 160.67, 86.55,
55.43, 369.39, 107.23, 93.39, 725.00, 725.00, 42.20,
138.52, 191.94, 266.05, 372.72, 52.36, 111.30)
## On utilise 'coverage' pour obtenir la densit des donnes
## modifies.
f <- coverage(dpareto, ppareto, deductible = 75, limit = 800)
## On utilise la fonction obtenue ci-dessus dans 'fitdistr' comme
## n'importe quelle autre fonction de densit.
fitdistr(x, f, start = list(shape = 1, scale = 1), method = "N")
## EXEMPLE 5.14
##
## On a des donnes groupes.
(x <- grouped.data(cj = c(0, 1, 2, 4, 8, 12, 16, 20, Inf),
nj = c(4, 3, 3, 4, 2, 0, 1, 3)))
##
##
##
##
##
##
##
##
ll
219
220
Modlisation paramtrique
##
## Le calcul d'un l'intervalle de confiance 95
## quelques lignes de programmation.
(m <- fit$estimate[2] * (fit$estimate[1] - 1)/2)
h <- c(fit$estimate[2], fit$estimate[1] - 1)/2
s <- sqrt(drop(crossprod(h, fit$vcov)) %*% h)
c(m - 1.96 * s, m + 1.96 *s)
221
mode estim
gradient estim
cart-type estim
intervalle
###
### ESTIMATION PAR DISTANCE MINIMALE
###
## L'interface de la fonction 'mde' de actuar est trs similaire
## celle de 'fitdistr'. Outre les donne, la fonction requiert en
## argument:
##
## - la mesure de distance utiliser ("CvM", "chi-square" ou "LAS");
## - une fonction pour calculer la fonction de rpartition thorique
##
ou l'esprance limite thorique;
## - des valeurs de dpart;
## - des poids, le cas chant.
##
## On va ajuster par la mthode de la distance minimale une
## distribution exponentielle aux donnes groupes d'assurance
## dentaire.
gdental
## La moyenne des donnes nous donnera une ide de la valeur de dpart
## fournir la fonction de minimisation de distance.
mean(gdental)
## Estimateur de la distance de Cramr-von Mises. On peut utiliser des
## donnes individuelles ou des donnes groupes avec cette distance.
## On fournit 'mde' la fonction de rpartition thorique; 'mde' se
## chargera de construire l'ogive et la fonction de distance
## minimiser, puis d'appler 'optim' pour effectuer la minimisation
## numrique.
mde(gdental, pexp, start = list(rate = 1/200), measure = "CvM")
## Pour viter l'avertissement lors d'une optimisation sur un seul
## paramtre, on peut utiliser la mthode de type quasi-Newton "BFGS".
mde(gdental, pexp, start = list(rate = 1/200),
measure = "CvM", method = "BFGS")
## Estimateur de la distance du khi carr modifi.
222
Modlisation paramtrique
5.8
Exercices
5.1 Soit X , une variable alatoire reprsentant le montant dun sinistre. On suppose
X | = Exponentielle()
Gamma(, ).
Les sinistres suivants ont t observs :
{1, 10, 200, 1 000, 5 000}.
Estimer et par la mthode des moments.
5.2 On dispose dun chantillon alatoire avec deux donnes infrieures 2 000
et quatre donnes entre 2 000 et 5 000. Les donnes suprieures 5 000 nont
pas t enregistres. crire la fonction de vraisemblance pour un modle de loi
exponentielle.
5.8. Exercices
223
0 < x < 1.
12 ( x )2
x > 0, > 0.
p
Lesprance de cette variable alatoire est donne par 2/2. On a observ
les cinq valeurs suivantes :
{4,9, 1,8, 3,4, 6,9, 4,0}.
Dterminer lestimateur de laide de la mthode des moments.
5.8 On suppose que la distribution du montant des sinistres obit une loi Weibull(, )
de paramtres inconnus.
a) Sachant que 50 % des sinistres sont suprieurs 1 000 $ et que 75 % des
sinistres sont suprieurs 500 $, estimer et par la mthode des quantiles.
b) partir des estimations obtenues en a), estimer le 80e centile.
224
Modlisation paramtrique
5.9 Soit X , la variable alatoire reprsentant le montant dun sinistre. On suppose que le montant dun sinistre pour un fix obit une distribution
Exponentielle() et que est une ralisation de la variable alatoire , o
Gamma(, ). la suite dune exprience, on observe que 0,1 % des sinistres sont suprieurs 450 et que 87,5 % des sinistres sont infrieurs 50.
Trouver lquation, uniquement fonction de , que lon doit rsoudre pour
estimer et qui, aprs avoir t rsolue, permet destimer le paramtre .
5.10 Pour des contrats en assurance automobile avec les modalits suivantes, on a
observ pour lanne 1999 :
1 , a <x <b
f X (x) = b a
0,
ailleurs.
5.12 On a enregistr n essais indpendants X 1 , . . . , X n de la variable alatoire X
Bernoulli(p). Trouver lestimateur du maximum de vraisemblance pour p.
5.13 Soit X 1 , . . . , X n , un chantillon alatoire provenant dune loi normale de paramtres et 2 inconnus.
a) Trouver les estimateurs du maximum de vraisemblance de et 2 .
2 ont approximativement une distribution normale
b) Dmontrer que et
conjointe avec moyennes et 2 et variances 2 /n et 24 /n.
2 ) de
c) Trouver lapproximation de la distribution de lestimateur h(,
c
h(, 2 ) = Pr[X c] =
.
5.8. Exercices
225
5.14 Soit X , une variable alatoire reprsentant les montants de sinistres dont on
possde un chantillon de taille n. La fonction de densit de probabilit de X
est
2
f (x) = 2xe x , x > 0.
Dterminer lestimateur du maximum de vraisemblance de .
5.15 Un assureur possde un chantillon alatoire x 1 , . . . , x n et il souhaite modliser la variable alatoire sous-jacente laide de la fonction
F (x) = x p ,
0 < x < 1.
x > 0.
f (x) = 2xe x ,
x > 0.
La seule information dont on dispose est quune des quatre observations est
infrieure 2. Calculer une estimation du maximum de vraisemblance de .
5.18 Un chantillon de taille 40 a t tir dune population dont la densit est
f (x) = (2)1/2 e x
/(2)
< x < .
226
Modlisation paramtrique
5.19 On suppose que X obit une distribution log-gamma :
f (x) =
2 ln(x)
x +1
x > 1.
5.21 Le tableau ci-dessous prsente les sinistres pays en 1999. On pose lhypothse que la svrit dun sinistre est distribue selon une loi de Pareto de
paramtres et 1. Dterminer lquation finale permettant de trouver lestimateur du maximum de vraisemblance de .
Montant
Nombre de sinistres
(0, 2]
(2, 5]
(5, 11]
(11, )
2
0
1
1
5.22 Le tableau ci-dessous prsente les sinistres pays par un assureur. On pose
que la distribution de X est une exponentielle de paramtre inconnu. Quel
est lestimateur du maximum de vraisemblance de ?
Montant
Nombre de sinistres
(0, 1]
(1, 2]
(2, )
1
0
1
f (x) = 2xe x ,
x > 0.
5.8. Exercices
227
a) Dterminer Pk .
b) Dterminer la variance de lestimateur trouv en a).
c) Si X 1 = X 2 = 10 et X 3 = 15, calculer Pr[P10 21 ].
5.24 Sachant quun chantillon alatoire X 1 , . . . , X 50 provenant dune distribution
de Pareto(, ) a conduit aux estimations = 1,5 et = 1 500 par la mthode
du maximum de vraisemblance, estimer les variances des estimateurs et
ainsi que leur covariance.
5.25 On suppose que le montant dun sinistre obit une loi de Pareto(, ). Pendant une anne, on observe 50 sinistres. laide des montants des 50 sinistres,
= 40. Si la covariance entre les es = 24 et Var[]
on obtient = 2, = 4, Var[]
228
Modlisation paramtrique
b) Calculer partir de lchantillon {2, 1, 2, 3, 3, 4} si la fonction de perte
choisie est lerreur quadratique.
5.30 On suppose que X |B = Exponentielle() et que la distribution a priori de
B est une Gamma(2, 3). On a lchantillon alatoire suivant :
{6, 11, 8, 13, 9}
a) Calculer lestimateur bayesien du paramtre si la fonction de perte est
lerreur quadratique.
b) Rpter la partie a) avec la fonction de perte valeur absolue de lerreur. On
fournit les valeurs
(7; 4,734) = 0,2
5.31 Au cours dune session, les tudiants en actuariat font des devoirs informatiques. En faisant ces devoirs, il leur arrive de rester bloqus. Le nombre
de fois o un tudiant reste bloqu dans un devoir suit une distribution
Binomiale(3, ), o lon suppose que est uniformment distribu sur lintervalle (0,25, 0,75). Deux tudiants sont rests bloqus chacun deux fois pendant
un certain devoir.
a) Trouver lestimateur bayesien de avec une fonction de perte quadratique.
b) Dterminer la probabilit a posteriori que se retrouve dans lintervalle
(0,6, 0,7).
5.32 Pour des contrats dassurance comportant une rtention de 1,5 millions, 40
catastrophes ont t dclares au rassureur. Le rassureur suppose que les
montants de sinistres obissent une loi de Pareto(, ). Soit Y la variable
alatoire reprsentant un montant de sinistre dclar au rassureur (en millions). laide des montants qui lui ont t dclars, le rassureur a estim les
paramtres et par la mthode du maximum de vraisemblance. Il a obtenu
= 5,084 et = 28,998.
a) Trouver, par la mthode du maximum de vraisemblance, lestimation de
Pr[Y > 29,5].
est
)
b) Si la matrice variance-covariance de (,
23,92
167,07
,
167,07 1 199,32
5.8. Exercices
229
f (x) =
(
e x , x > 0
0,
ailleurs.
230
Modlisation paramtrique
Rponses
5.1 = 3,45, = 3 048,87
5.2 L() = [(1 e 2 000 )2 (e 2 000 e 5 000 )4 ]/(1 e 5 000 )6
5.3 = 34,83, = 27,99
5.4 = 2, = 200
5.5 x/(1
x)
= 0,6368 et 0,056
5.6 = 7,40,
5.7 3,3511
5.8 a) = 1,2687, = 0,000747 b) 1 947
5.9 ( + 450)0,3010 = 0,3010 ( + 50)
5.10 = 0,48, = 0,01,
5.11 a = 10, b = 60
5.12 p = X
2 = S 2 b) h(,
2 ) N (h(, 2 ),V ), V = 2 ((c )/
)(1/n
5.13 a) = X ,
+ (c
2
2
/(2n
))
)
P
5.14 n/ ni=1 x i2
p
P
+ p)
+ p)
+
e) p/(1
1,96p(1
5.15 a) n/ ni=1 ln x i b) p 2 /n c) p 1,96p/ n d) p/(1
p
2
/ n
p)
5.16 3,8629
5.17
4
1
4 ln 3
5.18 0,20
p
p
P
5.19 a) X /( X 1) b) 2n/ ni=1 ln(X i )
5.20 a) 1/2 b) 1/64
5.21 L() = (1 (1/3) )2 ((1/6) (1/12) )(1/12)
5.22 ln(1,5)
P
2
2
5.23 a) 1 e k , = n/ ni=1 X i2 b) k 4 2 e 2k /n c) 0,4875
5.8. Exercices
5.27 a) 0,4 b) 0,3432
P
5.28 ( + ni=1 X i )/( + n)
P
5.29 a) Gamma n + 1, 3 + ni=1 ln(1 + x i ) b) 0,68
5.30 a) 0,14 b) 0,1334
5.31 a) 0,5668 b) 0,3055
5.32 a) 0,0365 b) 0,00057
5.33 1 302,50
5.34 0,0059
5.35 0,2286
231
6 Tests dadquation
Objectifs du chapitre
x valuer ladquation dun modle de svrit des sinistres aux donnes laide
de divers tests graphiques (comparaison des distributions empirique et thorique,
graphique quantile-quantile) et statistiques (KolmogorovSmirnov, khi carr, ratio
des vraisemblances).
x Calculer la valeur du critre bayesien de Schwarz et du critre dinformation de
Akaike.
x Choisir le modle le plus appropri partir des rsultats des tests et critres cidessus.
6.1
Tests graphiques
Lide ici est de comparer graphiquement une caractristique du modle postul avec la caractristique empirique correspondante.
6.1.1
234
Tests dadquation
empirique. videmment, il sagit ici dun critre trs subjectif qui permet seulement
une classification plutt grossire en bons, acceptables et mauvais modles.
Exemple 6.1. Soit lchantillon
{30, 38, 34, 29, 35, 34, 29, 39, 20, 37,
21, 27, 37, 18, 37, 34, 17, 28, 38, 32}.
On ajuste des distributions gamma et log-normale ces donnes par la mthode
du maximum de vraisemblance :
> x <- c(30, 38, 34, 29, 35, 34, 29, 39, 20, 37,
+
21, 27, 37, 18, 37, 34, 17, 28, 38, 32)
rate
17.2322211
0.5613102
sdlog
3.394966 0.252984
6.1.2
On trouvera le code R pour crer les graphiques des figures 6.1 et 6.2
la section 6.4.
Graphique quantile-quantile
235
1.0
ecdf(x)
0.8
Fn(x)
0.6
0.4
0.2
0.0
15
20
25
30
35
40
j +a
1
yj ,F
, j = 1, . . . , n,
n + 1 2a
236
Tests dadquation
0.04
0.00
0.02
Density
0.06
Histogram of x
10
20
30
40
237
25
20
15
Sample Quantiles
10
10
15
20
25
Theoretical Quantiles
> x <- c(10, 12, 13, 15, 18, 20, 22, 23, 29)
> (fit.e <- fitdistr(x, "exponential")$estimate)
rate
0.05555556
> (fit.ln <- fitdistr(x, "lognormal")$estimate)
meanlog
sdlog
2.837825 0.327125
238
Tests dadquation
0.06
0.00
0.02
0.04
Density
0.08
0.10
Histogram of x
10
15
20
25
30
6.2
On trouvera le code R pour crer les graphiques des figures 6.3 et 6.4
la section 6.4.
Tests statistiques
239
6.2.1
Test de Kolmogorov-Smirnov
Tests dadquation
0.8
240
0.7
Fn(x) F0(x)
0.6
0.4
0.5
Fn(x)
Fn(x) F0(x)
0.6
0.8
1.0
1.2
1.4
1.6
241
0,20
0,10
0,05
0,02
0,01
p
1,073/ n
p
1,223/ n
p
1,358/ n
p
1,518/ n
p
1,629/ n
x > 0,
242
Tests dadquation
xj
|F n (x j ) F 0 (x j )|
|F n (x j ) F 0 (x j )|
2
5
7
10
11
12
18
0,0514
0,1315
0,1019
0,0890
0,0191
0,1308
0,1431
0,1943
0,2744
0,2447
0,2318
0,1237
0,0121
0,0003
0,
1/7,
2/7,
3/7,
F 7 (x) =
4/7,
5/7,
6/7,
1,
x <2
x 2
x 5
x 7
x 10
x 11
x 12
x 18.
Le tableau 6.2 prsente les diffrences entre les fonctions de rpartition. La valeur de
la statistique D 7 est la plus grande diffrence de ce tableau. On a donc D 7 = 0,2744.
p
La valeur critique du test est c = 1,073/ 7 = 0,4056 > 0,2744. On ne rejette donc
pas le modle.
Ce rsultat est confirm par les valeurs p asymptotique et exacte telles que
calcules par la fonction ks.test :
> x <- c(2, 5, 7, 10, 11, 12, 18)
> ks.test(x, pexp, rate = 0.108, exact = TRUE)
One-sample Kolmogorov-Smirnov test
data:
243
1.0
0.8
0.6
Fn(x)
0.4
0.2
0.0
10
15
20
La figure 6.6 prsente une comparaison entre la fonction de rpartition empirique et la fonction de rpartition thorique. Les lignes pointilles indiquent
les bornes de lintervalle de confiance calcul avec la statistique de Kolmogorov
Smirnov. On constate visuellement que lhypothse nulle nest pas rejete puisque
la fonction de rpartition empirique se trouve lintrieur des bornes plus de
80 %.
244
Tests dadquation
6.2.2
Le test du khi carr est valide uniquement pour des donnes groupes ; on
regroupera donc arbitrairement les donnes individuelles. Soit n j le nombre de
P
donnes dans la classe (c j 1 , c j ], j = 1, . . . , r et rj =1 n j = n. La statistique est
Q=
r (E n )2
X
j
j
j =1
Ej
o
F (c j 1 ; )).
E j = n(F (c j ; )
Asymptotiquement, la statistique Q obit une loi 2 (r p 1), o p est le
nombre de paramtres estims. On rejette donc H0 si, et seulement si, la valeur de
la statistique est suprieure au 100(1 )e quantile dune telle loi.
On considre gnralement que ce test est valide si le nombre espr de donnes dans chacune des classes est suprieur 5. Si ce nest pas le cas, il faut regrouper des classes.
Exemple 6.4. Supposons que lon a les donnes groupes suivantes :
Intervalle
Frquence
(0 5]
(5 20]
(20 )
60
25
15
x
10
=
10 + x 10 + x
et
E 1 = 100(F (5) F (0)) = 33,33
E 2 = 100(F (20) F (5)) = 33,33
E 3 = 100(1 F (20)) = 33,33,
do la statistique du khi carr est
(60 33,33)2 (25 33,33)2 (15 33,33)2
+
+
33,33
33,33
33,33
= 33,50.
Q=
6.2.3
245
Le test du ratio des vraisemblances (likelihood ratio test ; LRT) est un test statistique bas sur la fonction de vraisemblance permettant de vrifier si un paramtre
dune distribution appartient ou non un sous-espace 0 de lespace des valeurs
possibles des paramtres. La forme gnrale du test est la suivante :
H0 : 0
H1 : c0 ,
o c0 est le complment de 0 . La statistique du test est
sup0 L(; X 1 , . . . , X n )
T = 2 ln
,
sup L(; X 1 , . . . , X n )
Q
o L(; X 1 , . . . , X n ) = ni=1 f (X i ; ) est la fonction de vraisemblance de lchantillon
X 1 , . . . , X n . Par les proprits des logarithmes, on peut aussi crire la statistique
ainsi :
"
#
T = 2 sup l (; X 1 , . . . , X n ) sup l (; X 1 , . . . , X n ) ,
0
L(; X 1 , . . . , X n ) L() =
e X i = n e n X .
i =1
Or, on sait (voir lexemple 5.3) que cette fonction trouve son maximum dans R+
lorsque = X 1 . La statistique du test LRT est donc, ici,
L(0 )
T = 2 ln
L( X 1 )
!
n0 e 0 n X
= 2 ln
X n e n
= 2n(ln 0 + ln X )n 2n + 20 n X .
246
Tests dadquation
Le 95e centile dune distribution 2 avec un degr de libert est 3,8415. On rejetterait
donc lhypothse H0 si T > 3,8415 pour un chantillon et une valeur de 0 donns.
Exemple 6.6. On a un chantillon alatoire de 20 valeurs dune distribution de Pareto de paramtre = 2 et inconnu. Lestimateur du maximum de vraisemblance
P
P20
de est 7,0. On sait galement que 20
i =1 ln(x i + 7,0) = 49,01 et i =1 ln(x i + 3,1) =
39,30. On souhaite tester
H0 : = 3,1
H1 : 6= 3,1.
On vrifie aisment que la fonction de log-vraisemblance est
l () = n ln + n ln ( + 1)
n
X
ln(x i + )
i =1
= 20 ln 2 + 40 ln 3
20
X
ln(x i + ).
i =1
maximum de vraisemblance. Si lon note 0 lestimateur du maximum de vraisemblance sur le sous-espace restreint 0 , la statistique du test peut scrire plus
simplement ainsi :
!
L(0 ; X 1 , . . . , X n )
T = 2 ln
X1, . . . , Xn )
L(;
X 1 , . . . , X n )].
= 2[l (0 ; X 1 , . . . , X n ) l (;
Dans notre contexte, il nest pas trs intressant de tester si un paramtre est
gal ou non une certaine valeur. Nous sommes davantage intresss comparer lajustement de deux distributions. Le test du ratio des vraisemblances nous
247
permet de faire cette comparaison, mais dans la mesure o une distribution est
un cas spcial de lautre. En effet, la maximisation de la fonction de vraisemblance
sur un sous-espace restreint sapparente alors fixer la valeur dun ou plusieurs
paramtres de la distribution gnrale. La distribution de la statistique est toujours
une 2 , mais le nombre de degrs de libert est gal la diffrence entre le nombre
de paramtres libres sous H1 et le nombre de paramtres libres sous H0 .
Exemple 6.7. La distribution de Pareto est un cas spcial de la distribution de
Pareto gnralise lorsque = 1. Ainsi, le test H0 : les donnes proviennent dune
distribution de Pareto(, ) versus H1 : les donnes proviennent dune distribution de Pareto gnralise(, , ) se ramne essentiellement au test du ratio des
vraisemblances
H0 : = 1
H1 : 6= 1,
avec l (, , ; X 1 , . . . , X n ) la fonction de log-vraisemblance de la distribution Pareto
gnralise. Soit 0 et 0 les estimateurs du maximum de vraisemblance des paramtres de la distribution de Pareto et soit 1 , 1 et 1 , ceux des paramtres de la
distribution de Pareto gnralise. La statistique du test ci-dessus est alors
T = 2[l ( 0 , 1, 0 ; X 1 , . . . , X n ) l ( 1 , 1 , 1 ; X 1 , . . . , X n )].
Sa distribution est une 2 avec 3 2 = 1 degr de libert.
6.3
6.3.1
p
ln n,
2
248
Tests dadquation
6.3.2
6.4
Code informatique
###
### TESTS GRAPHIQUES
###
## EXEMPLE 6.1
##
## Comparaison des courbes
## un petit chantillon de
x <- c(30, 38, 34, 29, 35,
21, 27, 37, 18, 37,
## EXEMPLE 6.2
##
## On ajuste des distributions exponentielle et log-normale un petit
## chantillon par la mthode du maximum de vraisemblance et on
## compare leur ajustement avec des graphiques quantile-quantile.
x <- c(10, 12, 13, 15, 18, 20, 22, 23, 29)
(fit.e <- fitdistr(x, "exponential")$estimate)
(fit.ln <- fitdistr(x, "lognormal")$estimate)
## La fonction 'ppoints' calcule les points o l'on valuera la
## fonction de quantile thorique. Par dfaut, a = 3/8 si le nombre de
## donnes est infrieur ou gal 10 (comme c'est le cas ici) et a =
## 1/2 sinon. On peut fixer la valeur de a nous-mmes.
(y <- ppoints(length(x)))
# avec 'ppoints'
(seq_along(x) - 3/8)/(length(x) + 0.25) # idem
## La fonction 'qqnorm' trace un graphique quantile-quantile pour des
## donnes normales. De manire plus gnrale, 'qqplot' trace un
## graphique quantile-quantile entre deux jeux de donnes. Pour nos
## besoins, il suffit de lui fournir les donnes et les quantiles
## thoriques.
##
## On commence avec le graphique de la log-normale.
qqplot(qlnorm(y, fit.ln[1], fit.ln[2]), x,
pch = 16, col = "darkblue",
xlab = "Theoretical Quantiles",
ylab = "Sample Quantiles")
## Pour ajouter le graphique de l'exponentielle au graphique existant,
## on doit utiliser la fonction 'points' et renverser l'ordre des
## arguments.
points(qexp(y, fit.e), x, pch = 16, col = "darkred")
## Enfin, ajout de la droite y = x.
abline(0, 1)
## L'ajustement de la log-normale est bien meilleur, comme en fait foi
## la comparaison des densits avec l'histogramme.
hist(x, breaks = c(10, 12, 18, 25, 30))
249
250
Tests dadquation
6.5
Exercices
Nombre de donnes
(0, 3]
(3, 7,5]
(7,5, 15]
(15, 40]
(40, )
180
180
235
255
150
Une loi de Pareto a t ajuste ces donnes et les estimateurs obtenus sont
= 3,5 et = 50. Quel est le seuil de signification le plus lev (parmi 5 %,
2,5 %, 1 % et 0,5 %) auquel on ne rejette pas ce modle avec le test du khi carr ?
6.3 On dispose de lchantillon alatoire {0,1, 0,4, 0,8, 0,8, 0,9} et on veut y ajuster
la distribution avec fonction de densit de probabilit
f (x) =
1 + 2x
,
2
0 x 1.
6.5. Exercices
251
Nombre de sinistres
(0, 2]
(0, 4]
(0, 8]
6
9
10
f (x) =
0 x 2,
et soit lchantillon tir de cette densit {0,5, 1, 1,25, 1,5}. Calculer la statistique
de KolmogorovSmirnov.
6.6 On veut tester si
(
f X (x) =
x
50 ,
0 < x < 10
0,
ailleurs
Frquence
10
5
5
252
Tests dadquation
pour guider le choix de la distribution. Voici quelques valeurs de la Gamma
incomplte : (3,5; 1,25) = 0,0729, (3,5; 5,51) = 0,8614, (3,5; 7) = 0,9488. De
plus, pour entier, on a
(; x) = 1
1
X
x j e x
.
j!
j =0
6.9 On a observ les sinistres du tableau ci-dessous en assurance mdicaments. Dterminer, laide de la statistique de Pearson, si lhypothse dune distribution
avec taux dchec constant
(x) = 0,01,
x >0
Nombre de sinistres
[0, 25)
[25, 40)
[40, 60)
[60, 80)
[80, )
10
5
10
5
20
Frquence
10
12
12
11
5
On hsite entre une loi de Pareto(1,5, 50) et une loi de Weibull(0,01, 1) pour la
distribution du montant dun sinistre.
a) Quel modle privilgier si on utilise la distance de Cramrvon Mises avec
poids unitaires pour guider le choix ?
b) Si la statistique de Pearson avait t utilise au lieu de la distance de
Cramrvon Mises, lhypothse de la loi Pareto(1,5, 50) aurait-elle t rejete un niveau de confiance = 0,05 ?
6.5. Exercices
253
Gains
29
19
18
25
17
10
15
11
Nombre de paramtres
Log-vraisemblance
3
3
2
2
1
219,1
219,2
221,2
221,4
224,4
Rponses
6.1 Q = 0,7740
6.2 0,5 %
254
Tests dadquation
6.3 D = 0,32
6.4 a) 0,3478 b) 0,0242
6.5 0,4375
6.6 D = 0,1329
6.7 1,1667
6.8 Gamma(3,5, 0,01)
6.9 Q = 1,8179
6.10 a) Weibull b) oui c) oui
6.11 a) D = 0,132
6.12 Pareto
7 Modles de frquence
Objectifs du chapitre
x Dfinir et connatre les principales caractristiques (fonction de masse de probabilit, fonction de rpartition, moments) des modles de frquence les plus usuels en
actuariat : loi de Poisson, loi binomiale ngative, loi gomtrique, loi binomiale.
x Estimer les paramtres des lois mentionnes ci-dessus et choisir quel modle utiliser
dans un contexte donn.
x Dfinir les familles de lois de probabilit discrtes (a, b, 0) et (a, b, 1) et en identifier
les principaux membres.
7.1
k = 0, 1, 2, . . .
pk t k .
k=0
On a, par exemple,
P N (e t ) = E [e t N ] = M N (t )
255
256
Modles de frquence
et
P N (1) = E [1N ] = 1.
Comme son nom lindique, la fonction gnratrice des probabilits permet de
gnrer les probabilits dune variable alatoire. En effet, la drive dordre k de la
fonction gnratrice des probabilits est
(k)
PN
(z) =
j!
z j k p j
(
j
k)!
j =k
et donc
(k)
PN
(0) = k!p k ,
do
pk =
1 (k)
P (0).
k! N
De plus, on a
(k)
PN
(1) = E [(N )(N 1) (N k + 1)]
et donc
(1)
E [N ] = P N
(1)
(2)
(1)
(1)
Var[N ] = P N
(1) + P N
(1) (P N
(1))2 .
7.2
On prsente dans cette section les principales caractristiques des lois de probabilit les plus souvent utilises pour le dnombrement des sinistres en assurance.
Lannexe A reprend bon nombre des rsultats ci-dessous dans une forme plus
compacte.
7.2.1
Loi de Poisson
k e
,
k!
k = 0, 1, . . . ,
257
k
X
pj
j =0
k j e
X
j!
j =0
= 1 (k + 1; ),
o () est la fonction gamma incomplte.
La fonction gnratrice des moments est
M N (t ) = e (e
1)
> 0.
k
Y
P Ni (z)
i =1
(z1)(1 ++k )
=e
do N Poisson(1 + . . . + k ).
258
Modles de frquence
De manire similaire, le thorme suivant tablit que la distribution de Poisson
est infiniment divisible.
Thorme 7.2. Si on dcompose lensemble des sinistres en m catgories possibles
P
de probabilits respectives q 1 , . . . , q m ( m
i =1 q i = 1) et si Ni est le nombre de sinistres
dans la catgorie i , alors si N Poisson(), les variables alatoires N1 , . . . , Nm sont
indpendantes et telles que Ni Poisson(q i ), i = 1, . . . , m.
Dmonstration. On trouve dabord la distribution des variables alatoires Ni , i =
1, . . . , m. Pour N = n fix, on a Ni |N = n Binomiale(n, q i ). Ainsi, la distribution
conjointe est
N1 , . . . , Nm |N = n Multinomiale(n, q 1 , . . . , q m ).
On a alors
Pr[Ni = n i ] =
Pr[Ni = n i |N = n]Pr[N = n]
n=n i
n=n i
!
n
n ni
nn i e
q (1 q i )
ni i
n!
ni
e q i
ni !
X
n=n i
n i
1
(1 q i )nni n
(n n i )!
nn i (1 q )nn i
X
i
(n
n
i )!
n=n i
(q i ) e
ni !
[(1 q )] j
(q i )ni e X
i
ni !
j
!
j =0
(q i )ni e (1qi )
e
ni !
(q i )ni e qi
,
ni !
do Ni Poisson(q i ).
De plus,
Pr[N1 = n 1 , . . . , Nm = n m ] = Pr[N1 = n 1 , . . . , Nm = n m |N = n]Pr[N = n]
n
n!
nm e
n1
q qm
=
n ! nm ! 1
n!
1
!
!
n 1 q 1
(q 1 ) e
(q m )nm ! e p m
=
.
n1 !
nm !
259
Ainsi, la distribution conjointe est le produit des distributions marginales et, par
consquent, les variables alatoires N1 , . . . , Nm sont stochastiquement indpendantes.
Exemple 7.1. Soit un portefeuille dassurance automobile pour lequel la frquence
des sinistres a une distribution de Poisson de moyenne 1 000. On value que 65 %
des accidents sont causs par des jeunes conducteurs et que 35 % des accidents
sont causs par les conducteurs plus ags. (Il faut noter, ici, que cela ne signifie pas
que le portefeuille contient 65 % de jeunes conducteurs et 35 % de conducteurs
plus vieux.)
Par le thorme prcdent, on a donc les modles suivants pour la distributions des sinistres du portefeuille global (N ), des jeunes conducteurs (N J ) et des
conducteurs plus ags (NV ) :
N Poisson(1 000)
N J Poisson(650)
NV Poisson(350).
7.2.2
Binomiale ngative
La fonction de masse de probabilit de la loi binomiale ngative avec paramtres r > 0 et 0 1 est
Pr[N = k] =
(r + k)
r (1 )k ,
(k + 1)(r )
k = 0, 1, . . .
!
k +r 1
(r + k)
=
.
(k + 1)(r )
r 1
r (1 ) r (1 )
= Var[N ] ,
260
Modles de frquence
avec galit dans le cas peu intressant o = 1.
La distribution binomiale ngative cultive des liens troits avec la distribution
gamma. En effet, la premire est lanalogue discret de la seconde (voir la distribution gomtrique, plus bas). De plus, la binomiale ngative rsulte dun mlange
(continu) entre une distribution de Poisson et la gamma.
7.2.3
Gomtrique
= F X (x + 1) F X (x)
= e x e (x+1)
= e x (1 e )
= (1 )x ,
avec = 1 e . Ainsi, la partie entire dune loi exponentielle est une loi gomtrique. On peut en dduire le lien semblable entre les lois gamma et binomiale
ngative.
7.2.4
Binomiale
261
une normale lorsque n augmente. Enfin, on note que son esprance est suprieure
sa variance. En effet,
E [N ] = n > n(1 ) = Var[N ].
7.3
Estimation et tests
Lestimation des paramtres des lois discrtes nest pas trs diffrent de celle des
lois continues. On utile principalement les mthodes des moments et du maximum
de vraisemblance. ce propos, pour les lois discrtes, la fonction de vraisemblance
reprsente rellement la probabilit dobserver un chantillon donn.
Pour lestimation, on dispose dun chantillon de la variable alatoire N qui
peut tre soit le nombre de sinistres dun portefeuille pour plusieurs priodes, soit
les nombres de sinistres de plusieurs portefeuilles sur une priode. On pose que n
est la taille de lchantillon. Il faut faire attention, car n ne reprsente pas le nombre
total de sinistres.
Soit n k le nombre de fois o la valeur k apparat. On a donc
n=
nk
k=0
P
et le nombre total de sinistres est
kn k .
k=0
On notera que pour la distribution binomiale, le paramtre n est gnralement
considr connu. Il est donc rare quil fasse lobjet dune procdure destimation.
Quant ladquation entre des modles de frquence et les donnes, elle est
gnralement mesure avec le test du khi carr de Pearson. En effet, celui-ci se
prte tout naturellement lexercice puisquil compare des frquences observes
et modlises. De plus, le test de Kolmogorov-Smirnov nest valide que pour des
distributions continues.
Y
k=0
pk k
!
k e n k
Y
=
k!
k=0
262
Modles de frquence
Anne
Frquence
Anne
Frquence
Anne
Frquence
1970
1973
1974
1976
1979
1
1
2
1
1
1980
1983
1988
1989
1990
1
2
2
3
3
1991
1992
1993
1994
1995
3
2
3
1
4
n k ln(p k )
k=0
n k (k ln() ln(k!) ) .
k=0
1 X
=
kn k .
n k=0
Ici, on a n = 26 observations (les annes de 1970 1995, inclusivement) rparties ainsi : n 0 = 11, n 1 = 6, n 2 = 4, n 3 = 4, n 4 = 1 et n k = 0 pour tout k 5. Par
consquent,
30
= .
26
Notre modle pour la frquence annuelle des catastrophes coteuses est donc
N Poisson(30/26).
Pour vrifier ladquation entre ce modle et les donnes, on utilise le test du
khi carr. On rappelle que pour que ce test soit valide, la frquence de chaque
classe devrait tre suprieure 5. Pour ce faire, nous regroupons les valeurs 2, 3 et
4.
Soit
k e
pk =
k!
la probabilit dobserver la valeur k selon le modle et n pk la frquence espre de
cette valeur dans un chantillon de taille n. Le tableau 7.2 prsente les rsultats des
calculs ncessaires pour effectuer le test du khi carr. On a alors
263
pk
n pk
nk
0
1
2+
0,3154
0,3639
0,3297
8,20
9,46
8,34
11
6
9
TAB. 7.2: Quantits ncessaires pour le test du khi carr de lexemple 7.2
Loi
Poisson
Binomiale
Binomiale Ngative
Gomtrique
p0
0
/(1 )
1
1
(n + 1)/(1 )
(r 1)(1 )
0
e
(1 )n
r
Q=
7.4
Soit N une variable alatoire discrte avec une distribution quelconque et,
comme prcdemment, Pr[N = k] = p k . On dit que la distribution de N est membre
de la famille (ou classe) (a, b, 0) sil est possible de trouver deux constantes a et b
telles que
pk
b
= a + , k = 1, 2, 3, . . . ,
p k1
k
avec
p0 = 1
pk .
k=1
On peut dmontrer que les seuls membres de cette famille sont les distributions
binomiale, binomiale ngative, gomtrique et de Poisson. Le tableau 7.3 prsente
les valeurs des paramtres a, b et p 0 pour ces lois.
264
Modles de frquence
On peut rcrire la relation liant les probabilits successives des membres de la
famille (a, b, 0) sous la forme
k
pk
= ak + b,
p k1
k = 1, 2, . . .
On remarque que la partie de droite de cette galit est lquation dune droite de
pente a et dordonne lorigine b. Or, du tableau 7.3, cette pente est nulle pour
une distribution de Poisson, ngative pour une distribution binomiale et positive
pour une distribution binomiale ngative. Cela suggre une procdure graphique
pour choisir quelle distribution devrait tre utilise pour modliser la frquence : il
sagit de dterminer la pente dun graphique de
k
pk
nk
=k
pk1
n k1
en fonction de k.
Un peu plus gnralement, on dit que la distribution de N est membre de la
famille (ou classe) (a, b, 1) sil est possible de trouver deux constantes a et b telles
que
pk
b
= a + , k = 2, 3, . . .
p k1
k
On remarque que la relation rcursive dbute cette fois k = 2 plutt qu k = 1. La
valeur de p 0 est donc quelconque.
Les rapports entre les probabilits successives des membres de cette famille
sont les mmes que ceux des membres de la famille (a, b, 0). Le comportement des
distributions est donc similaire. En revanche, la probabilit zro diffre et permet
une plus grande flexibilit. Si p 0 = 0, on dit que la distribution est zro tronque.
Il peut tre trs utile dans certaines applications de pouvoir exclure la possibilit
de ne pas avoir de sinistres. Lorsque p 0 > 0 (et que la probabilit est diffrente de la
valeur correspondante du tableau 7.3), on dit que la distribution est modifie
zro.
Les sections 6.5 et 6.6 de Klugman et collab. (2012a) offrent plus de dtails sur
les familles (a, b, 0) et (a, b, 1) ainsi que des exemples utiles pour mieux comprendre
le mcanisme de rcursivit.
7.5
Exercices
7.1 Un assureur dcide de modliser la frquence des sinistres par une distribution
N Binomiale(m, ) dont le paramtre m est connu.
a) Dmontrer que lestimateur du maximum de vraisemblance de est sans
biais.
7.5. Exercices
265
Nombre de risques
0
1
2
3
4+
9 048
905
45
2
0
k
( + 1)k+1
k = 0, 1, . . .
!
k +r 1
k
Pr[N = k] =
,
r 1 ( + 1)k+r
k = 0, 1, . . .
266
Modles de frquence
Frquence
Hommes
Femmes
0
1
2
3
4
5+
901
92
5
1
1
0
947
50
2
1
0
0
Nombre de risques
0
1
2
3
4
5
6
7+
861
121
13
3
1
0
1
0
a) Ajuster une distribution Binomiale(7, ) ces donnes en estimant le paramtre par la mthode du maximum de vraisemblance.
b) Ajuster plutt une distribution binomiale ngative aux donnes par la mthode des moments. Utiliser la paramtrisation de lexercice 7.2 c)).
c) Rpter la partie b) en estimant plutt par la mthode du maximum de
vraisemblance.
7.5 Dmontrer que la distribution Binomiale ngative(r, ( + 1)1 ) est le rsultat
du mlange continu de distributions de Poisson suivant
N | = Poisson()
Gamma(r, ).
7.5. Exercices
267
7.6 Un assureur modlise la frquence des sinistres par une distribution Binomiale
ngative(3, 1/6). La svrit des sinistres est modlise par une distribution
Exponentielle(0,01). Si une franchise de 20 $ est ajoute au contrat, calculer
E [N ], lesprance de la frquence modifie.
7.7 Un portefeuille dassurance compte 1 000 contrats. Le tableau ci-dessous rsume linformation connue propos de la frquence des sinistres.
Nombre de sinistres
Nombre de contrats
0
1
2
3
4
5
6
7+
100
267
311
208
87
23
4
0
Nombre de jours
0
1
2
3
4
5
209
111
33
7
3
2
268
Modles de frquence
Rponses
q
)/(mn)
t 1 e t d t
1
(a, b)
Z
0
avec
(a, b) =
t a1 (1 t )b1 d t ,
Z
0
t a1 (1 t )b1 d t =
269
(a)(b)
(a + b)
270
A.1
A.1.1
u (1 u)
,
x(, )
u=
v
,
1+v
v=
F (x) = (, ; u)
k ( + k/)( k/)
, < k <
()()
( + 1/)( 1/)
E [X ; x] =
( + 1/, 1/; u) + x(1 F (x))
()()
E [X k ] =
A.1.2
Burr (, , )
Racine : burr
Paramtres : shape1 (), shape2 (), rate ( = 1/), scale ()
u (1 u)
,
x
F (x) = 1 u
f (x) =
u=
1
,
1+v
v=
k (1 + k/)( k/)
, < k <
()
(1 + 1/)( 1/)
(1 + 1/, 1/; u) + xu
E [X ; x] =
()
E [X k ] =
A.1.3
Burr inverse (, , )
Racine : invburr
Paramtres : shape1 (), shape2 (), rate ( = 1/), scale ()
u (1 u)
,
x
F (x) = u
f (x) =
u=
v
,
1+v
v=
271
k ( + k/)(1 k/)
, < k <
()
( + 1/)(1 1/)
E [X ; x] =
( + 1/, 1 1/; u) + x(1 u )
()
E [X k ] =
A.1.4
Pareto gnralise (, , )
Racine : genpareto
Paramtres : shape1 (), shape2 (), rate ( = 1/), scale ()
f (x) =
u (1 u)
,
x(, )
u=
v
,
1+v
v=
F (x) = (, ; u)
k ( + k)( k)
, < k <
()()
E [X ] =
, >1
1
2 ( + 1)
Var[X ] =
, >2
( 1)2 ( 2)
E [X ; x] =
( + 1, 1; u) + x(1 F (x))
1
E [X k ] =
A.1.5
Pareto (, )
u=
1
,
1+v
v=
k (k + 1)( k)
, 1 < k <
()
E [X ] =
, >1
1
2
Var[X ] =
, >2
( 1)2 ( 2)
1
, 6= 1
1
x +
E [X ; x] =
ln
,
=1
x +
E [X k ] =
272
A.1.6
Pareto inverse (, )
Racine : invpareto
Paramtres : shape (), scale ()
u (1 u)
,
x
F (x) = u
u=
f (x) =
v
,
1+v
v=
k ( + k)(1 k)
, < k < 1
()
Z u
y
E [X ; x] = k
d y + x(1 u )
0 1 y
E [X k ] =
A.1.7
Log-logistique (, )
Racine : llogis
Paramtres : shape (), rate ( = 1/), scale ()
u(1 u)
,
x
F (x) = u
f (x) =
u=
v
,
1+v
E [X k ] = k (1 + k/)(1 k/),
v=
< k <
A.1.8
Paralogistique (, )
Racine : paralogis
Paramtres : shape (), rate ( = 1/), scale ()
2 u (1 u)
,
x
F (x) = 1 u
f (x) =
u=
1
,
1+v
v=
k (1 + k/)( k/)
, 2 < k < 2
()
(1 + 1/)( 1/)
E [X ; x] =
(1 + 1/, 1/; u) + xu
()
E [X k ] =
A.1.9
Paralogistique inverse (, )
Racine : invparalogis
Paramtres : shape (), rate ( = 1/), scale ()
273
u=
v
,
1+v
v=
k ( + k/)(1 k/)
, 2 < k <
()
( + 1/)(1 1/)
E [X ; x] =
( + 1/, 1 1/; u) + x(1 u )
()
E [X k ] =
A.2
A.2.1
Racine : trgamma
Paramtres : shape1 (), shape2 (), rate (), scale ( = 1/)
u e u
,
x()
F (x) = (; u)
f (x) =
E [X k ] =
u = (x)
( + k/)
, k >
k ()
( + 1/)
E [X ; x] =
( + 1/; u) + x(1 (; u))
()
A.2.2
Racine : invtrgamma
Paramtres : shape1 (), shape2 (), rate (), scale ( = 1/)
u e u
,
x()
F (x) = 1 (; u)
f (x) =
E [X k ] =
u = (x)
( k/)
, k <
k ()
( 1/)
E [X ; x] =
(1 ( 1/; u)) + x(; u)
()
A.2.3
Gamma (, )
Racine : gamma
Paramtres : shape (), rate (), scale ( = 1/)
274
E [X k ] =
u = x
( + k)
, k >
k ()
E [X ] =
Var[X ] = 2
( + 1)
E [X ; x] =
( + 1; u) + x(1 (; u))
()
M (t ) =
t
A.2.4
Gamma inverse (, )
Racine : invgamma
Paramtres : shape (), rate (), scale ( = 1/)
u e u
,
u = (x)1
x()
F (x) = 1 (; u)
f (x) =
E [X k ] =
( k)
, k <
k ()
( 1)
E [X ; x] =
(1 ( + 1; u)) + x(; u)
()
A.2.5
Weibull (, )
Racine : weibull
Paramtres : shape (), scale ( = 1/)
ue u
,
x
F (x) = 1 e u
f (x) =
E [X k ] =
(1 + k/)
u = (x)
, k >
k
(1 + 1/)
E [X ; x] =
(1 + 1/; u) + xe u
A.2.6
275
Weibull inverse (, )
u = (x)
f (x) =
E [X k ] =
(1 k/)
, k <
k
(1 1/)
E [X ; x] =
(1 (1 1/; u)) + x(1 e u )
A.2.7
Exponentielle ()
Racine : exp
Paramtre : rate ()
ue u
,
x
F (x) = 1 e u
u = x
f (x) =
E [X k ] =
(k + 1)
k
k > 1
1
Var[X ] = 2
1 e u
E [X ; x] =
M (t ) =
t
E [X ] =
A.2.8
Exponentielle inverse ()
Racine : invexp
Paramtres : rate (), scale ( = 1/)
ue u
,
x
F (x) = e u
u = (x)1
f (x) =
E [X k ] =
(1 k)
k
k <1
276
A.3
A.3.1
Racine : norm
Paramtres : mean ( < < ), sd ()
n 1 x 2 o
1
f (x) = p
, < x <
exp
2
2
Z x
x
2
1
F (x) =
e y d y
, (x) = p
E [X ] =
Var[X ] = 2
M (t ) = e t +
2 2
A.3.2
t /2
Log-normale (, 2 )
Racine : lnorm
Paramtres : meanlog (), sdlog ()
n 1 ln x 2 o
1
f (x) = p
exp
2
2 x
ln x
F (x) =
E [X k ] = e k+k
E [X ] = e +
2 /2
/2
Var[X ] = e 2+ (e 1)
2
A.3.3
Log-gamma (, )
Racine : lgamma
Paramtres : shapelog (), ratelog ()
f (x) =
(ln x)1
x +1 ()
F (x) = (; ln x),
x >1
x >1
277
E [X ] =
k
E [X ] =
1
2
Var[X ] =
2
1
E [X ; x] =
(; ( 1) ln x) + x(1 (; ln x))
1
k
A.3.4
Pareto translate (, )
Racine : pareto1
Paramtres : shape (), min ()
,
x >
x +1
F (x) = 1
,
x >
x
f (x) =
k
, k <
k
E [X ; x] =
1 ( 1)x 1
E [X k ] =
Cette loi est galement appele Pareto un paramtre. Seul est considr comme
un vritable paramtre de la distribution. Le paramtre est la borne infrieure du
support de la distribution et est en gnral considr connu.
A.3.5
Bta gnralise (, , , )
Racine : genbeta
Paramtres : shape1 (), shape2 (), shape3 (), rate ( = 1/), scale ()
f (x) =
u (1 u)1
,
x(, )
u=
0<x <
F (x) = (, ; u)
k ( + )( + k/)
, k >
()( + + k/)
( + )( + 1/)
E [X ; x] =
( + 1/, ; u) + x(1 (, ; u))
()( + + 1/)
E [X k ] =
278
A.3.6
Bta (, )
Racine : beta
Paramtres : shape1 (), shape2 ()
f (x) =
( + ) 1
x
(1 x)1 ,
()()
0<x <1
F (x) = (, ; x)
( + )( + k)
, k >
()( + + k)
E [X ] =
+
Var[X ] =
2
( + ) ( + + 1)
( + )( + 1)
( + 1, ; u) + x(1 (, ; x))
E [X ; x] =
()( + + 1)
E [X k ] =
A.4
A.4.1
Racine : binom
Paramtres : size (n), prob ()
!
n x
Pr[X = x] =
(1 )nx ,
x
E [X ] = n
Var[X ] = n(1 )
M (t ) = (1 + e t )n
P (z) = (1 (z 1))n
A.4.2
Racine : nbinom
Paramtres : size (r ), prob (), mu ( = r (1 )/)
!
x +r 1 r
Pr[X = x] =
(1 )x ,
r 1
0 < < 1, x = 0, 1, . . .
r (1 )
Var[X ] =
2
r
M (t ) =
1 (1 )e t
E [X ] =
P (z) = (1 (1 )z)r
A.4.3
Gomtrique ()
Racine : nbinom
Paramtre : prob ()
Pr[X = x] = (1 )x ,
0 < < 1, x = 0, 1, . . .
1
Var[X ] = 2
E [X ] =
1 (1 )e t
P (z) = (1 (1 )z)1
M (t ) =
A.4.4
Poisson ()
Racine : pois
Paramtre : lambda ()
x e
,
x!
E [X ] =
Pr[X = x] =
Var[X ] =
M (t ) = e (e
1)
P (z) = e (z1)
x = 0, 1, . . .
279
281
282
Si dsir, remplacer la valeur de loption repos par lURL dun autre site miroir
de CRAN.
Les utilisateurs de GNU Emacs voudront ajouter une autre option. Le code
entrer dans le fichier ~/.Rprofile sera plutt
options(repos = "http://cran.ca.r-project.org",
menu.graphics = FALSE)
Consulter la rubriques daide de Startup pour les dtails sur la syntaxe et lemplacement des fichiers de configuration, celles de library et .libPaths pour la gestion
des bibliothques et celle de options pour les diffrentes options reconnues par R.
Aprs un redmarrage de R, la bibliothque personnelle aura prsance sur la
bibliothque principale et il ne sera plus ncessaire de prciser le site miroir de
CRAN lors de linstallation de packages. Ainsi, la simple commande
> install.packages("actuar")
C Table de quantiles
de la loi normale
x
2
1
p e y /2 d y
2
(x) = 1 (x)
Pr[X x] = (x) =
(x)
(x)
(x)
0,00
0,05
0,10
0,15
0,20
0,25
0,30
0,35
0,40
0,45
0,50
0,55
0,60
0,65
0,70
0,75
0,80
0,85
0,90
0,95
1,00
1,05
0,500
0,520
0,540
0,560
0,579
0,599
0,618
0,637
0,655
0,674
0,691
0,709
0,726
0,742
0,758
0,773
0,788
0,802
0,816
0,829
0,841
0,853
1,10
1,15
1,20
1,25
1,282
1,30
1,35
1,40
1,45
1,50
1,55
1,60
1,645
1,65
1,70
1,75
1,80
1,85
1,90
1,95
1,96
2,00
0,864
0,875
0,885
0,894
0,900
0,903
0,911
0,919
0,926
0,933
0,939
0,945
0,950
0,951
0,955
0,960
0,964
0,968
0,971
0,974
0,975
0,977
2,05
2,10
2,15
2,20
2,25
2,30
2,326
2,35
2,40
2,45
2,50
2,55
2,576
2,60
2,65
2,70
2,75
2,80
2,85
2,90
2,95
3,00
0,980
0,982
0,984
0,986
0,988
0,989
0,990
0,991
0,992
0,993
0,994
0,995
0,995
0,995
0,996
0,997
0,997
0,997
0,998
0,998
0,998
0,999
283
D Table de quantiles
de la loi khi carr
x
Pr[X x] =
1
(r /2)2r /2
y r /21 e r /2 d x
Pr[X x]
r
0,01
0,025
0,05
0,95
0,975
0,99
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
0,000
0,020
0,115
0,297
0,554
0,872
1,239
1,646
2,088
2,558
3,053
3,571
4,107
4,660
5,229
5,812
6,408
7,015
7,633
0,001
0,051
0,216
0,484
0,831
1,237
1,690
2,180
2,700
3,247
3,816
4,404
5,009
5,629
6,262
6,908
7,564
8,231
8,907
0,004
0,103
0,352
0,711
1,145
1,635
2,167
2,733
3,325
3,940
4,575
5,226
5,892
6,571
7,261
7,962
8,672
9,390
10,117
3,841
5,991
7,815
9,488
11,070
12,592
14,067
15,507
16,919
18,307
19,675
21,026
22,362
23,685
24,996
26,296
27,587
28,869
30,144
5,024
7,378
9,348
11,143
12,833
14,449
16,013
17,535
19,023
20,483
21,920
23,337
24,736
26,119
27,488
28,845
30,191
31,526
32,852
6,635
9,210
11,345
13,277
15,086
16,812
18,475
20,090
21,666
23,209
24,725
26,217
27,688
29,141
30,578
32,000
33,409
34,805
36,191
285
286
Pr[X x]
r
0,01
0,025
0,05
0,95
0,975
0,99
20
21
22
23
24
25
26
27
28
29
30
8,260
8,897
9,542
10,196
10,856
11,524
12,198
12,879
13,565
14,256
14,953
9,591
10,283
10,982
11,689
12,401
13,120
13,844
14,573
15,308
16,047
16,791
10,851
11,591
12,338
13,091
13,848
14,611
15,379
16,151
16,928
17,708
18,493
31,410
32,671
33,924
35,172
36,415
37,652
38,885
40,113
41,337
42,557
43,773
34,170
35,479
36,781
38,076
39,364
40,646
41,923
43,195
44,461
45,722
46,979
37,566
38,932
40,289
41,638
42,980
44,314
45,642
46,963
48,278
49,588
50,892
Chapitre 1
1.1 On a
lim
x0 2
1
2
et
1 x2 1
= .
x0 2
24 2
lim
La figure E.1 prsente le graphique de la fonction et des deux bornes, ainsi que
le code R pour crer ce graphique.
1.2 Il suffit dappliquer la rgle de lHpital :
lim
x0 ln(x + 1)
= lim
d x/d x
x0 d ln(x + 1)/d x
1
x0 1/(x + 1)
= 1.
= lim
287
288
0.35
0.40
f(x)
0.45
0.50
289
1.3 Il faut faire quelques modifications avant de pouvoir utiliser la rgle de lHpital. On passe dabord la forme logarithmique
y = (1 + x)1/x
ln(y) = ln(1 + x)1/x
=
ln(1 + x)
,
x
x0
x0
x0
x0
1
=e
= e.
1.4 a) On utilise la rgle de lHpital pour valuer
lim
x ln(x)
d x/d x
= lim
x d ln(x)/d x
= lim
x 1/x
= lim x
x
= .
Il est donc possible de conclure que le numrateur tend plus rapidement
vers linfini que le dnominateur, cest--dire que x tend plus rapidement
vers linfini que ln(x).
b) De manire similaire,
lim
x e x
= lim
d x/d x
x d e x /d x
= lim
x e x
= 0,
290
x2 x4
+
...
2!
4!
x3 x5
+
...
3!
5!
c) On obtient
i 2x2 i 3x3 i 4x4 i 5x5
+
+
+
+...
2!
3!
4!
5!
x2
x3 x4
x5
= 1+ix
i
+
+i
....
2!
3!
4!
5!
ei x = 1 + i x +
x3 x5
x2 x4
+
... +i x
+
...
ei x = 1
2!
4!
3!
5!
e x
> 0,
(1 + e x )2
291
1.7 La fonction g (x) est clairement positive. Il faut dmontrer que lintgrale sur la
totalit du domaine de cette fonction est 1 :
Z
Z
f (x)
g (x) d x =
dx
x0
x 0 1 F (x 0 )
R
x f (x) d x
= 0
1 F (x 0 )
1 F (x 0 )
=
1 F (x 0 )
= 1.
1.8 On a
S(x) = Pr[X > x]
Z
dt
=
+1
x (t + )
=
.
x +
La figure E.2 prsente le graphique de cette fonction.
1.9 On a que Y = n X si, et seulement si, X = n Y . Ainsi,
Pr[Y = y] = Pr[X = n y]
!
n
=
p ny (1 p)n(ny)
ny
!
n
=
(1 p) y p ny , y = 0, 1, . . . ,
y
do Y Binomiale(n, 1 p).
1.10 a) On a Y = e X o X N (, 2 ). Par consquent,
F Y (x) = Pr[Y x]
= Pr[e X x]
= Pr[X ln x]
= F X (ln x)
et
f Y (x) = F Y0 (x)
1
= f X (ln x).
x
292
> library(actuar)
0.6
0.2
0.4
S(x)
0.8
1.0
1000
2000
3000
4000
5000
2 2
Var[Y ] = E [Y 2 ] E [Y ]2
= E [e 2X ] E [e X ]2
= M X (2) M X2 (1)
2
= e 2+2 e 2+
= e 2+ (e 1).
2
t /2
. On a
293
1.11 On a
|x| 1
dx
2
1 + x
Z
Z 0
x 1
x 1
dx +
dx
=
2
2
0 1+x
1 + x
Z
2
x
=
dx
0 1 + x2
Z
2 a x
dx
= lim
a 0 1 + x 2
= lim ln(1 + a 2 )
Z
E [|X |] =
= .
E [g (X )] =
g (x)
x=0
e x
x!
(x + 1)g (x)
e x+1
.
(x + 1)!
e x+1
g (x)
=
x!
x=0
X
x=0
x +1
x +1
X
e x
E g (X ) =
xg (x 1)
x!
x=1
X
x=0
xg (x 1)
e x
x!
= E [X g (X 1)].
1.13 Il suffit de remarquer que M + m = X + Y . Le rsultat dcoule ensuite directement par linarit de lesprance : E [M ] + E [m] = E [X ] + E [Y ].
294
y 3
= Pr X
4
y 3
= FX
4
= 1e
y3
4
7 7 (y3)
e 4
,
4
y > 3.
= Pr X 3 y
h
i
1
= Pr X y 3
Z 1
1 y3 2
=
x dx
9 0
y
= .
27
1
,
27
0 y 27.
1.16 Selon lnonc, X N (0,2 ) et Y = X 2 . Il faut voir que Y = X 2 nest pas une
transformation bijective ( une valeur de Y correspond plus dune valeur de
X ). On pose W = |X | et on trouve la densit de W laide de la technique de
la fonction de rpartition :
FW (w) = Pr[|X | w]
= Pr[w X w]
= F X (w) F X (w)
295
et donc
f W (w) = f X (w) + f X (w)
2
2
2
= p e x /(2 ) .
2
On pose maintenant Y = W 2 = |X |2 = X 2 et on trouve la densit de Y par la
technique du changement de variable :
d 1/2
f Y (y) = f W (y 1/2 )
y
dy
1
= f W (y 1/2 ) p
2 y
2
1
y/(22 )
= p e
p
2 y
2
=
(22 )1/2
( 21 )
y 1/2 e y/(2
p
puisque ( 12 ). On a donc que Y Gamma( 12 , 12 2 ). De manire quivalente, on peut aussi poser X = Z , o Z N (0, 1), et utiliser le rsultat connu
que Z 2 2 (1) Gamma( 12 , 12 ).
1.17 Si X est une variable alatoire dont la distribution est symtrique autour du
point a, alors E [X ] = a. On a donc
3 = E [(X a)3 ]
Z
(x a)3 f (x) d x
=
Z a
Z
=
(x a)3 f (x) d x +
(x a)3 f (x) d x.
y 3 f (y + a) d y +
3
y 3 f (y + a) d y
y f (y + a) d y +
Z
0
y 3 f (y + a) d y
= 0,
puisque f (y + a) = f (y + a) par symtrie autour du point a. Par consquent,
1 = 3 /3/2
2 = 0.
296
h
/5
/5
i
X 1 X 2
Pr | X 1 X 2 | <
= Pr
< p
< p
p
5
/ n/2 / n/2 / n/2
r
r
1 n
1 n
= Pr
<Z <
,
5 2
5 2
p
p
o Z N (0, 1). On doit donc trouver une valeur de n tel que Pr[Z n/(5p 2)]
p
0,005. On trouve dans une table de quantiles de la loi normale que n/(5 2) =
2,576, et donc que n 332.
1.22 a) On a X i Gamma(25, 12 ). Or, une somme de n lois gamma indpendantes
P
de paramtres i et est une loi gamma de paramtres ni=1 i et . Par
P
consquent, ni=1 X i Gamma(2 500, 12 ) et X Gamma(2 500, 50).
b) On obtient avec R
297
49 50 X 50 51 50
<
<
1
1
1
= 749 500
2(1 000)2
1 000 2
(2)(1)
2
= 500.
Lerreur quadratique moyenne est
= Var[]
+ b ()2
MSE()
= 750 + (500)2
= 250 750.
1.24 a) Par linarit de lesprance,
"
n
X
ai X i =
i =1
n
X
i =1
n
X
a i E [X i ]
ai
i =1
n
X
i =1
= .
ai
298
Var
n
X
ai X i =
n
X
i =1
i =1
= 2
a i2 Var[X i ]
n
X
i =1
Pn
2
i =1 a i
n
X
i =1
a i2
a i2 .
sous la contrainte
Pn
i =1 a i
= 1. Or,
1
1 2
=
ai
+
n
n
i =1
2
n
X
1
1
=
ai
+ ,
n
n
i =1
n
X
P
tant donn que le produit crois vaut 0. Ainsi, lexpression ni=1 a i2 est
minimise en choisissant a i = 1/n pour tout i . Par consquent,
X =
n 1
X
Xi
i =1 n
possde la plus petite variance parmi tous les estimateurs sans biais linaires.
1.25 On a
#
n
n
1X
1X
2
(X i ) =
E [(X i )2 ]
E
n i =1
n i =1
"
n
1X
2
n i =1
= 2 .
1.26 En utilisant la dfinition de lesprance, on obtient
1
E [T (X )] = 0 + 0 + (2)
= .
2
X
X
=n
1
.
n
n
299
On a
2
= E [X ] E [X ]
E []
n
np(1 p) + (np)2
= np
n
= np p(1 p) np 2
= np(1 p) p(1 p)
= p(1 p).
Par consquent, est un estimateur de avec un biais de p(1 p).
1.28 On sait que
Var[ X ] =
Var[X ]
= .
n
n
De plus,
X 2
ln f (X ; )
=E
1
= 2 Var[X ]
1
= .
= Var X .
n
300
1 1/1
x
(1 x)11 ,
n
X
1
1
ln x i n ln .
i =1
Par consquent,
d
`() =
d
Pn
i =1 ln x i
2
P
et = n 1 ni=1 ln x i .
c) On a
n Z 1 1
1X
E [] =
(ln x i )x i1/1 d x i
n i =1 0
Z 1
n
1X
1/
1
1/1
=
x ln x|0
xi
d xi
n i =1
0
n
1X
=
()
n i =1
= .
Chapitre 2
2.1 La franchise permet lassureur dconomiser au plus 250 $ par contrat. Lassureur conomise donc, pour les 12 contrats de son portefeuille,
250, 110, 250, 213, 98, 250, 250, 162, 131, 250, 250, 250,
pour un total de 2 464 $. Le montant total des sinistres sans la franchise est de
4 982 $. Le rapport dlimination de perte est donc
LER =
2 464
= 0,4946.
4 982
301
2.2 La limite permet lassureur dconomiser lexcdent de 100 000 $ par contrat.
Lassureur conomise donc, pour les huit contrats de son portefeuille,
0, 23 000, 323 000, 0, 113 000, 0, 0, 78 000,
pour un total de 537 000 $. Le montant total des sinistres sans la limite est de
1 146 000 $. Le rapport dlimination de perte est donc
LER =
537 000
= 0,4686.
1 146 000
E [Y ]
E [X ]
R 10
0 x f X (x) d x
50
0,87616
=
50
= 0,0175.
b) Avec une franchise forfaitaire, on a plutt
(
X , X 10
Y =
10, X > 10.
Le rapport dlimination de perte est
LER =
=
E [Y ]
E [X ]
R 10
R
0 x f X (x) d x + 10 10 f X (x) d x
E [X ; 10]
=
E [X ]
9,0634
=
50
= 0,1813.
E [X ]
302
(
0,
X 100
X 100,
X > 100.
E [Y ]
E [X ]
R
100 (x 100) f X (x) d x
E [X ]
0,138
=
40
= 0,0034.
Il est galement possible de rcrire la variable alatoire comme tant
(
Y =
X X,
X 100
X 100,
X > 100.
LER =
1
X
j =0
y j e y
j!
pour obtenir
LER =
40 39,862
= 0,0034.
40
303
2.5 Il est dit dans la question que lassureur limite ses paiements 200, la limite est donc de 270. En introduisant dabord la limite, lassureur conomise,
respectivement,
0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 90, 130,
pour un total de 220. En introduisant ensuite la franchise, lassureur conomise
en plus, respectivement,
20, 50, 70, 70, 70, 70, 70, 70, 70, 70, 70, 70,
pour un total de 770. Le montant total des sinistres sans la limite et la franchise
est de 1 745. Le rapport dlimination de perte est donc
LER =
770 + 220
= 0,567.
1 745
LERd =100 =
x f X (x) d x +
E [X ]
E [X ; 100]
E [X ]
E [X ; 100]
=
2 000
= 0,0465
X 100
X ,
Y =
100,
100 < X 30 000
X 100
X X + X ,
=
Ainsi,
LER =
X X + 100,
100 < X 30 000
304
e x
2
1
1
= (2e 2x ) + e x .
2
2
Lesprance limite est donc, en utilisant les formules pour lesprance
limite dune exponentielle,
1 1
1
2d
E [X ; d ] =
(1 e
)+
(1)(1 e d )
2 2
2
1 e 2d 1 e d
+
.
4
2
Y =
0,
X 0,25
ou encore
Y = max(X 0,25, 0)
= X min(X , 0,25)
(
X X,
X 0,25
=
X 0,25, X > 0,25.
partir de cette reprsentation, il est facile de voir que
E [Y ] = E [X ] E [X ; 0,25]
1 1
3 e 0,5 2e 0,25
=
+
4 2
4
= 0,541.
Lesprance de la svrit est de un sinistre tous les dix ans, donc de 0,1.
Ainsi, la prime pure est
= (0,541)(0,1) = 0,0541.
305
(
0,
X 50 000,
X 50 000
X > 50 000,
ou encore
Y = max(X 50 000, 0)
= X min(X , 50 000)
(
X X,
X 50 000
=
X 50 000, X > 50 000
partir de cette reprsentation, il est facile de voir que
E [Y ] = E [X ] E [X ; 50 000]
= 1 091,09.
b) Soit Y la variable alatoire du montant conomis par le rassureur. On
dfinit
(
X,
X 100 000
Y =
100 000, X > 100 000.
306
1
2 500
=
1,5 1
= 5 000.
E [X ] =
4 219,13
= 0,8438.
5 000
f X (x + d )
,
1 F X (d )
x > 0.
On a donc
E [Y P ] =
=
=
=
=
Z
1
x f X (x + d ) d x
1 F X (d ) 0
Z
1
(y d ) f X (y) d y
1 F X (d ) d
Z
1
y f X (y) d y d (1 F (d ))
1 F X (d ) d
Z
Z d
1
y f X (y) d y
y f X (y) d y d (1 F (d ))
1 F X (d ) 0
0
E [X ] E [X ; d ]
1 F X (d )
307
0.00
0.00
0.04
0.10
0.08
0.20
10
20
30
40
50
10
20
30
40
50
10
20
30
40
50
(c) Coassurance de 80 %
franchise = TRUE)
lwd = 3)
308
E [X 1995
] = E [X 1995 ] E [X 1995 ; 500]
= 3 939,3 421,3
= 3 518.
Enfin, le rapport dlimination de perte est
LER =
3 939,3 3 518
= 0,1069.
3 939,3
= 0,4107.
c) La nouvelle variable alatoire est
0,
Y =
X 1995 500
4 000 500, X
1995 > 4 000,
ou encore,
Y = max(min(X 1995 , 4 000) 500, 0)
= min(X 1995 , 4 000) min(X 1995 , 500)
X 1995 500,
4 000 500,
do
E [Y ] = E [X 1995 ; 4 000] E [X 1995 ; 500]
= 1 255,23.
309
+2 /2
ln(u) 2
ln(u)
+u 1
f Y P (x) =
do Y P Pareto(, + ).
310
Pr[X > x + d ]
Pr[X > d ]
(x+d )
e
e d
=e
,
x
E A [Y P ] =
311
X X , X < d
= X d,
u d ,
d X <u
X u.
X d
0,
(Y S )2 =
2 (X 2 2d X + d 2 ), d < X < u
2 (u 2 2ud + d 2 ), X u
2
2
2
(X X 2d X + 2d X ), X d
2 (X 2 d 2 2d X + 2d d ),
2 (u 2 d 2 2d u + 2d d ),
d <X <u
X u.
On a alors
E [(Y S )2 ] = 2 (E [X 2 ; u 2 ] E [X 2 ; d 2 ] 2d E [X ; u] + 2d E [X ; d ]).
La variance est donc
Var[Y S ] = E [(Y S )2 ] E [Y S ]2
= 2 (E [X 2 ; u 2 ] E [X 2 ; d 2 ] 2d E [X ; u] + 2d E [X ; d ])
2 (E [X ; u] E [X ; d ])2 .
c) Suite une inflation de 100r %, la dfinition de la variable alatoire Y S
quivalente celle utilise en a) est
d
u
min X ,
Y S = (1 + r ) min X ,
1+r
1+r
X < d /(1 + r )
X X ,
= (1 + r ) X d /(1 + r ),
d /(1 + r ) X < u/(1 + r )
u i
d
E [Y S ] = (1 + r ) E X ;
E X;
.
1+r
1+r
312
X <d
0,
YS=
X d, d X < u
u d , X u.
Par consquent,
E [Y S ] = (0)Pr[X < d ] +
ud
ud
d
u
0
(x d ) f X (x) d x
(x d ) f X (x) d x
+ (u d )(1 F X (u))
Z u
Z
=
x f X (x) d x d F X (u)
0
d
0
x f X (x) d x + d F X (d )
+ (u d )(1 F X (u))
Z
Z u
x f X (x) d x + u(1 F X (u))
=
0
d
0
x f X (x) d x d (1 F X (d ))
= E [X ; u] E [X ; d ].
2.19 Lorsquil y a bonus, son montant est
1
450 000 S
B = max 0,
= 150 000 min(S, 450 000),
3
3
do
1
E [S; 450 000]
3
1
= 150 000 (220 321,36)
3
= 76 559,55.
E [B ] = 150 000
313
d
(X d ), d < X d
P
Y = d d
X,
X > d .
b) Pour pouvoir valuer lesprance, il est plus facile de rcrire la variable
sous la forme
d
d
P
Y =X+
min(X , d )
min(X , d ) X > d
d d
d d
d
d
X +
X
d, d < X d
d d
d d
=
d
d
d
d , X > d .
X +
d d
d d
Par la dfinition de lesprance limite ou en utilisant le rsultat de lexercice 2.9,
on obtient directement
E [Y P ] =
E [X ] + d E [X ; d ]/(d d ) d E [X ; d ]/(d d )
.
1 F X (d )
Chapitre 3
3.1 a) On peut calculer puis tracer la fonction de rpartition empirique aisment
avec la fonction ecdf de R ; voir la figure E.4. Quant la fonction de masse
de probabilit empirique, la faon la plus simple de la calculer est partir
de la fonction table ; voir la figure E.5.
b) Il faut dabord dterminer le nombre de donnes dans chacune des classes.
On a n 1 = 4, n 2 = 10, n 3 = 2 et n 4 = 4. Lquation de logive est alors
0,
x 2
(x 2)/25,
2<x 7
(x 5)/10,
7 < x 12
F20 (x) =
(x + 58)/100, 12 < x 22
(x + 42)/80, 22 < x 38
1,
x > 38
Les fonctions grouped.data et ogive de actuar permettent, dans lordre,
de dfinir un objet de donnes groupes et de calculer son ogive ; voir la
figure E.6.
314
1.0
ecdf(x)
0.8
0.6
Fn(x)
0.4
0.2
0.0
10
20
30
40
315
> table(x)
x
3
9 10 11 16 21 23 26 29 36
fn
10
15
20
25
30
35
unique(x)
316
1.0
ogive(xg)
0.8
0.0
0.2
0.4
F(x)
0.6
10
15
20
25
30
35
317
> hist(x)
4
0
Frequency
Histogram of x
10
20
30
40
0,
1/25,
1/10,
f20 (x) =
1/100,
1/80,
0,
x 2
2<x 7
7 < x 12
12 < x 22
22 < x 38
x > 38.
Le package actuar dfinit une mthode de la fonction hist pour les donnes groupes ; voir la figure E.7.
318
310 + 69 379 + y
,
= (0,5)
+
500
500
do lon trouve que y = 81.
3.4 Les donnes sont entres dans R avec
> (x <- grouped.data(Group = 1000 * c(0, 1, 3, 5, 10, 25, 50, 100, Inf),
+
Frequency = c(16, 22, 25, 18, 10, 5, 3, 1)))
Group Frequency
1
(0,
1000]
16
(1000,
3000]
22
(3000,
5000]
25
(5000, 10000]
18
(10000, 25000]
10
(25000, 50000]
7 (50000, 100000]
8 (100000,
Inf]
Pour calculer logive de ces donnes, la borne infinie de la dernire classe doit
tre remplace par une valeur trs grande par rapport aux autres bornes. Il
ne faut pas que cette valeur soit trop grande si on veut avoir un graphique
intressant. Il ne faut pas supprimer la dernire classe. La figure E.8 prsente
les ogives avec 200 000 et 2 000 000 comme dernire borne. On cherche
Pr[2 000 X 6 000] = F 100 (6 000) F 100 (2 000).
Or,
319
> plot(Gn)
> plot(Gn)
1.0
0.0
0.4
0.2
0.4
F(x)
50000
150000
x
0.8
0.2
0.6
0.0
F(x)
0.8
ogive(x)
0.6
1.0
ogive(x)
500000
1500000
x
F IG . E.8: Ogive des donnes de lexercice 3.4 avec diffrentes dernires bornes
Density
0.010
0.000
Density
100
200
300
0.020
320
400
100
N = 1 Bandwidth = 20
150
200
250
300
350
N = 4 Bandwidth = 20
5
X
f 5 (y j )t j (150)
j =1
300 150
= 225.
2
b) La figure E.9(a) prsente les quatre noyaux (quatre densits) sur le mme
graphique et la figure E.9(b) prsente leur somme pondre, cest--dire
f(x).
3.7 La figure E.10 prsente la distribution empirique des donnes.
a) On voit que pour une largeur de bande de 0,5, aucune donne ne va contribuer la densit au point 6,2.
321
fn
10
unique(x)
322
1
6,2 7 + 1
f (6,2) =
= 0,02.
10
(1)2
g) Pour une largeur de bande de 2, il y a trois valeurs, 5, 7 et 8 qui vont contribuer la densit au point 6,2 :
1
6,2 7 + 2
1
6,2 8 + 2
f(6,2) =
+
10
(2)2
10
(2)2
3
6,2 + 5 + 2
+
10
(2)2
= 0,095.
3.8 On utilise lquation dun estimateur avec noyaux triangulaires :
1/5
1/5
f (5) =
(5 (6 a)) +
(5 (4 + a))
a2
a2
= 0,01.
En simplifiant, on trouve 0,05a 2 2a + 2 = 0 et, en choisissant la bonne racine,
a = 1,0263.
3.9 On entre les donnes individuelles de lexercice dans R avec
> x <- c(2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 2, 3, 3, 3, 3,
+
4, 4, 4, 5, 5, 5, 5, 6, 6, 6, 6, 8, 8, 9, 15,
323
1.0
ecdf(x)
0.8
0.6
Fn(x)
0.4
0.0
0.2
10
20
30
40
324
0.15
0.00
0.05
0.10
Density
0.20
0.25
0.30
Histogram of xg
10
20
30
40
50
xg
une classe. Pour les deux bornes extrmes, 0 aurait peut-tre t un peu
plus intuitif comme choix de borne infrieure que 1,5. La borne suprieure
est logique, car suprieure la valeur maximale de lchantillon, mais est
totalement arbitraire sinon (on aurait pu choisir, par exemple, 50).
c) La figure E.12 prsente lhistogramme des donnes cr partir de lobjet
de donnes groupes.
d) On a simplement
> mean(x)
[1] 9.225
> sd(x)
[1] 10.23691
325
3.10 Tout dabord, on a F 10 (4) = 0,40, E [F 10 (4)] = F (4), et Var[F 10 (4)] = F (4)(1
F (4))/10. En utilisant le Thorme central limite, on peut poser que
F 10 (4) E [F 10 (4)]
Pr 1,96 p
1,96 0,95,
Var[F 10 (4)]
soit
#
p
10(F 10 (4) F (4))
1,96 0,95.
Pr 1,96 p
F (4)(1 F (4))
"
p
p
F 10 (4)(1 F 10 (4)) = 0,24 puis en isolant
F
(4)(1
F
(4))
10
10
r
0,24
0,4 (1,96)
10
s
(0,0964, 0,7036).
3.11 En utilisant lquation de lestimateur de Nelsonalen, on obtient
3 10 11
+
+
52 40 19
= 0,8866.
Hn (1 200) =
326
yi
Hn
1 000
3 400
4 500
7 500
15 000
17 500
1,5000
1,1000
0,9167
0,4167
0,8333
1,2500
0,0
0,0
0,5
0,5
0,5
0,0
0,0500
0,0526
0,0556
0,0588
0,0625
0,0667
TAB. E.1: Rsultats intermdiaires du calcul de lestimation par noyaux pour les
donnes de lexercice 3.12
3.12 a) En utilisant lquation de lestimateur de Nelsonalen, on trouve
1
20
= 0,0500
1
1
+
Hn (3 400) =
20 19
= 0,1026
1
1
1
Hn (4 500) =
+
+
20 19 18
= 0,1582
1
1
1
1
Hn (7 500) =
+
+
+
20 19 18 17
= 0,2170
1
1
1
1
1
Hn (15 000) =
+
+
+
+
20 19 18 17 16
= 0,2795
1
1
1
1
1
1
Hn (17 500) =
+
+
+
+
+
20 19 18 17 16 15
= 0,3462.
Hn (1 000) =
h(10
000) =
(0,5)(0,0556 + 0,0588 + 0,0625)
6 000
= 0,00001449.
3.13 On a
1 =
3
.
3
327
On a donc
1 =
3 200 0003/2
= 0,559.
On a donc
4
4
2 =
328
E [min(X , 320)] =
50 + 0
20
100 + 50
34
+
2
100
2
100
22
320 + 200 9,6
200 + 100
+
+
2
100
2
100
24 9,6
+ 320
100
329
3.17 tant donn que lintervalle est petit, on peut en calculer le niveau de confiance
exactement. En utilisant la loi binomiale avec paramtres n = 5 et p = 0,5, on
obtient
!
3
X
5
Pr[X (2) 0,5 < X (4) ] =
(0,5)k (0,5)5k
k
k=2
= 0,625.
3.18 tant donn que lintervalle est grand, on va utiliser lapproximation normale
avec correction pour la continuit pour dterminer le niveau de confiance.
On a, avec Y N (250, 125),
Pr[240 0,50 < 260] Pr[239,5 K < 259,5]
= (0,85) (0,94)
= 0,6287.
3.19 La valeur 10 est la 10e statistique dordre et la valeur 20 est la 14e statistique
dordre. Comme lintervalle est petit, on peut en calculer le degr de confiance
exactement. Soit N Binomiale(20, 0,55),
Pr[X (10) 0,55 < X (14) ] = Pr[N = 10, 11, 12, 13]
!
13 20
X
=
(0,55)k (0,45)20k
k
k=10
= 0,1593 + 0,1771 + 0,1623 + 0,1221
= 0,6208.
3.20 a) On obtient aisment
1
f Y (ln(x))
x
=
(ln(x))1 e ln(x) ,
()x
f X (x) =
x > 1.
On remarque que comme Y est dfinie sur [0, ), X = e Y est dfinie sur
[1, ). Cette distribution est la log-gamma de paramtres et .
b) La fonction R de la figure E.13 calcule le biais empirique pour des valeurs
de , n et r donnes. On remarquera que cette fonction dfinit une fonction interne qui se charge des tapes 2 et 3 de lalgorithme prsent dans
lexpos de lexercice. Cette fonction est ensuite passe replicate pour
raliser efficacement ltape 4 de lalgorithme.
330
331
332
0.0
0.2
0.4
Fn(x)
0.6
0.8
1.0
ecdf(x)
4.8
5.0
5.2
5.4
et
> 0.30 * xs[70] + 0.70 * xs[71]
[1] 5.068952
333
> plot(Gn)
ogive(xg)
1.0
3.0
Histogram of x
0.8
0.6
0.4
F(x)
2.0
0.2
1.0
Density
0.0
0.0
4.6
4.8
5.0
5.2
5.4
4.6
4.8
5.0
5.2
5.4
(a) Histogramme
(b) Ogive
45%
70%
4.972414 5.068000
Chapitre 4
4.1 En utilisant la technique de la fonction de rpartition, on a
F Y (y) = Pr[Y y]
= Pr[c X y]
= Pr[X y/c]
= F X (y/c)
= 1
+ y/c
c
= 1
c + y
et donc, Y Pareto(, c).
334
F X (ln(y)) =
1 x
e dx
2
1
= e ln(y)/ .
2
Pour 0 < x < , on a 1 < y < , et donc
0
1 x/
e dx +
2
1 ln(y)
= 1 e .
2
Z
F X (ln(y)) =
ln(y)
Z
0
1 x/
dx
e
2
Par consquent,
(
F (y) =
1 ln(y)/
,
2e
1 ln(y)/
,
1 2 e
0<y <1
y 1.
x Z x
1
1
t /
2
t /
=
t
e
( 1)t
e
dt
+
0
( 1)!
0
Z x
x 1 e x/
1
=
+
t 2 e t / d t
( 1)!1 ( 2)!1 0
Z x
1
= Pr[Y = 1] +
t 2 e t / d t
( 2)!1 0
Pr[X x] =
335
y
F Y (y) = F X
1 y
y/(1 y)
= , ;
y/(1 y) +
= (, ; y),
o (a, b; x) est la fonction de rpartition dune distribution Bta(a, b) value
au point x. On a donc que Y Bta(, ).
4.5 En utilisant la technique de la fonction de rpartition, on obtient
F Y (y) = Pr[Y y]
= Pr[5X 1/4 y]
y 4
= Pr X >
5
y 4
= 1 FX
5
1
=
1 + (5/y)4
qui est la fonction de rpartition dune variable alatoire avec distribution Burr
inverse de paramtres = , = 4 et = 5.
4.6 a) Par dfinition, la distribution de Y est nomme log-gamma. On remarque
que comme X est dfinie sur [0, ), Y = e X est dfinie sur [1, ). On a donc
1
f X (ln(y))
y
=
(ln(y))1 e ln(y) ,
()y
f Y (y) =
y 1.
=
,
1
>1
E [Y 2 ] = E [e 2X ]
= M X (2)
=
,
2
> 2,
336
2
=
,
2
1
> 2.
=
,
k
> k.
(1 + i )
= 1
(1 + i ) + y
F Y (y) = F X
f Y (y) =
()y +1
= 1
+ y
et donc, Y Burr(, , 1/ ).
337
= 1
+ (e y 1)
= 1 e y ,
y 0.
Ainsi, Y Exponentielle().
4.10 Soit Y , la variable alatoire du montant des sinistres en 2007. On dfinit
Y = (1,04)(1,045)(1,16)X
= 1,260688X .
En se reportant lexercice 4.7 b), on a que Y Burr( = 0,5, = 2, = 3,7821)
et donc que
Pr[Y > 4] = 1 F Y (4) = 0,6870.
4.11 a) On observe que la variable alatoire X obit une distribution Pareto
translate(3, 1). En utilisant la technique de la fonction de rpartition, on
trouve
y
F Y (y) = F X
1,10
1,10 3
= 1
.
y
!
10 x
Pr[X = x] =
(1 )10x d
x
0
!Z
1
10
=
x (1 )10x d
x 0
338
!
10 (x + 1)(11 x)
Pr[X = x] =
x
(12)
10!
x!(10 x)!
(10 x)!x!
11!
1
= .
11
=
Par consquent,
Pr[X > 6] =
10
X
Pr[X = i ]
i =7
4
.
11
d
f X (x) =
()
()
0
Z
x 1
+1 e (x+) d
=
()() 0
qui devient, en reconnaissant sous lintgrale la forme fonctionnelle dune
distribution Gamma( + , x + )
f X (x) =
x 1 ( + )
()()(x + )+
Pr[X = x] =
(1 )x1
( + )
=
()()
Z
0
( + ) 1
(1 )1 d
()()
(1 )x+2 d
339
f X (x) =
x 1 e x
x 1 e x
d
()
Z
x 1 (x +)
=
e
d
()
0
x 1
,
=
( + x )+1
E [X ] = E [E [X |]]
(4)(2)
=E
(5)
1
= E []
4
5
=
4
et
(3)(3)
(4)(2) 2
(4)(2)
=E
+ Var
(5)
(5)
(5)
5
1
=
E [2 ] +
Var[]
48
16
145
=
.
48
4.17 Pour commencer, on utilise le lien entre le taux dchec et la fonction de
340
Rx
0
(x|) d t
Z x
3
= exp
dt
+t
0
= exp 3 ln
+x
3
= exp ln
+x
=
.
+x
F (x|) = 1
+x
=E
2
= 500
et
Var[X ] = E [Var[X |]] + Var[E [X |]]
3
+ Var
=E
4
2
= 850 000.
4.18 Soit f () la fonction de densit de probabilit dune Log-normale(, 2 ) et g ()
celle dune Gamma(, ). Pour comparer les queues de ces deux distributions,
il faut valuer
f (x)
lim
.
x g (x)
En liminant les termes qui ne dpendent pas de x, on obtient
lim
x 1 e (ln(x))
x 1 e x/
/22
= lim e (ln(x))
x
/22 ln(x)+x/
Or, de lexercice 1.4 on sait que x tend plus rapidement vers linfini que ln(x).
Lexposant tend donc vers , do la distribution log-normale a une queue
plus paisse que la distribution gamma.
341
4.19 Une fonction desprance de vie rsiduelle linaire en x indique une distribution de Pareto telle que
1
x+
.
e(x) =
1
1
partir de e(x) = 2 000 + 2x, on trouve que = 1,5 et = 1 000. En utilisant les
formules de lannexe A pour lesprance limite dune loi de Pareto, on a que
le LER est
E [X ] E [X ; x]
E [X ]
= 0,30115.
LER =
0,5e x
,
0<x 2
(1; 2)
25 251 5x
0,2
5 x
e
f X (x) =
,
2<x 8
(25; 40) (25; 10)
(25)
0,3
12144 x 1441 e 12x
, 8 < x 16.
(144; 192) (144; 96)
(144)
342
f X (x) =
p/10,
0 < x < 10
3
(3)(100 )
1
, x 10
(1 p)
4
(100 + x) (100/110)3
0 < x < 10
p/10,
= (3)(1103 )
(1 p), x 10.
(100 + x)4
p
,
10
soit
p
3
(1 p) = .
110
10
En rsolvant pour p, on trouve p = 3/14.
4.23 a) On pose X Weibull(1 , 1 ) et Y = X 1 Weibull inverse(2 , 2 ). On sait
de lannexe A que tous les moments positifs de la distribution Weibull
existent, alors que ceux de la distribution Weibull inverse nexistent que
pour k < 2 . Par ce critre, on voit que la distribution Weibull Inverse
possde une queue plus lourde.
Dautre part, on a
2
1 ( x)
f Y (x) 2 2 2 x 2 e 2
=
1 1 1 (1 x)1
f X (x)
1 1 x
e
x 1 2 e (2 x)
+(1 x)1
do
f Y (x)
(1 x)1 (2 x)2 (1 + 2 ) ln(x).
ln
f X (x)
Lorsque x , le terme central tend vers 0. Comme x tend plus rapidement vers que ln(x), on a que
f Y (x)
lim ln
= lim (1 x)1 (2 x)2 (1 + 2 ) ln(x)
x
x
f X (x)
= ,
343
0.000
0.005
0.010
0.015
0.020
16
18
20
22
24
F IG . E.16: Comparaison des queues des distributions Weibull (trait mince) et Weibull inverse (trait pais)
do
f Y (x)
S Y (x)
= lim
= .
x f X (x)
x S X (x)
lim
344
S X (y)
et y
dy
E [X ]
0
Z t y
ty
S X (y)
f X (y)
e
e
dy
=
+
t
E [X ] 0
t
E [X ]
0
1
M X (t )
=
+
t E [X ] t E [X ]
M X (t ) 1
=
.
t E [X ]
Z
M Y (t ) =
Ce rsultat suppose que lim y e t y S X (y) = 0. En appliquant la rgle de lHpital, on voit quil sagit de la mme limite que t 1 lim y e t y f X (y) qui doit
tre 0 sinon lintgrale dfinissant M X (t ) ne convergerait pas.
4.25 a) Par dfinition de la fonction de survie :
Z
S(x) =
(1 + 2t 2 )e 2t d t
x
= (1 + x + x 2 )e 2x ,
x 0.
h(x) =
c) On a dabord
Z
x
S(t ) d t =
(1 + t + t 2 )e 2t d t
= (1 + x + 0,5x 2 )e 2x
et donc
e(x) =
S(x) d x
S(x)
1 + x + 0,5x 2
.
1 + x + x2
d) On a
= 2.
1 + 2x
1 + x + x2
345
e) On a
lim e(x) =
limx h(x)
1
= .
2
Chapitre 5
P
P
5.1 On a 5i =1 x i = 6 211 et 5i =1 x i2 = 26 040 101. Pour trouver les estimateurs des
moments de et , on pose
E [X ] = E [E [X |]] = E [1 ] =
6 211
=
1
5
et
E [X 2 ] = E [E [X |]] = E [2 ] =
22
26 040 101
=
.
( 1)( 2)
5
P
P
5.3 On a 5i =1 x i = 5 850 et 5i =1 x i2 = 5 867 500. Pour trouver les estimateurs des moments de et il suffit de poser gaux les deux premiers moments empiriques
et thoriques :
5 850
=
6
346
5 867 500
.
6
F (x) =
E [X ] =
xpx p1 d x
p
.
p +1
p =
P
P
5.6 On a 5i =1 x i = 10 000 et 5i =1 x i2 = 30 000 000. On galise les deux premiers
moments thoriques et empiriques :
e +
2
2
et
2
e 2+2 =
10 000
5
30 000 000
.
5
347
p
2
= x = 4,2,
2
do = 3,3511.
5.8 a) Posons = . On a
= 0,25
et
1 e ( 0,80 ) = 0,80.
On trouve 0,80 = ( ln 0,20)1/ / = 1 947.
5.9 La distribution marginale de la variable X est une loi de Pareto(, ). Ainsi,
pour estimer les paramtres et par la mthode des quantiles, on pose
S X (450) =
+ 450
= 0,001
et
S X (50) =
+ 50
= 0,125.
ln +50
ln(0,125)
ln(0,001) ln
+450
0,3010
ln
= ln
+ 450
+ 50
0,3010 ( + 50) = ( + 450)0,3010 .
348
E [X ; d ]
E [X ]
E [X ; d ] d (1 F (d ))
.
E [X ]
x a
,
ba
on a
50 a
= 0,80
ba
55 a
= 0,90,
ba
do on obtient a = 10 et b = 60.
5.12 Pour X Bernoulli(p), on a la fonction de vraisemblance
L(p; x 1 , . . . , x n ) =
n
Y
p xi (1 p)1xi
i =1
Pn
=p
i =1 x i
Pn
(1 p)n
i =1 x i
349
et la fonction de log-vraisemblance
l (p; x 1 , . . . , x n ) =
n
X
x i ln(p) + n
n
X
x i ln(1 p),
i =1
i =1
do
Pn
0
l (p; x 1 , . . . , x n ) =
i =1 x i
Pn
i =1 x i
p
1p
n x n n x
=
.
p
1p
On trouve donc p = X .
5.13 a) On a la fonction de log-vraisemblance
n (x )2
1X
i
l (, ) = n ln p
2 i =1
2
2
i
l (, 2 ) =
2
i =1
Pn
2
n
i =1 (x i )
2
)
=
l
(,
+
.
2
22
24
(X i X )2 = S 2 .
n i =1
Pn
2
2
n
i =1 (x i )
2
l
(,
)
=
4
4
2
6
Pn
2
xi
l (, 2 ) = i =1 4
.
2
350
Pn
i =1 (X i )
=0
E
4
et
"
n
E
24
Pn
i =1 (X i )
6
n
24
/n
0
.
=
0
24 /n
Or, on sait que la distribution asymptotique conjointe des estimateurs du
maximum de vraisemblance est une normale multivarie sans biais et de
matrice variance-covariance .
c) On rappelle que () et () sont, dans lordre, les fonctions de rpartition
et de densit de probabilit dune loi N (0, 1). Par consquent,
A=
1 c
h(, 2 ) =
B=
1 c c
2
h(,
)
=
,
2
2 3
et
do
2 /n
0
A
)] = A B
Var[h(,
4
0
2 /n B
c 2 1 (c )2
+
.
=
n
2n2
2
351
!
n
Pn
Y
2
n n
L() = 2
x i e i =1 xi
i =1
l () = n ln(2) + n ln() +
n
X
ln x i
n
X
i =1
i =1
x i2
n
n X
l 0 () =
x 2.
i =1 i
P
On trouve alors = n/ ni=1 x i2 . En calculant la drive seconde de la fonction
de log-vraisemblance, on voit quil sagit bien dun maximum.
5.15 a) On a f (x) = px p1 , do
L(p) = p n
n
Y
i =1
p1
xi
l (p) = n ln(p) + (p 1)
n
X
ln(x i )
i =1
l 0 (p) =
n
n X
+
ln(x i ).
p i =1
P
On trouve alors p = n/ ni=1 ln(x i ). En calculant la drive seconde, on
voit quil sagit bien dun maximum.
n
p2
do
I (p) = nE [p 2 ] =
n
p2
et
=
Var[p]
1
p2
=
.
I (p)
n
c) On sait que
p p 1,96
De a) et b), on a donc
Var[p].
p
p p 1,96 p .
n
352
E [X ] =
=
xpx p1 d x
p
.
p +1
p
,
1+p
h 0 (p) =
1
.
(1 + p)2
do
p
1
=
1+p
n
et
d E [X ] =
Var
et donc
1
1 + p
p 2
n
p
p
1,96
E [X ]
p .
1 + p
2 n
(1 + p)
5.16 On a
L() = 5 5
5
Y
( + x i )1
i =1
l () = 5 ln() + 5 ln ( + 1)
5
X
i =1
l 0 () =
5
X
5
+ 5 ln()
ln( + x i ).
i =1
ln( + x i )
353
P
On obtient alors = 5/( 5i =1 ln(+x i )5 ln ) = 3,8629. En calculant la drive
seconde de la fonction de log-vraisemblance, on vrifie quil sagit bien dun
maximum.
X2
1
I () = nE
2 2 3
2n
= 2
)
d ]
MSE(
2 2
n
= 0,20.
=
2
E [X ] =
.
1
354
X
= p
.
X 1
b) On a
Q
2n ni=1 ln(x i )
L() = Qn
+1
i =1 x i
n
n
X
X
l () = 2n ln() +
ln(ln(x i )) ( + 1) ln(x i )
i =1
i =1
n
2n X
l 0 () =
ln(x i ).
i =1
Pn
i =1 ln(x i ).
5.20 a) On a
L() =
P5
Q
(5 )5 ( 5i =1 x i4 )e i =1 xi
((5))5
5
5
X
X
x i 5 ln (5)
l () = 25 ln() + 4 ln x i
i =1
i =1
5
25 X
l 0 () =
xi ,
i =1
P
do = 25/ 5i =1 x i = 1/2.
b) On a
l 00 () =
25
2
25
I () = E 2
25
=
(5/8)2
= 64.
=
Par consquent, Var[]
1
64 .
355
2
1
1
1
1
= 1
.
3
6
12
12
tant donn quil faudra faire appel des mthodes numriques pour rsoudre ce problme, on peut tout aussi bien minimiser la fonction de vraisemblance au lieu de la fonction de log-vraisemblance.
5.22 On a
L() = Pr[0 X 1]Pr[X 2]
= (1 e )e 2
l () = ln(1 e ) 2
l 0 () =
e
1 e
2.
On obtient = ln(1,5).
5.23 a) Par la mthode du maximum de vraisemblance habituelle, on trouve
n
= Pn
2
i =1 x i
P k
Var[]
= (k 2 e k )2 Var[].
n
X
i =1
n
n X
l 0 () =
x2
i =1 i
n
l 00 () = 2
x i2 + . . .
356
h ni
n
E 2 = 2
et
=
Var[]
2
.
n
Par consquent,
2
Var[Pk ] =
k 4 2 e 2k
.
n
= 0,0405.
3
n
X
ln( + x i )
i =1
et
n
2
l (, ; x) = 2
2
2
n
X
2
n
1
l
(,
;
x)
=
+
(
+
1)
2
2
i =1 + x i
2
n
n X
1
l (, ; x) =
.
i =1 + x i
Pour la suite, on aura besoin des rsultats intermdiaires
1
dx
2
+1
0 ( + x) ( + x)
= 2
( + 2)
Z
1
1
E
=
dx
+1
+ X
0 + x (x + )
1
=
.
+1
1
+ X
357
Ainsi,
h n i
2
E
l
(,
;
X
)
=
E
2
2
n
= 2
"
2 #
2
n
X
1
n
l (, ; X ) = E 2 + ( + 1)
E
2
i =1 + X i
n
= 2
( + 2)
"
2
#
n
1
n X
l (, ; X ) = E
i =1 + X i
n
=
( + 1)
I (, ) =
n
( + 1)
n
( + 1)
.
n
2 ( + 2)
2 ( + 1)2
=
( + 1)( + 2)
n
( + 1)( + 2)
.
2
2
( + 1) ( + 2)
De l, on obtient
2 ( + 1)2
= 0,28125
50
2
2
= ( + 1) ( + 2) = 656 250
d ]
Var[
50
+ 1)( + 2)
= (
d ,
)
Cov(
= 393,75.
50
d ]
=
Var[
358
=
+ 10
h(, )
ln
=
+ 10
+ 10
h(, )
10
=
.
+ 10
( + 10)2
Or,
= 0,0816
)
h(,
h(, )
= 0,1023
(,
h(, )
= 0,0292
(,
)
et donc
24 10 0,1023
d
)] = 0,1023 0,0292
Var[h(,
0,0292
10 40
= 0,2254.
p
Lintervalle de confiance est donc 0,0816 (1,44) 0,2254. tant donn quil
sagit dun intervalle pour une probabilit, la borne infrieure ne peut tre
plus petite que 0 (et la borne suprieure ne peut tre plus grande que 1).
Lintervalle de confiance est donc (0, 0,7653).
5.26 On sait que lestimateur du maximum de vraisemblance de est = X 1 et il
= 2 /n. Ici, on a = 0,0187.
est simple dtablir que Var[]
a) On a
h() = E [X ; 50]
1 e 50
d h() e 50 (50 + 1) 1
=
,
d
2
=
do
50
(50 + 1) 1
= e
Var[h()]
2
!2
2
n
359
et
= (686,57)2 (0,000 043 88)
d )]
Var[h(
= 20,68.
b) On procde comme en a) avec
h() = 0,95
ln(0,05)
d h()
ln(0,05)
.
=
d
2
=
= 3 196.
d )]
On obtient Var[h(
5.27 a) On a que X |A = obit une loi de Bernoulli() et que A obit une loi
U (0, 1). On cherche
P3
f (|x 1 , x 2 , x 3 )
i =1 x i
P3
(1 )3
i =1 x i
(1)
= (1 ) .
On reconnat ici la forme fonctionnelle dune distribution Bta(2, 3). On
sait que, si la fonction de perte choisie est lerreur quadratique, lestimateur bayesien est lesprance de la distribution a posteriori. On a donc
=
2
2
= = 0,4.
2+3 5
b) On a
Z
0,4
0,2
(5)
(1 2 + 2 ) d
(2)(3)
= 0,3432.
5.28 On a que X | = Poisson() et que Gamma(, ). On a donc
!
!
Pn
e n i =1 xi e 1
f (|x 1 , . . . , x n )
Qn
()
i =1 x i !
Pn
= e (+n) +
i =1 x i 1
360
n
3
f (|x 1 , . . . , x n ) 3e
Qn
+1
i =1 (1 + x i )
n e 3
+1
i =1 (1 + x i )
e 3
n
= Qn
i =1 (1 + x i )
= Qn
= n e
P
avec = 3 + ni=1 ln(1 + x i ). On reconnat alors la forme fonctionnelle
dune loi Gamma. On a donc, comme densit a posteriori, une loi Gamma
de paramtres = n + 1 et .
b) On sait que, si la fonction de perte choisie est lerreur quadratique, lestimateur bayesien est lesprance de la distribution a posteriori. On a
donc
n +1
7
=
=
= 0,68.
Pn
3 + i =1 ln(1 + x i ) 3 + 7,27
5.30 a) On a que X |B = obit une loi Exponentielle() et que B obit une loi
Gamma(2, 3). On a
P5
f (|x 1 , . . . , x 5 ) 6 e (3+
i =1 x i )
P
Puisque 5i =1 x i = 47, on reconnat ici la forme fonctionnelle dune loi
Gamma(7, 50). Avec une fonction de perte quadratique, lestimateur bayesien est lesprance de la distribution a posteriori. On a donc
7
=
= 0,14.
50
b) Avec une fonction de perte valeur absolue lestimateur bayesien est la
mdiane de la distribution a posteriori. Il faut donc choisir tel que
= x] = (7; 50)
= 1.
Pr[B |X
2
361
2 2
2 (1 )
f (|x 1 = 2, x 2 = 2) = R 0,75
3 2
2 2 d
0,25 2 2 (1 )
4 (1 )2
= R 0,75
4
2
0,25 (1 ) d
= 141,22 4 (1 )2 .
Avec une fonction de perte quadratique, lestimateur bayesien est lesprance de la distribution a posteriori. Ainsi,
Z 0,75
= 141,22
5 (1 )2 d = 0,5668.
0,25
b) On a
Pr[0,6 < < 0,7|X 1 = 2, X 2 = 2] = 1441,22
0,7
0,6
4 (1 )2 d p
= 0,3055.
5.32 a) Soit X la variable alatoire du montant dun sinistre en millions. On a
Y = X |X > 1,5. Par consquent,
Pr[Y > 29,5] =
F X (29,5) F X (1,5)
.
1 F X (1,5)
+1,5
+29,5
b
Pr(Y
> 29,5) = 1
+ 1,5
+ 29,5
= 0,0365.
+1,5
!
362
+ 1,5
h(, ) =
,
+ 29,5
h
+ 1,5
+ 1,5
=
ln
,
+ 29,5
+ 29,5
h
( + 1,5)1
= 28
( + 29,5)+1
et
h(, )
= 0,0238
(,
h(, )
= 0,0029,
)
(,
do
167,07 0,0238
= 0,0238 0,0029 23,92
d ,
)]
Var[h(
167,07 1 199,32
0,0029
= 0,00057.
5.33 Le montant pay par lassureur est Y = min(X , 3 000) 100|X > 100, do
f (y + 100)
X
, 0 y < 2 900
1 F X (100)
S X (3 000)
f Y (y) =
, y = 2 900
1 F X (100)
0,
y > 2 900,
,
0 y < 2 900
e 2 900 ,
0,
y = 2 900
y > 2 900.
n
Y
f Y (y i )
i =1
8 (100++1 500)
= e
= 8 e 10 420 .
(e 2 900 )2
363
1 10 420
= 1 302,50.
=
8
y < 150
f X (y),
f Y (y) =
1 F X (150), y = 150
0,
y > 150,
, y < 150
e
e 150 ,
0,
y = 150
y > 150.
n
Y
f Y (y i )
i =1
5 (10++110)
= e
(e 150 )3
= 5 e 845 .
Par la technique usuelle, on trouve = 0,0059.
5.35 On a la fonction de rpartition empirique
0, x 0
2, 0 < x 2
F 9 (x) = 96
9, 2 < x 5
1, 5 < x 8.
2
2 2
6 2
2
5
= 1e
+ 1e
+ 1 e 8 1 .
9
9
364
Chapitre 6
6.1 La fonction de rpartition thorique est
F X (x) = 1
+x
= 1
1 000
1 000 + x
4 (n E )2
X
j
j
j =1
Ej
Soit 23, 0,10 = 6,2514 le 90e centile dune distribution khi carr avec trois degrs
de libert. Puisque 0,7740 < 6,2514, on ne rejette pas le modle.
6.2 La fonction de rpartition du modle est
F X (x) = 1
+x
= 1
50
50 + x
3,5
365
Q=
Or, Pr[22 > 9,5046] = 0,0086 (o 22 est une variable alatoire avec distribution
khi carr avec deux degrs de libert). Par consquent, on ne rejette pas le
modle avec un seuil de signification de 0,86 %. Des seuils proposs, seul 0,5 %
est donc valide.
6.3 La fonction de rpartition empirique est
0,
0,2,
F 5 (x) = 0,4,
0,8,
1,
x < 0,1
0,1 x < 0,4
0,4 x < 0,8
0,8 x < 0,9
x 0,9.
F (x) =
1 + 2y
x
d y = (1 + x),
2
2
0 x 1.
366
0, x < 1
10 , 1 x < 2
6 , 2x <3
F 10 (x) = 10
8
10 , 3 x < 4
10 , 4 x < 8
1, x 8.
b) La fonction de rpartition thorique est
F X (x) = 1
On a donc F (1) = 59 , F (2) = 34 , F (3) =
Cramrvon Mises est
Q CvM =
10
X
2
2+x
.=
21
8
25 , F (4) = 9
et F (8) =
24
25
La distance de
(F (x i ) F 10 (x i ))2
i =1
3 6 2
21 8 2
5 2 2
+4
+2
= 2
9 10
4 10
25 10
2
2
8 9
24
+
+
1
9 10
25
= 0,3478.
c) On a cette fois une fonction de rpartition empirique telle que F 10 (2) =
9
F 10 (4) = 10
et F 10 (8) = 1. La valeur de la distance est donc
Q
CvM
3 6
=
4 10
8 9
+
9 10
24
+
1
25
6
10 ,
= 0,0242.
6.5 On trouve dabord la fonction de rpartition thorique :
x
F (x) =
y
x2
dy = ,
2
4
0 x 2.
On a ensuite F 4 (0,5) = 1/4, F 4 (1) = 2/4, F 4 (1,25) = 3/4, et F 4 (1,5) = 1. Le tableau E.2 prsente les diffrences entre les fonctions de rpartition. La statistique D 4 est donc 7/16 = 0,4375.
367
|F n (x i ) F (x i )|
|F n (x i 1 ) F (x i )|
1
2
3
4
3/16
1/4
0,359375
7/16
1/16
0
0,109375
3/16
TAB. E.2: Diffrences entre les fonctions de rpartition thorique et empirique pour
les donnes de lexercice 6.5
xi
|F n (x i ) F (x i )|
|F n (x i 1 ) F (x i )|
1
4
6
7
8
9
9,5
0,1329
0,1257
0,0686
0,0814
0,0743
0,0471
0,0975
0,0100
0,0171
0,0743
0,0614
0,0686
0,0957
0,0454
TAB. E.3: Diffrences entre les fonctions de rpartition thorique et empirique pour
les donnes de lexercice 6.6
F (x) =
x2
y
dy =
,
50
100
0 x 10.
On a ensuite F 7 (1) = 1/7, F 7 (4) = 2/7, F 7 (6) = 3/7, F 7 (7) = 4/7, F 7 (8) = 5/7,
F 7 (9) = 6/7 et F 7 (9,5) = 1. Le tableau E.3 prsente les diffrences entre les
fonctions de rpartition. La statistique de KolmogorovSmirnov vaut donc
p
D = 0,1329. Puisque la valeur critique du test est c = 1,36/ 7 = 0,5140 > D, on
ne rejette pas le modle.
6.7 On a
F X (x) = 1
+x
= 1
8
.
8+x
368
Q=
(1,25)2
(3; 1,25) = 1 e 1,25 1 + 1,25 +
= 0,1315
2
(5,5)2
5,5
= 0,9116
(3; 5,5) = 1 e
1 + 5,5 +
2
(7)2
(3; 7) = 1 e 7 1 + 7 +
= 0,9704.
2
Le tableau E.4 prsente les calculs pour les deux distributions postules. Pour la
Gamma(3, 0,01), la statistique de KolmogorovSmirnov est D = 0,6616 et pour
la Gamma(3,5, 0,01), la statistique est D = 0,6114. On choisit donc la deuxime
distribution pour la modlisation les donnes.
6.9 Lhypothse de taux dchec constant correspond une distribution exponentielle de paramtre = 0,01. On a donc F X (x) = 1 e x/100 et
E 1 = 50 (F (25) F (0)) = 11,0600
E 2 = 50 (F (40) F (25)) = 5,4240
E 3 = 50 (F (60) F (40)) = 6,0754
E 4 = 50 (F (80) F (60)) = 4,9741
E 5 = 50 (F () F (80)) = 22,4664.
369
|F n (x i ) F (x i )|
|F n (x i 1 ) F (x i )|
xi
=3
= 3,5
=3
= 3,5
125
550
700
0,1185
0,1616
0,0296
0,1771
0,1114
0,0512
0,1315
0,6616
0,2204
0,0729
0,6114
0,1988
TAB. E.4: Diffrences entre les fonctions de rpartition thorique et empirique pour
les donnes de lexercice 6.8
Q=
Puisque Pr[23 > 1,8179] = 0,61 > 0.05, on ne rejette pas lhypothse.
6.10 a) On les valeurs suivantes de la fonction de rpartition empirique : F 50 (25) =
0,20, F 50 (50) = 0,44, F 50 (100) = 0,68 et F 50 (200) = 0,90. Pour la distribution
de Pareto, on a F (25) = 0,4557, F (50) = 0,6464, F (100) = 0,8075 et F (200) =
0,9106. La distance de Cramrvon Mises est alors
Q CvM = (0,4557 0,2)2 + (0,6464 0,44)2
+ (0,8075 0,68)2 + (0,9106 0,9)2
= 0,1244.
Pour la distribution de Weibull, on a F (25) = 0,2212, F (50) = 0,3935, F (100) =
0,6321 et F (200) = 0,8647. La distance est alors
Q CvM = (0,2212 0,2)2 + (0,3935 0,44)2
+ (0,6321 0,68)2 + (0,8647 0,9)2
= 0,0062.
Comme 0,0062 < 0,1244, la distribution de Weibull est un meilleur modle.
370
Q=
Or, 23, 0,05 = 7,815 < Q. On rejette donc le modle avec distribution de
Pareto.
c) Comme 0,10 < 0,1244, le choix de la distribution log-normale serait meilleur.
6.11 a) On a
H0 : numros de dpart quiprobables
H1 : numros de dpart non quiprobables.
b) Pour un total de 144 courses et une probabilit uniforme de victoire de
1
8 , le nombre de victoires espr pour chaque numro est 144/8 = 18.
Les rsultats cumuls observs et esprs sont prsents dans le tableau
suivant :
Numro
Gains observs
Gains thoriques
29
18
48
36
66
54
91
72
108
90
118
108
133
126
144
144
cart absolu
11
12
12
19
18
10
371
Chapitre 7
7.1 a) On trouve dabord lestimateur du maximum de vraisemblance du para k
nk
mtre . On a p k = Pr[N = k] = m
, k = 0, . . . , m et donc
k (1 )
L() =
l () =
m
Y
k=0
m
X
(p k )nk
n k ln p k
k=0
!
!
m
=
n k ln
+ k ln() + (m k) ln(1 )
k
k=0
m
X
k m k
0
.
l () =
nk
1
k=0
m
X
372
E [] = E
m
E [N ]
=
m
m
=
m
= .
b) On a
N
Var[] = Var
m
Var[N ]
=
nm 2
m(1 )
=
nm 2
(1 )
=
.
nm
c) De la partie a), on a
d2
k
m k
ln p k = 2
d 2
(1 )2
do
2
= E [n d ln p N ]
I ()
d 2
n
n(m k)
= 2 E [N ] +
E [m N ]
(1 )2
nm mn(1 )
=
+
(1 )2
nm
=
(1 )
et donc
= I 1 ()
Var[]
=
(1 )
.
nm
373
Var[]
soit
s
z /2
(1 )
.
mn
(1
z /2
.
mn
7.2 a) On a Pr[N = k] = k e /k!, k = 0, 1, . . . , et donc les fonctions de vraisemblance
!
k e n k
Y
L() =
k!
k=0
et de log-vraisemblance
l () =
n k (k ln ln k!).
k=0
= N = Pk=0
= 0,1001
n
k=0 k
= et Var[]
= Var[N ] = /n. On a donc N (, /n). Par conspuis E []
quent, un intervalle de confiance approximatif 95 % pour le paramtre
est
q
1,96
d ],
Var[
= /n.
d ]
avec Var[
Lintervalle de confiance est donc
r
0,1001 1,96
0,1001
.
10 000
374
!n k
Y
k
L() =
k+1
k=0 ( + 1)
X
n k (k ln (k + 1) ln( + 1))
l () =
k=0
et donc
P
kn k
= N = Pk=0
= 0,1001.
n
k=0 k
=
d ]
On trouve ensuite que E []
do Var[
0,1001(1,1001)/10 000. Lintervalle de confiance est donc
s
0,1001(1,1001)
0,1001 1,96
.
10 000
c) En posant = (+1)1 dans les formules de lannexe A, on trouve E [N ] = r
et Var[N ] = r ( + 1). Les estimateurs des moments de r et sont donc les
solutions des quations
P
kn k
r = Pk=0
= 0,1001
n
k=0 k
et
k 2 nk
k=0
P
n
k=0 k
r (1 + ) =
kn k
k=0
P
n
k=0 k
!2
= 0,10028
mu
1.001000e-01
(3.797344e+02) (3.166543e-03)
375
P7
k=0
kn k
7n
= 0,0237.
r = Pk=0
kn k
k=0
nk
= 0,166
et
k 2 nk
k=0
P
n
k=0 k
r (1 + ) =
kn k
k=0
P
n
k=0 k
!2
= 0,2244
mu
0.16600239
(0.21012471) (0.01442188)
7.5 On a
Pr[N = k| = ] f () d
Z
r
=
r +k1 e (+1) d .
(r )k! 0
Pr[N = k] =
376
(r + k)r
(r )k!( + 1)r +k
(r + k)
r
1 k
=
(r )(k 1) + 1
+1
(r + k)
=
r (1 )k ,
(r )(k 1)
7.8 On regroupe les trois dernires classes pour obtenir une frquence significative
pour le calcul de la statistique. Si N Poisson(0,6), on a
E 0 = 365Pr[N = 0] = 200,32
E 1 = 365Pr[N = 1] = 120,19
E 2 = 365Pr[N = 2] = 36,06
E 3+ = 365Pr[N 3]
= 365 E 0 E 1 E 2 = 8,43.
377
2.5
2.0
knk/nk1
1.5
1.0
3 (n E )2
X
j
j
j =0
Ej
Bibliographie
Akaike, H. 1974, A new look at the statistical model identification, IEEE Transactions on Automatic Control, vol. 19, no 6, p. 716723.
Bguin, L.-P. 1990, Lexique gnral des assurances : lexique anglais-franais et
franais-anglais, Cahiers de lOffice de langue franaise, Publications du Qubec,
ISBN 2-55114107-9.
Charbonnier, J. 2004, Dictionnaire de la gestion des risques et des assurances, La
Maison Du Dictionnaire, Paris, ISBN 978-2-85608178-5.
Cossette, H., V. Goulet, M. Jacques et M. Pigeon. 2009, Modlisation des distributions
de sinistres : Exercices et solutions, Document libre publi sous contrat Creative
Commons, ISBN 978-2-9811416-1-3.
Dutang, C., V. Goulet et M. Pigeon. 2008, actuar: An R package for actuarial science, Journal of Statistical Software, vol. 25, no 7. URL http://www.jstatsoft.
org/v25/i07.
Goulet, V. 2012, Introduction la programmation en R, 3e d., Document libre
publi sous contrat Creative Commons, ISBN 978-2-9811416-2-0. URL http:
//cran.r-project.org/other-docs.html.
Goulet, V. 2013, ACT 2002 Mthodes numriques en actuariat, Partie III : Analyse
numrique, cole dactuariat, Universit Laval. URL http://libre.act.ulaval.
ca/cours/act_2002.
Hogg, R. V., A. T. Craig et J. W. McKean. 2005, Introduction to Mathematical Statistics,
6e d., Prentice Hall, Upper Saddle River, NJ, ISBN 0-13008507-3.
Hogg, R. V. et S. A. Klugman. 1984, Loss Distributions, Wiley, New York, ISBN 04718792-9-0.
Klugman, S. A., H. H. Panjer et G. Willmot. 2004, Loss Models: From Data to Decisions, 2e d., Wiley, New York, ISBN 0-4712157-7-5.
379
380
Bibliographie
Klugman, S. A., H. H. Panjer et G. Willmot. 2008, Loss Models: From Data to Decisions, 3e d., Wiley, New York, ISBN 978-0-470-18781-4.
Klugman, S. A., H. H. Panjer et G. Willmot. 2012a, Loss Models: From Data to
Decisions, 4e d., Wiley, New York, ISBN 978-1-118-31532-3.
Klugman, S. A., H. H. Panjer et G. Willmot. 2012b, Solutions Manual to Accompany
Loss Models: From Data to Decisions, 4e d., Wiley, New York, ISBN 978-1-11831531-6.
R Development Core Team. 2013, R: A Language and Environment for Statistical
Computing, R Foundation for Statistical Computing, Vienna, Austria. URL http:
//www.r-project.org.
Silverman, B. W. 1986, Density Estimation for Statistics and Data Analysis, Chapman
& Hall, London, ISBN 0-41224620-1.
Therneau, T. 2013, survival : A Package for Survival Analysis in S. URL http://cran.
r-project.org/package=survival, R package version 2.37-4.
Venables, W. N. et B. D. Ripley. 2002, Modern Applied Statistics with S, 4e d., Springer,
New York, ISBN 0-38795457-0.
Wheeler, B. 2013, SuppDists : Supplementary Distributions. URL http://cran.
r-project.org/package=SuppDists, R package version 1.1-9.
Ce document a t produit avec le systme de mise en page LATEX. Le texte principal est en
Adobe Utopia 11 points, le code informatique en Bera Mono (un clone de Bitstream Vera)
et les titres en Adobe Myriad Pro. Les graphiques ont t raliss avec R.