solution originale

This commit is contained in:
François Pelletier 2014-06-02 19:49:28 -04:00
parent 513b72ef67
commit 6c108afe25
7 changed files with 33273 additions and 0 deletions

109
agents.csv Normal file
View file

@ -0,0 +1,109 @@
Numero_Agent:CodePostal_Agent:Latitude_Agent:Longitude_Agent
1:T0K 0K0:49.1950663:-113.301798
2:T1K 8A3:49.6622436:-112.7929136
3:T1J 0T8:49.6926999:-112.8398053
4:T1G 1S8:49.7845994:-112.1454271
5:T1B 0G4:49.9964255:-110.6439871
6:T1B 4K3:50.0095016:-110.6650225
7:T0L 0T0:50.024926:-113.5836063
8:T1A 6P1:50.0499509:-110.6888123
9:T1V 1Y1:50.5732269:-113.8721543
10:T1R 1E8:50.5775589:-111.8972558
11:T1S 0B8:50.7322168:-113.9509506
12:T2X 0R2:50.8777235:-114.024696
13:T2Y 2Z7:50.905906:-114.066468
14:T2Z 0H2:50.9342099:-113.9600274
15:T2H 0M4:50.975786:-114.0688266
16:T3E 7C4:51.0177079:-114.1619435
17:T2T 1Z5:51.0240451:-114.1094093
18:T1X 1L6:51.0375722:-113.8401517
19:T1P 1M6:51.0388164:-113.4019482
20:T3H 2Z9:51.0486527:-114.1847879
21:T3H 5Z9:51.0594686:-114.2124593
22:T2E 2P6:51.0600361:-114.0627838
23:T2A 4T8:51.0632759:-113.9563593
24:T1W 2B2:51.0892322:-115.3597398
25:T2L 2L4:51.0958997:-114.1406025
26:T3J 0C9:51.124554:-113.9479397
27:T3G 4P5:51.1255584:-114.2104369
28:T3K 1J7:51.127039:-114.0710999
29:T3A 6A7:51.1417566:-114.1320584
30:T3G 5Y6:51.146564:-114.256096
31:T4C 2G4:51.1864754:-114.4746373
32:T4B 2R3:51.2928628:-114.0235277
33:T0J 2J0:51.3550584:-110.4835044
34:T0M 0G0:51.3846891:-113.537125
35:T0J 0Y0:51.4638385:-112.7111381
36:T0J 1P0:51.6440582:-111.92872
37:T4H 1P8:51.7916103:-114.1073735
38:T4H 1P8:51.7916406:-114.1081947
39:T4G 1M1:52.0265077:-113.9512284
40:T4N 3Y6:52.2558949:-113.8141449
41:T4N 1V1:52.2679058:-113.8117233
42:T4P 0M9:52.2869687:-113.7692526
43:T4N 4C9:52.2922951:-113.8138239
44:T0C 2L0:52.3239321:-112.7050294
45:T4T 1B1:52.3730587:-114.9155283
46:T4L 2A6:52.4645222:-113.7334275
47:T0B 0H0:52.5834633:-112.9713008
48:T0B 1N0:52.5840409:-112.0638005
49:T4J 1J3:52.6741235:-113.5795563
50:T0B 4C0:52.7738023:-111.692669
51:T9W 1P5:52.8352589:-110.8591752
52:T9A 3M5:52.9655093:-113.3903635
53:T4V 0K1:53.0180589:-112.8395784
54:T0B 4N0:53.0936131:-111.7784371
55:T0C 0P0:53.10999:-114.4713948
56:T7A 1R8:53.2231787:-114.9834766
57:T0C 2P0:53.2254891:-114.0468373
58:T9V 1H5:53.2703129:-110.0060112
59:T9X 1A8:53.3539758:-110.8580358
60:T4X 1T8:53.3615951:-113.4152885
61:T4X 1T8:53.3616849:-113.4154891
62:T0B 4J0:53.370743:-112.6690703
63:T7V 2A4:53.4236739:-117.5459319
64:T6J 7B9:53.4538133:-113.5098197
65:T6R 3L6:53.4549697:-113.5653848
66:T6H 0K5:53.4880747:-113.5006307
67:T9C 1R7:53.4947266:-112.0527531
68:T6E 6T9:53.504064:-113.4428277
69:T6E 6T9:53.5041123:-113.4428274
70:T5T 1L6:53.5261737:-113.6171847
71:T8A 4W2:53.5290354:-113.290455
72:T7Z 2L6:53.5315832:-114.0006747
73:T8A 4E3:53.537056:-113.320885
74:T7X 4P9:53.5407198:-113.9286684
75:T5N 4C2:53.5471783:-113.5322974
76:T7E 1V1:53.5835146:-116.433637
77:T7E 1V1:53.5846528:-116.4332401
78:T5A 5A1:53.6005516:-113.4218937
79:T5X 5V2:53.6142423:-113.4930047
80:T8N 7C8:53.6202572:-113.6117086
81:T8N 2W5:53.629014:-113.6306592
82:T8N 6G5:53.652502:-113.6287467
83:T8L 1Y8:53.7127179:-113.2127324
84:T0B 2R0:53.7599293:-112.7859936
85:T8R 1A5:53.8021407:-113.6526042
86:T0B 0C0:53.8781325:-112.3359786
87:T0E 1N0:53.953004:-115.1390042
88:T0A 3A3:53.9913448:-111.2791062
89:T7N 1A5:54.1213867:-114.4043008
90:T7S 1P2:54.1402499:-115.6867717
91:T7P 2C7:54.1522529:-113.8490956
92:T9N 2G9:54.2686094:-110.7385026
93:T9M 1P1:54.4142732:-110.2096606
94:T0A 0M0:54.5889981:-112.8077637
95:T9S 1C2:54.7207172:-113.2843934
96:T0H 3N0:55.06916:-117.2811059
97:T8V 7K2:55.1680612:-118.7912393
98:T8V 0T8:55.1707559:-118.7888646
99:T8V 2M7:55.1783861:-118.7945169
100:T0G 2A0:55.2786965:-114.7716474
101:T0H 1M0:55.7378786:-117.2021459
102:T0H 1L0:56.0665277:-118.3841124
103:T0H 1W0:56.1899457:-117.6064166
104:T8S 1S1:56.2323403:-117.2898995
105:T8S 1S1:56.2365884:-117.2900498
106:T9H 3Z7:56.7238513:-111.3717021
107:T9K 1Y1:56.7336707:-111.4259718
108:T0H 1Z0:58.512194:-117.1461249
1 Numero_Agent:CodePostal_Agent:Latitude_Agent:Longitude_Agent
2 1:T0K 0K0:49.1950663:-113.301798
3 2:T1K 8A3:49.6622436:-112.7929136
4 3:T1J 0T8:49.6926999:-112.8398053
5 4:T1G 1S8:49.7845994:-112.1454271
6 5:T1B 0G4:49.9964255:-110.6439871
7 6:T1B 4K3:50.0095016:-110.6650225
8 7:T0L 0T0:50.024926:-113.5836063
9 8:T1A 6P1:50.0499509:-110.6888123
10 9:T1V 1Y1:50.5732269:-113.8721543
11 10:T1R 1E8:50.5775589:-111.8972558
12 11:T1S 0B8:50.7322168:-113.9509506
13 12:T2X 0R2:50.8777235:-114.024696
14 13:T2Y 2Z7:50.905906:-114.066468
15 14:T2Z 0H2:50.9342099:-113.9600274
16 15:T2H 0M4:50.975786:-114.0688266
17 16:T3E 7C4:51.0177079:-114.1619435
18 17:T2T 1Z5:51.0240451:-114.1094093
19 18:T1X 1L6:51.0375722:-113.8401517
20 19:T1P 1M6:51.0388164:-113.4019482
21 20:T3H 2Z9:51.0486527:-114.1847879
22 21:T3H 5Z9:51.0594686:-114.2124593
23 22:T2E 2P6:51.0600361:-114.0627838
24 23:T2A 4T8:51.0632759:-113.9563593
25 24:T1W 2B2:51.0892322:-115.3597398
26 25:T2L 2L4:51.0958997:-114.1406025
27 26:T3J 0C9:51.124554:-113.9479397
28 27:T3G 4P5:51.1255584:-114.2104369
29 28:T3K 1J7:51.127039:-114.0710999
30 29:T3A 6A7:51.1417566:-114.1320584
31 30:T3G 5Y6:51.146564:-114.256096
32 31:T4C 2G4:51.1864754:-114.4746373
33 32:T4B 2R3:51.2928628:-114.0235277
34 33:T0J 2J0:51.3550584:-110.4835044
35 34:T0M 0G0:51.3846891:-113.537125
36 35:T0J 0Y0:51.4638385:-112.7111381
37 36:T0J 1P0:51.6440582:-111.92872
38 37:T4H 1P8:51.7916103:-114.1073735
39 38:T4H 1P8:51.7916406:-114.1081947
40 39:T4G 1M1:52.0265077:-113.9512284
41 40:T4N 3Y6:52.2558949:-113.8141449
42 41:T4N 1V1:52.2679058:-113.8117233
43 42:T4P 0M9:52.2869687:-113.7692526
44 43:T4N 4C9:52.2922951:-113.8138239
45 44:T0C 2L0:52.3239321:-112.7050294
46 45:T4T 1B1:52.3730587:-114.9155283
47 46:T4L 2A6:52.4645222:-113.7334275
48 47:T0B 0H0:52.5834633:-112.9713008
49 48:T0B 1N0:52.5840409:-112.0638005
50 49:T4J 1J3:52.6741235:-113.5795563
51 50:T0B 4C0:52.7738023:-111.692669
52 51:T9W 1P5:52.8352589:-110.8591752
53 52:T9A 3M5:52.9655093:-113.3903635
54 53:T4V 0K1:53.0180589:-112.8395784
55 54:T0B 4N0:53.0936131:-111.7784371
56 55:T0C 0P0:53.10999:-114.4713948
57 56:T7A 1R8:53.2231787:-114.9834766
58 57:T0C 2P0:53.2254891:-114.0468373
59 58:T9V 1H5:53.2703129:-110.0060112
60 59:T9X 1A8:53.3539758:-110.8580358
61 60:T4X 1T8:53.3615951:-113.4152885
62 61:T4X 1T8:53.3616849:-113.4154891
63 62:T0B 4J0:53.370743:-112.6690703
64 63:T7V 2A4:53.4236739:-117.5459319
65 64:T6J 7B9:53.4538133:-113.5098197
66 65:T6R 3L6:53.4549697:-113.5653848
67 66:T6H 0K5:53.4880747:-113.5006307
68 67:T9C 1R7:53.4947266:-112.0527531
69 68:T6E 6T9:53.504064:-113.4428277
70 69:T6E 6T9:53.5041123:-113.4428274
71 70:T5T 1L6:53.5261737:-113.6171847
72 71:T8A 4W2:53.5290354:-113.290455
73 72:T7Z 2L6:53.5315832:-114.0006747
74 73:T8A 4E3:53.537056:-113.320885
75 74:T7X 4P9:53.5407198:-113.9286684
76 75:T5N 4C2:53.5471783:-113.5322974
77 76:T7E 1V1:53.5835146:-116.433637
78 77:T7E 1V1:53.5846528:-116.4332401
79 78:T5A 5A1:53.6005516:-113.4218937
80 79:T5X 5V2:53.6142423:-113.4930047
81 80:T8N 7C8:53.6202572:-113.6117086
82 81:T8N 2W5:53.629014:-113.6306592
83 82:T8N 6G5:53.652502:-113.6287467
84 83:T8L 1Y8:53.7127179:-113.2127324
85 84:T0B 2R0:53.7599293:-112.7859936
86 85:T8R 1A5:53.8021407:-113.6526042
87 86:T0B 0C0:53.8781325:-112.3359786
88 87:T0E 1N0:53.953004:-115.1390042
89 88:T0A 3A3:53.9913448:-111.2791062
90 89:T7N 1A5:54.1213867:-114.4043008
91 90:T7S 1P2:54.1402499:-115.6867717
92 91:T7P 2C7:54.1522529:-113.8490956
93 92:T9N 2G9:54.2686094:-110.7385026
94 93:T9M 1P1:54.4142732:-110.2096606
95 94:T0A 0M0:54.5889981:-112.8077637
96 95:T9S 1C2:54.7207172:-113.2843934
97 96:T0H 3N0:55.06916:-117.2811059
98 97:T8V 7K2:55.1680612:-118.7912393
99 98:T8V 0T8:55.1707559:-118.7888646
100 99:T8V 2M7:55.1783861:-118.7945169
101 100:T0G 2A0:55.2786965:-114.7716474
102 101:T0H 1M0:55.7378786:-117.2021459
103 102:T0H 1L0:56.0665277:-118.3841124
104 103:T0H 1W0:56.1899457:-117.6064166
105 104:T8S 1S1:56.2323403:-117.2898995
106 105:T8S 1S1:56.2365884:-117.2900498
107 106:T9H 3Z7:56.7238513:-111.3717021
108 107:T9K 1Y1:56.7336707:-111.4259718
109 108:T0H 1Z0:58.512194:-117.1461249

32670
assures.csv Normal file

File diff suppressed because it is too large Load diff

20
census_division.csv Normal file
View file

@ -0,0 +1,20 @@
Division,Population,Superficie
1,78694,20526
2,156536,17662
3,38566,13866
4,10078,21467
5,53263,16775
6,1311022,12646
7,40232,1921
8,189243,9909
9,2129,18922
10,93039,2045
11,1203115,15768
12,63427,30047
13,68919,24374
14,28584,26965
15,35983,284
16,67516,97255
17,61504,192116
18,14534,33205
19,109712,2052
1 Division Population Superficie
2 1 78694 20526
3 2 156536 17662
4 3 38566 13866
5 4 10078 21467
6 5 53263 16775
7 6 1311022 12646
8 7 40232 1921
9 8 189243 9909
10 9 2129 18922
11 10 93039 2045
12 11 1203115 15768
13 12 63427 30047
14 13 68919 24374
15 14 28584 26965
16 15 35983 284
17 16 67516 97255
18 17 61504 192116
19 18 14534 33205
20 19 109712 2052

366
cooperators.Rnw Normal file
View file

@ -0,0 +1,366 @@
%
\documentclass[letter]{article}
\usepackage[margin=0.75in]{geometry}
\usepackage{Sweave}
\usepackage{tikz} \usepackage{multirow}
\usepackage[utf8]{inputenc}
\usepackage{amsmath}
\numberwithin{equation}{section}
\usepackage{amsfonts} \usepackage{amsthm} \usepackage{thmtools}
\usepackage{hyperref}
\usepackage{cleveref}
\usepackage[hypcap]{caption}
\usepackage[off]{auto-pst-pdf}
\title{Risque de concentration en Alberta}
\author{François Pelletier}
\begin{document}
\maketitle
<<librairie,echo=FALSE>>=
set.seed(908144032)
library("actuar")
library("multicore")
source("mcsapply.R")
nb.conc <- 50
@
\section{Risque associé aux tempêtes de grêle}
Une tempête de grêle est un évènement de météo extrême qui peut créer beaucoup
de dommages dans une zone concentrée en peu de temps. Plusieurs facteurs
influencent la sévérité des dommages causés par ces tempêtes, dont le principal
est le diamètre moyen des grêlons. Le facteur influençant le nombre d'assurés
qui encourent un dommage est principalement la taille de l'orage de grêle et la
direction des vents.
J'ai choisi de modéliser les orages de grêle selon une force de 1 à 4, qui
influencera à la fois la taille moyenne de l'orage et la taille des grêlons.
Afin de pouvoir obtenir facilement une formule analytique pour déterminer si un
assuré est touché par l'orage, j'ai choisi de lui donner une forme elliptique.
Ainsi, l'orage de force 4 sera formée de 4 ellipses concentriques, chacune ayant
un niveau de dommages allant de 1 à 4.
Je n'ai pas modélisé le lieu de formation des orages. Ainsi, un orage a la même
possibilité de se former à n'importe quel point à l'intérieur d'un rectangle
regroupant les assurés.
Les coordonnées en longitude et en latitude des assurés seront les seules
données utilisées ici.
La position des agents n'est pas nécessaire.
L'ensemble des montants et des paramètres des distributions statistiques ont été
attribués arbitraitement afin de donner des montants raisonnables à mon avis.
Enfin, afin de représenter le risque d'un assuré moins concentré, j'ai
simplement décidé de former un sous-groupe d'assurés formé d'un assuré sur
\Sexpr{nb.conc} dans l'ordre fourni, ce qui a pour effet de diminuer la concentration
du risque tout en conservant l'effet de la densité de population, qui est déjà
concentrée sur le territoire.
\subsection{Définition des ellipses}
L'ellipse est une section conique qui peut être représentée selon une formule
paramétrique \ref{eq:ellipse} qui sert à définir les
coordonnées en fonction des paramêtres. L'ellipse est alors définie par son
centre $X_c,Y_c$ et la mesure de son plus grand rayon $a$ et de son plus petit
rayon $b$, qui sont perpendiculaires. L'inclinaison du plus grand rayon par
rapport à l'abscisse est donnée par l'angle $\varphi$ en radian.
\begin{eqnarray}
\label{eq:ellipse}
X(t)&=X_c + a\,\cos t\,\cos \varphi - b\,\sin t\,\sin\varphi\\
Y(t)&=Y_c + a\,\cos t\,\sin \varphi + b\,\sin t\,\cos\varphi
\end{eqnarray}
La formule canonique \ref{eq:canonique} est une seconde forme qui sert à
déterminer si un point $(x,y)$ fait partie de l'ellipse.
\begin{eqnarray}
\label{eq:canonique}
\frac{\left((x-X_c)\cos \varphi-(y-Y_c)\sin \varphi\right)^2}{a^2} +
\frac{\left((x-X_c)\sin \varphi+(y-Y_c)\cos \varphi\right)^2}{b^2} = 1
\end{eqnarray}
Ces deux dernières équations nous permettrons de simuler des ellipses,
représentant les orages, et de déterminer si des assurés sont touchés par les
orages. Les functions suivantes seront utilisées à cette fin:
<<ellipse1>>=
ellipse_param <- function(t,X,Y,a,b,phi)
{
cbind(X + a * cos(t)*cos(phi) - b*sin(t)*sin(phi),Y + a *
cos(t)*sin(phi) + b*sin(t)*cos(phi))
}
dans_ellipse_param <- function(coord,param)
{
x <- coord[,1]
y <- coord[,2]
xc <- param[1]
yc <- param[2]
a <- param[3]
b <- param[4]
phi <- param[5]
(((x-xc)*cos(phi)-(y-yc)*sin(phi))^2/a^2 +
((x-xc)*sin(phi)+(y-yc)*cos(phi))^2/b^2) <= 1
}
@
\subsection{Importation des données}
On importe des données dans une structure
<<>>=
agents <- read.csv("agents.csv",sep=":")
assures <- read.csv("assures.csv",sep=":")
n.assures <- nrow(assures)
@
On visualise les données sur une carte. Les agents sont représentés par un
triangle rouge et les assurés par une croix noire.
<<fig=TRUE>>=
plot(assures$Longitude_Assure,assures$Latitude_Assure,pch=3,
main="Carte des assurés et des agents",xlab="Longitude",
ylab="Latitude")
points(agents$Longitude_Agent,agents$Latitude_Agent,pch=25,
col="red")
@
On définit la zone où se situent les assurés par un rectangle:
<<>>=
coord_assures <- cbind(assures$Longitude_Assure,
assures$Latitude_Assure)
range_longitude <- c(min(coord_assures[,1]),
max(coord_assures[,1]))
range_latitude <- c(min(coord_assures[,2]),
max(coord_assures[,2]))
@
On définit le groupe restreint pount l'évaluation avec moins de concentration
<<>>=
coord_assures2 <- coord_assures[seq(1,32669,by=nb.conc),]
@
La longitude est située entre \Sexpr{range_longitude[1]} et \Sexpr{range_longitude[2]}.
La latitude est située entre \Sexpr{range_latitude[1]} et \Sexpr{range_latitude[2]}.
\subsection{Hypothèses}
On pose plusieurs hypothèses de modélisation:
<<>>=
franchise <- 500
n.annees <- 100
n.moyen.orages <- 10
longueur.orage <- 1
largeur.orage <- .2
angle.dominant <- 1/2
ecart.angle <- 1
prob_forces <- c(.8, .15, .04, .01)
facteur_forces <- c(1,1.5,2,2.5)
param_dommages <- c(50,1000,20000,400000)
param_forme <- 2
@
\begin{enumerate}
\item La franchise est de \Sexpr{franchise}\$.
\item Le nombre d'années simulées est de \Sexpr{n.annees} années.
\item Le nombre d'orages par année suit une loi Poisson de moyenne
\Sexpr{n.moyen.orages}.
\item La longitude et la latitude du centre de l'orage sont distribués
uniformément sur le rectangle défini précédemment.
\item La longueur et la largeur suivent une loi exponentielle de paramètres
\Sexpr{1/longueur.orage} et \Sexpr{1/largeur.orage}. L'unité est en degrés de
latitude.
\item Les vents dominants sont de direction \Sexpr{angle.dominant} radians.
\item L'angle varie dans un invervalle de \Sexpr{ecart.angle} radian de part et
d'autre selon la distribution $\beta\left(2,2\right)$.
\item La force entre 1 et 4 suit une distribution discrète où les probabilités
sont respectivement \Sexpr{prob_forces[1]},
\Sexpr{prob_forces[2]}, \Sexpr{prob_forces[3]} et \Sexpr{prob_forces[4]}.
\item La force multiplie les dimensions de l'orage respectivement par
\Sexpr{facteur_forces[1]}, \Sexpr{facteur_forces[2]}, \Sexpr{facteur_forces[3]}
et \Sexpr{facteur_forces[4]}.
\item Les dommages suivent une distribution Pareto de paramètre de forme
\Sexpr{param_forme} et de paramètre d'échelle dépendant de la force, prenant les
valeurs
\Sexpr{param_dommages[1]}, \Sexpr{param_dommages[2]}, \Sexpr{param_dommages[3]}
et \Sexpr{param_dommages[4]}.
\end{enumerate}
\subsection{Simulation des orages}
En utilisant les hypothèse définies précédemment, on simule des orages
<<>>=
orages_annees <- rpois(n.annees,n.moyen.orages)
n.orages <- sum(orages_annees)
sim_xc <- runif(n.orages,range_longitude[1],range_longitude[2])
sim_yc <- runif(n.orages,range_latitude[1],range_latitude[2])
sim_a <- rexp(n.orages,1/longueur.orage)
sim_b <- rexp(n.orages,1/largeur.orage)
sim_angle <- angle.dominant - ecart.angle + 2*ecart.angle*
rbeta(n.orages,2,2)
sim_force <- sample(x=1:4,size=n.orages,replace=TRUE,
prob=prob_forces)
sim_param_orages_1 <- cbind(sim_xc,sim_yc,
facteur_forces[1]*sim_a,
facteur_forces[1]*sim_b,sim_angle)
sim_param_orages_2 <- cbind(sim_xc,sim_yc,
facteur_forces[2]*sim_a,
facteur_forces[2]*sim_b,sim_angle)
sim_param_orages_3 <- cbind(sim_xc,sim_yc,
facteur_forces[3]*sim_a,
facteur_forces[3]*sim_b,sim_angle)
sim_param_orages_4 <- cbind(sim_xc,sim_yc,
facteur_forces[4]*sim_a,
facteur_forces[4]*sim_b,sim_angle)
@
Pour chaque orage, on vérifie le nombre d'assurés touchés. On utilise une boucle
ici pour éviter de surcharger la mémoire vive. C'est plus rapide ainsi.
<<>>=
touches1 <- numeric(n.orages)
touches2 <- numeric(n.orages)
touches3 <- numeric(n.orages)
touches4 <- numeric(n.orages)
for(i in 1:n.orages)
{
touches1[i] <- sum((sim_force[i] >= 1) *
dans_ellipse_param(coord_assures,
sim_param_orages_1[i,]))
touches2[i] <- sum((sim_force[i] >= 2) *
dans_ellipse_param(coord_assures,
sim_param_orages_2[i,]))
touches3[i] <- sum((sim_force[i] >= 3) *
dans_ellipse_param(coord_assures,
sim_param_orages_3[i,]))
touches4[i] <- sum((sim_force[i] >= 4) *
dans_ellipse_param(coord_assures,
sim_param_orages_4[i,]))
}
rpareto_tronque <- function(n,shape,scale,deductible)
{
pmax(rpareto(n,shape,scale)-deductible,0)
}
dommages <- mcsapply(mcsapply(as.list(touches1),
rpareto_tronque,param_forme,
param_dommages[1],franchise),sum)+
mcsapply(mcsapply(as.list(touches2),
rpareto_tronque,param_forme,
param_dommages[2],franchise),sum)+
mcsapply(mcsapply(as.list(touches3),
rpareto_tronque,param_forme,
param_dommages[3],franchise),sum)+
mcsapply(mcsapply(as.list(touches4),
rpareto_tronque,param_forme,
param_dommages[4],franchise),sum)
dommages_cum <- cumsum(dommages)
orages_annees_cum <- cumsum(orages_annees)
dommages_annuels <- diff(dommages_cum[orages_annees_cum])
@
On obtient une moyenne de \Sexpr{mean(dommages_annuels)}, un écart-type de
\Sexpr{sd(dommages_annuels)} et les quantiles suivants:
<<>>=
quantile(dommages_annuels,c(.5,.75,.9,.95,.99))
@
\subsection{Comparaison avec moins de concentration}
On sélectionne un assuré sur \Sexpr{nb.conc} afin de réduire l'exposition et on
peut comparar les résultats avec les mêmes orages.
<<>>=
touches21 <- numeric(n.orages)
touches22 <- numeric(n.orages)
touches23 <- numeric(n.orages)
touches24 <- numeric(n.orages)
for(i in 1:n.orages)
{
touches21[i] <- sum((sim_force[i] >= 1) *
dans_ellipse_param(coord_assures2,
sim_param_orages_1[i,]))
touches22[i] <- sum((sim_force[i] >= 2) *
dans_ellipse_param(coord_assures2,
sim_param_orages_2[i,]))
touches23[i] <- sum((sim_force[i] >= 3) *
dans_ellipse_param(coord_assures2,
sim_param_orages_3[i,]))
touches24[i] <- sum((sim_force[i] >= 4) *
dans_ellipse_param(coord_assures2,
sim_param_orages_4[i,]))
}
dommages2 <- mcsapply(mcsapply(as.list(touches21),
rpareto_tronque,param_forme,
param_dommages[1],franchise),sum)+
mcsapply(mcsapply(as.list(touches22),
rpareto_tronque,param_forme,
param_dommages[2],franchise),sum)+
mcsapply(mcsapply(as.list(touches23),
rpareto_tronque,param_forme,
param_dommages[3],franchise),sum)+
mcsapply(mcsapply(as.list(touches24),
rpareto_tronque,param_forme,
param_dommages[4],franchise),sum)
dommages2_cum <- cumsum(dommages2)
dommages2_annuels <- diff(dommages2_cum[orages_annees_cum])
@
On obtient une moyenne de \Sexpr{mean(dommages2_annuels)}, un écart-type de
\Sexpr{sd(dommages2_annuels)} et les quantiles suivants:
<<>>=
quantile(dommages2_annuels,c(.5,.75,.9,.95,.99))
@
\section{Risque associé au vol}
<<>>=
taux_vol_petits <- 1600/100000
taux_vol_grands <- 60/100000
@
\subsection{Hypothèses}
Le risque associé au vol peut être modélisé à l'aide d'une approche
fréquence-sévérité régulière. On pose donc comme hypothèse une distrubution
poisson pour la fréquence et une distribution gamma pour la sévérité.
On divise les vols en petits vols et en grands vols.
<<>>=
vols_petits <- taux_vol_petits * n.assures
freq_vols_petits <- rpois(n.annees,vols_petits)
vols_petits_totaux <- sum(freq_vols_petits)
sev_vols_petits <- rgamma(vols_petits_totaux,2,.005)
sev_vols_petits_cum <- cumsum(sev_vols_petits)
freq_petits_vols_cum <- cumsum(freq_vols_petits)
sev_vols_petits_annuels <- diff(sev_vols_petits_cum[freq_petits_vols_cum])
vols_grands <- taux_vol_grands * n.assures
freq_vols_grands <- rpois(n.annees,vols_grands)
vols_grands_totaux <- sum(freq_vols_grands)
sev_vols_grands <- rgamma(vols_grands_totaux,2,.0001)
sev_vols_grands_cum <- cumsum(sev_vols_grands)
freq_grands_vols_cum <- cumsum(freq_vols_grands)
sev_vols_grands_annuels <- diff(sev_vols_grands_cum[freq_grands_vols_cum])
vols_totaux <- rowSums(cbind(sev_vols_petits_annuels,sev_vols_grands_annuels))
@
On obtient une moyenne de vols annuels de \Sexpr{mean(vols_totaux)} avec un
écart-type de \Sexpr{sd(vols_totaux)}. Les quantiles de la distirubtion sont les
suivants:
<<>>=
quantile(vols_totaux,c(.5,.75,.9,.95,.99))
@
\end{document}

6
ellipse.R Normal file
View file

@ -0,0 +1,6 @@
# TODO: Add comment
#
# Author: francois
###############################################################################

12
mcsapply.R Normal file
View file

@ -0,0 +1,12 @@
require("multicore")
mcsapply <- function (X, FUN, ..., simplify = TRUE, USE.NAMES = TRUE)
{
FUN <- match.fun(FUN)
answer <- mclapply(X = X, FUN = FUN, ...)
if (USE.NAMES && is.character(X) && is.null(names(answer)))
names(answer) <- X
if (!identical(simplify, FALSE) && length(answer))
simplify2array(answer, higher = (simplify == "array"))
else answer
}

90
projet.R Normal file
View file

@ -0,0 +1,90 @@
# TODO: Add comment
#
# Author: francois
###############################################################################
## Importer les données dans un data.frame
agents <- read.csv("agents.csv",sep=":")
assures <- read.csv("assures.csv",sep=":")
## Assurés plus uniformes: en prendre le 1/8, multiplier la fréquence par 8
assures8 <- assures[seq(from=1,to=32669,by=8),]
## Visualiser les agents
plot(agents$Longitude_Agent,agents$Latitude_Agent,pch=6,col="red")
## Visualiser les assurés
points(assures$Longitude_Assure,assures$Latitude_Assure,pch=20)
#fonction orage
sim_dommages_orage <- function(assures)
{
franchise <- 0
## Définir un orage
sim_centre <- cbind(runif(1,range_longitude[1],range_longitude[2]),runif(1,range_latitude[1],range_latitude[2]))
sim_taille <- cbind(rexp(1,.25),rexp(1,4))
sim_angle <- runif(1,0,2*pi)
sim_force <- sample(x=1:4,size=1,replace=TRUE,prob=c(.95,.03,.016,.004))
if(sim_force==1)
{
ellipse1 <- ellipse_param(seq(0,2*pi,by=0.1),sim_centre[1],sim_centre[2],sim_taille[1],sim_taille[2],sim_angle)
}
if(sim_force==2)
{
ellipse2 <- ellipse_param(seq(0,2*pi,by=0.1),sim_centre[1],sim_centre[2],sim_taille[1],sim_taille[2],sim_angle)
ellipse1 <- ellipse_param(seq(0,2*pi,by=0.1),sim_centre[1],sim_centre[2],2*sim_taille[1],2*sim_taille[2],sim_angle)
}
if(sim_force==3)
{
ellipse3 <- ellipse_param(seq(0,2*pi,by=0.1),sim_centre[1],sim_centre[2],sim_taille[1],sim_taille[2],sim_angle)
ellipse2 <- ellipse_param(seq(0,2*pi,by=0.1),sim_centre[1],sim_centre[2],2*sim_taille[1],2*sim_taille[2],sim_angle)
ellipse1 <- ellipse_param(seq(0,2*pi,by=0.1),sim_centre[1],sim_centre[2],3*sim_taille[1],3*sim_taille[2],sim_angle)
}
if(sim_force==4)
{
ellipse4 <- ellipse_param(seq(0,2*pi,by=0.1),sim_centre[1],sim_centre[2],sim_taille[1],sim_taille[2],sim_angle)
ellipse3 <- ellipse_param(seq(0,2*pi,by=0.1),sim_centre[1],sim_centre[2],2*sim_taille[1],2*sim_taille[2],sim_angle)
ellipse2 <- ellipse_param(seq(0,2*pi,by=0.1),sim_centre[1],sim_centre[2],3*sim_taille[1],3*sim_taille[2],sim_angle)
ellipse1 <- ellipse_param(seq(0,2*pi,by=0.1),sim_centre[1],sim_centre[2],4*sim_taille[1],4*sim_taille[2],sim_angle)
}
## Définir le nombre d'assurés touchés
touches1 <- FALSE
touches2 <- FALSE
touches3 <- FALSE
touches4 <- FALSE
touches1 <- dans_ellipse_param(assures$Longitude_Assure,assures$Latitude_Assure,ellipse1)
if(sim_force>1)
{
touches2 <- dans_ellipse_param(assures$Longitude_Assure,assures$Latitude_Assure,ellipse2)
if(sim_force>2)
{
touches3 <- dans_ellipse_param(assures$Longitude_Assure,assures$Latitude_Assure,ellipse3)
if(sim_force>3)
{
touches4 <- dans_ellipse_param(assures$Longitude_Assure,assures$Latitude_Assure,ellipse4)
}
}
}
sum(max(rpareto(sum(touches1*!touches2),2,50)-franchise,0)) +
sum(max(rpareto(sum(touches2*!touches3),2,1000)-franchise,0)) +
sum(max(rpareto(sum(touches3*!touches4),2,20000)-franchise,0)) +
sum(max(rpareto(sum(touches4),2,400000)-franchise,0))
}
set.seed(123)
dommages <- replicate(1000,sim_dommages_orage(assures),simplify=TRUE)
set.seed(123)
dommages8 <- replicate(1000,sim_dommages_orage(assures8),simplify=TRUE)
mean(dommages)
var(dommages)
quantile(dommages,c(.25,.5,.75,.9,.95,.99,.995))
mean(dommages8)
var(dommages8)
quantile(dommages8,c(.25,.5,.75,.9,.95,.99,.995))