Category Archives: Actuarial science

Les tables de mortalité

Chose promise, chose dûe, un court billet expliquant les principales tables utilisées,

  • Les tables TD et TV 88-90

Ces tables datent un peu, et si je continue à en parler en cours, c’est parce qu’elles sont simples à utiliser (et pour continuer à me faire croire que je n’ai pas vieilli depuis mes études). Cette table est d’ailleurs tellement sérieuse qu’on la retrouve dans la loi (ici), dans un arrêté d’avril 1993. La table dite TD 88-90 (pour Décès) a été établie par l’INSEE suite aux observations réalisées entre 1988 et 1990 sur une population d’hommes. Elle était appliquée pour le calcul des primes des contrats d’assurance décès. La table dite TV 88-90 (pour Vie) a été établie par l’INSEE suite aux observations réalisées entre 1988 et 1990 sur une population de femmes. Elle était appliquée pour le calcul des primes des contrats d’assurance en cas de vie. Ces tables peuvent se récupérer à l’aide des codes suivants,

> TD=read.table("https://perso.univ-rennes1.fr/arthur.charpentier/TD8890.csv",
+ sep=";",header=TRUE)
> TV=read.table("https://perso.univ-rennes1.fr/arthur.charpentier/TV8890.csv",
+ sep=";",header=TRUE)

Ces tables ont été remplacées par les tables dites TH et TF, respectivement.

  • Les tables TH et TF 00-02

Ces tables ont été établies à partir des données INSEE de la population française entre 2000 et 2002 et ont été lissées. Ce sont des tables générationnelles, qui nécessitent un correctif d’âge pour tenir compte des écarts de mortalité entre générations. Elles sont applicables à partir du 1er janvier 2006. L’institut des actuaires a proposé une “notice d’utilisation” en ligne ici, et Cimon en avait parlé sur son blog ().

  • Les tables TPRV 95

La table TPRV 93 (pour Table Prospective de Rente Viagère) est un extrait de la table dite plancher pour la tarification des contrats de rente viagère. Elle a été publiée par l’arrêté du 28 juillet 1993 (ici sans les annexes), et correspond à une table prospective qui retrace la mortalité des générations 1887 à 1993 (les tables prospectives sont au programme de Master 2).

La TPRV 93 représente la table complète de la génération 1950. La table est en ligne ici (en csv) lisible sous R avec le code suivant,

> TD=read.table("https://perso.univ-rennes1.fr/arthur.charpentier/TPRV.csv",
+ sep=";",header=TRUE)

Mais ce n’est que la première étape. On utilise ensuite un décallage d’âge.

On a parlé de ce point en TD (mais sur les taux de mortalité), c’est aussi ce qui s’appelle l’hypothèse de Rueff, et traduit la diminution des taux de mortalité par un rajeunissement (i.e. un gain sur les âges), i.e.

Notons qu’en fait,  dépend du niveau des taux d’intérêt. L’ensemble des tables de décallage sont en ligne ici.

  • Les tables TGH et TGF 05

Il s’agit ici, là encore, d’une population de rentiers, comme pour la (ici pour les 38 pages du code). Aussi, ces tables ont été construites sur une population différente des tables TH et TF (construites sur l’ensemble de la population française). La méthodologie est décrite ici.

Actuariat, calcul(s) de PM (partie 1)

Suite aux deux billets sur la modélisation des sinistres en assurance auto (ici et ), un billet pour décrire le contrat d’assurance-vie sur lequel on travaillera en salle machine mercredi.
Pour s’entraîner, on va considérer un contrat d’assurance sur deux têtes, signé par un couple d’âge 50 ans (pour l’homme) et 52 ans (pour la femme). On supposera leur durées de vies résiduelles modélisées par les tables TV88-90 et TD88-90 respectivement. On supposera dans un premier temps que les durées de vie sont indépendantes (hypothèse que l’on relâchera éventuellement par la suite). Il s’agit d’un contrat d’assurance-vie proposant simplement en garantie en cas de vie (ça n’a rien de plus compliqué si l’on rajoute une garantie décès, on utilise l’additivité des VAP). Le couple cotise (annuellement) pendant 20 ans, soit en payant une prime https://perso.univ-rennes1.fr/arthur.charpentier/latex/tp-sum-21.png si les deux époux sont en vie, soit en payant une prime https://perso.univ-rennes1.fr/arthur.charpentier/latex/tp-sum-22.pngsi l’un des époux décède. En contrepartie, l’assureur verse 10 000 euros si les deux époux sont en vie, et 60% de ce montant au conjoint survivant.
Comme dans le cours, le plus simple pour comprendre ce contrat est sûrement de faire un petit dessin. On met en rouge les flux positifs pour l’assureur (il encaisse de l’argent) et en vert les flux négatifs (il dépense de l’argent). Dans un premier temps, les deux conjoints décèdent pendant leur retraite,

Dans un second temps, un des conjoint décède avant le début de la période de retraite,

  • Calculs des valeurs actuelles probables

La valeur actuelle probable de l’assureur s’écrit comme l’espérance de la somme actualisée des engagements de l’assureur,

https://perso.univ-rennes1.fr/arthur.charpentier/latex/tp-sum-23.png

que l’en sépare entre les engagements versés si une ou deux personnes sont en vie. Si l’on regarde les engagements de l’assuré, alors

https://perso.univ-rennes1.fr/arthur.charpentier/latex/tp-sum-24.png

En notant que https://perso.univ-rennes1.fr/arthur.charpentier/latex/tp-sum-25.png, et https://perso.univ-rennes1.fr/arthur.charpentier/latex/tp-sum-26.png (i.e. les deux membres du couple sont en vie à la signature du contrat).
On calcule la prime annuelle pi en notant que les deux VAP doivent être égaux,

https://perso.univ-rennes1.fr/arthur.charpentier/latex/tp-sum-27.png

Pour aller plus loin, il faut expliciter les probabilités jointes, en particulier https://perso.univ-rennes1.fr/arthur.charpentier/latex/tp-sum-28.png (les autres découlant de cette dernière). Pour importer les tables de mortalité, on pourra utiliser

> Lxv=read.table("https://perso.univ-rennes1.fr/arthur.charpentier/TV8890.csv",header=TRUE,sep=";")
> Lxd=read.table("https://perso.univ-rennes1.fr/arthur.charpentier/TD8890.csv",header=TRUE,sep=";")

On poursuivra en salle machine pour calculer tous ces coefficients, mais aussi calculer les provisions mathématiques….. à suivre donc.
PS et je ferais un jour un billet sur les autres tables utilisées, car je continue à utiliser des tables qui sont dépassées (mais qui, à mon avis, m’enlève rien à la modélisation)

Les Pays-Bas inondés, février 1953

Après un week-end de trois jours en Vendée, je vais faire un rapide billet sur les inondations costales (après la catastrophe improbable du barrage de Vajont, ici), que j’aurais du mal à qualifier d’improbables….

Lorsque l’on parle des l’histoire la théorie des valeurs extrêmes (ou plutôt la second partie de l’histoire, la première avait été racontée ici), on évoque souvent la crue de février 1953 qui a inondé les Pays-Bas (je renvoie aux slides de Paul Embrechts par exemple, ici). Pour résumer, pendant la nuit du 31 janvier au 1er février 1953, une très forte tempête traversa la mer du nord. Dans des villes où la marée haute atteignait en moyenne 1 mètre au dessus de niveau de la mer à Amsterdam (ce que l’on appelle NAP, Normaal Amsterdams Peil), on atteignait presque les 4 mètres. Comme souvent (je pense par exemple à la canicule de 2003 en particulier, à laquelle je consacrerais bientôt un billet), le problème n’est pas seulement les hauteurs extrêmes, mais la persistance, i.e. le fait que l’évènement ait duré longtemps. En effet, le 1er février 1953, la marée basse n’avait pas causé de véritable baisse du niveau de l’eau, si bien que quand la marée haute arriva, l’eau avait déjà atteint un niveau très élevé. Mais comme le notait de nombreux commentaires de l’époque, les grandes marées étaient prévues 2 semaines plus tard, autrement dit la catastrophe aurait pu être pire encore si elle était survenue 15jours plus tard !

Les digues subirent plus de 400 brèches, et plus de 500 km de digues furent à reconstruire. Mais surtout, 1836 personnes trouvèrent la mort, ce qui, ramené à la population des Pays-Bas de l’époque (de l’ordre de 10 millions d’habitants, ici), correspondrait à une catastrophe causant la mort de plus de 10000 personnes en France (ou plus de 60000 aux États Unis) aujourd’hui. La carte ci-dessus montre la région qui s’est retrouvée sous l’eau soit plus de 1600 km2 (sinon je renvois ici pour des photos). Comme toujours dans les catastrophes, l’intérêt ex-post sont les leçons que l’on en a tirées. Aux Pays-Bas, la catastrophe a donné naissance au Deltawerken, qui initia (et finança) davantage de recherche sur le niveau des digues, ce que nous appellerions en gestion des risques un problème de calcul de Value-at-Risk: à quelle hauteur mettre une digue pour éviter qu’un évènement se produise dans 99,999% des cas ? Je renvoie ici ou pour plus d’information sur la catastrophe…

Provisionnement et règle de trois

Longtemps, dans les formations d’actuariat, le provisionnement était un sujet traité en 2 heures, car malgré son importance stratégique dans les résultats d’une entreprise d’assurance, il ne s’agissait que d’une “règle de trois“.

  • De l’importance des provisions en assurance

D’un point de vue comptable, les provisions pour sinistres à payer sont une quantification des engagements pris par l’assureur envers ses assurés. En reprenant des vieux tableaux dans l’Argus, avec une compagnie d’assurance spécialisée en grands risques (AXA Corporate Solutions), une mutuelle (MAIF), et une autre société anonyme (Groupama), on obtient

Pour les grands risques, on notera que le montant de provision (plus de 5 milliards) est plus de 3 fois supérieur au montant des primes encaissées (1,5 milliards). Aussi, si l’on pense qu’il faut augmenter de 5% le montant des provisions (ce qui pourrait paraître acceptable à un statisticien) il faudrait augmenter le chiffre d’affaire de 15%.

  • Histoire de la “règle de trois

Le Petit Robert prétendait que l’expression « règle de trois » remontrait à 1538, même si le résultat sous-jacent est tout simplement le théorème de Thalès… Les manuels d’arithmétiques expliquaient ainsi que “la règle de trois est la méthode servant à résoudre les problèmes dans lesquels, étant données des grandeurs pouvant former les trois termes d’une proportion, on demande de trouver le quatrième terme“. On peut retrouver beaucoup plus d’infos sur le site images des mathématiques (ici) ou (), où on trouve d’ailleurs cette illustration,

La “règle de trois” est aussi le premier souvenir de raisonnement mathématique exposée par Paul Lévy dans sa biographie (“quelques aspects de la pensée d’un mathématicien“), lorsqu’il raconte qu’à 7 ans, il tenait la main d’une fille de 12 ans, il était surpris de lever la main. En effet, certes la jeune fille était plus grande que lui, mais son bras devait également être plus grand. Il compris alors le principe des règles de proportionnalité.
Bref, comme toujours, les actuaires ont inventé un nouveau mot pour la “règle de trois” (je renvoie à un vieux papier de Claude Bébéar sur la novlangue utilisée par les actuaires, ici, mais je reviendrais longuement sur le jargon une autre fois).

  • Les méthodes de provisionnement: Chain Ladder

J’avais évoqué il y a quelques temps les travaux d’Eugène Astesan (ici), précurseur de la méthode Chain Ladder, même si l’estimateur qu’il avait alors proposé n’était pas celui de variance minimale dans la classe des estimateurs. Les premiers travaux sur les modèles stochastiques de provisionnement ont chercher à justifier la règle de trois, par exemple via le papier de Thomas Mack “which stochastic model is underlying the chain ladder method” (ici).
Pour venir deux minutes sur la méthode “chain ladder“, par exemple à l’aide du triangle suivant, on peut estimer les montants à payer (dans la partie inférieure) via une simple règle de trois. Considérons le triangle suivant, tiré du livre de Partar, Lecoeur, Nessi,  Nisipasu et Reiz

Pour la seconde colonne, on fait une simple règle de trois, en calculant les sommes sur chacune des colonnes,

Aussi, ici,

i.e.

117,945 = 56,762 x 440,438 / 211,961

On peut ensuite passer à la troisième colonne, puis la quatrième, etc. La règle de trois donne ici

Voici pour la méthode dite “chain ladder“. Et Thomas Mack a été le premier à essayer de quantifier l’erreur de prédiction commise en utilisant cet estimateur. En particulier, les travaux théoriques des années 90 ont essayer d’étudier l’incertitude associée à l’estimation de la charge ultime, calculée aujourd’hui. Comme le notait un article paru dans l’Argus, les compagnies d’assurance qui utilisaient ces techniques étaient encore rares il y a 5 ans (et celles qui le faisaient cherchaient également à le faire savoir, ici). Mais les derniers documents (ici ou ) partent sur un objectif radicalement différent, que j’exposerais d’ici quelques jours, en revenant sur la novlangue des récents papiers sur le provisionnement..

Les modèles en réassurance

Publication d’un papier sur les modèles en réassurance dans la revue Risques (ici), suite à un questionnement sur la pertinence des modèles classiques utilisés par les réassureurs. La question initiale était partie de la constatation – que l’on retrouve ici ou – sur l’utilisation (ou la mauvaise utilisation) de modèles sophistiques en finance de marché, en essayant d’expliquer que – d’un point de vue épistémologique au moins – les modèles des réassureurs étaient plus robustes. En particulier, on notera que les modèles les plus anciens utilisés par les réassureurs (en particulier la loi de Pareto comme je l’avais évoqué ici) ont eu une légitimité pratique durant plusieurs décennies avant d’être justifiés par la théorie des valeurs extrêmes. Je ferais d’ailleurs bientôt un billet sur l’histoire des valeurs extrêmes en statistiques, en revenant en particulier sur les travaux de Gumbel ou de Fréchet.

Sinon le code utilisé dans le papier est en ligne ici, et la base

> sinpe = read.table("https://perso.univ-rennes1.fr/arthur.charpentier/sinpe.csv",header=TRUE,sep=";")
> head(sinpe)
      DSUR  MNTPE
1 19850206 240439
2 19851228 125674
3 19850504 488331
4 19851118 457347
5 19850220 990919
6 19851214 182939
> annee=as.numeric(substr(as.character(sinpe$DSUR),1,4))
> sinistres=sinpe$MNTPE[annee>1992]
> XS=sinistres/100000

On se limite ici aux sinistres survenus après 1992. Si l’on visualise ces montants de sinistres, on obtient

> datesur=as.Date(as.character(sinpe$DSUR),"%Y%m%d")
> jour=datesur[annee>1992]
> plot(jour,sinistres/100000,xlab="",ylab="Coût individuel",cex=.5,ylim=c(0,600))
> ded=50
> abline(h=ded)

On peut aussi faire le graphique Pareto-log-log,

> library(evir)
> n=length(X)
> plot(log(sort(X)),log((n:1)/(n+1)),
+ xlab="Coûts des sinistres (logarithme)",ylab="Fonction de survie (logarithme)",cex=.8)
> out <- gpd(XS, 15)
> XI=as.numeric(out$par.ests[1]); BETA=as.numeric(out$par.ests[2]) 
> x0=seq(2,8,.01)
> lines(x,-1/XI*(x-log(15)),col="red")

(l’ajustement de la loi de Pareto permettant de tracer la droite rouge) ou encore visualiser l’estimateur de Hill de l’indice de queue,

> hill(X)

Pour finir, le code suivant permet de calculer la prime pure (ou plus généralement une prime de Wang) soit de manière non-paramétrique (burning cost) ou bien en utilisant l’ajustement d’une loi de Pareto. Le papier étant un papier de vulgarisation, j’ai pris le seuil de manière arbitraire, sans aucune recherche de valeur “optimale”

> DEDUC = seq(10,50,by=5)
> lambda=0;  seuil=5
> WG1=WG2=rep(NA,length(DEDUC))
> for(k in 1:length(DEDUC)){
+ deductible=DEDUC[k]
+ out <- gpd(XS, seuil)
+ XI=as.numeric(out$par.ests[1]); BETA=as.numeric(out$par.ests[2]) 
+ G0=function(x){1-pgpd(x+seuil, xi = XI, mu = seuil, beta = BETA)}
+ G=function(x){(G0(x+deductible-seuil))/(G0(deductible-seuil))}
+ F=function(x){pnorm(qnorm(G(x))+lambda)}
+ (wang1=integrate(F, 0, Inf))
+ X=XS[XS>deductible]
+ n=length(X)
+ FS= function(z){
+ m=rep(NA,length(z))
+ for(i in 1:length(m)){
+ m[i]=sum(X>z[i]+deductible)/n}
+ return(m)
+ }
+ G=function(x){pnorm(qnorm(FS(x))+lambda)}
+ (wang2=sum(G(seq(0,800,.01))*.01))
+ WG1[k]=as.numeric(wang1$value)
+ WG2[k]=wang2
+ }
> plot(DEDUC,WG2-DEDUC,type='b',xlab="Niveau de la priorité ('00 000 euros)",ylab="Prime pure par sinistres réassuré",ylim=c(0,50))
> lines(DEDUC,WG1-DEDUC,type='b',col="red",pch=4)
> legend(10,50,c("Aujustement d'une loi de Pareto","'Burning cost'"),
+ col=c("red","black"),lty=1,cex=.8,pch=c(4,1))

Market valuation des garanties des produits d’assurance

Pour le dernier cours d’actuariat 2, on parlera de market valuation, i.e. de valorisation financière. Je vais aussi en profiter pour éclairer Bruno (pour qui j’avais déjà commencé à poster des billets sur le contrôle optimal) sur “comment introduit-t-on la mesure risque neutre dans un passif ?” , sous-entendu en assurance-vie (en assurance non-vie c’est à mon avis plus compliqué, mais c’est la question que l’on se pose quand on valoriser des cat-bonds par exemple). Notons que ce problème s’est posé assez tôt dans la littérature économique, je renverrais par exemple au papier de Brennan et Schwartz en 1976.

  • Valorisation de produits financiers (en marché complet)

Je passerais un peu rapidement (en renvoyant ici ou là pour les grands principes), mais les techniques principales pour valoriser les actifs financiers sont les prix d’Arrow-Debreu, les probabilités risque-neutre, et les déflateurs. Mais auparavant, rappelons que le marché est composé de n actifs financiers (bien distincts) dont les prix et les payoffs sont parfaitement connus, et de n états de la natures, de telle sorte que la matrice des payoffs soit de plein rang. Introduisons le concept de portefeuillle de réplication. Compte tenu du fait que la matrice des payoffs est de plein rang, on a l’existence et l’unicité d’un tel portefeuille. Considérons un ensemble de payoff contingents aux états de la nature

https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-valuation-17.png

Il est alors possible de créer un portefeuille (dit de réplication) qui permette de toucher précisément ces payoffs à la date 1. Le prix (à la date 0) d’un tel actif est alors le prix du portefeuille de réplication.
Considérons les actifs rapportant 1 dans un état de la nature, et 0 sinon, à la date 1, i.e.

https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-valuation-16.png

appelés actifs d’Arrow-Debreu. On note https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-valuation-15.png le prix aujourd’hui d’un tel actif. Notons que le prix d’un actif sans risque, rapportant 1 demain dans tous les états du monde est https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-valuation-1.png. Aussi, si on note https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-valuation-14.png le taux sans risque, alors

https://perso.univ-rennes1.fr/arthur.charpentier/latex/correc.png

Le prix de l’actif rapportant https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-valuation-12.png à la date 1 est alors https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-valuation-10.png.
Posons https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-valuation-09.png tel que https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-valuation-08.png soit une distribution de probabilité. Alors

https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-valuation-C-09.png

où l’espérance est prix sous la mesure de probabilité https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-pst.pngappelée probabilité risque-neutre.
Je renvoie ici (par exemple) pour l’introdution de cette probabilité via les arbres binomiaux. Bref, c’est ainsi que l’on arrive au théorème fondamental de la valorisation d’actif en marché complet: le prix d’un actif contingent l’espérance sous probabilité risque neutre, des payoff actualisés.
Aussi, prenons l’exemple d’un put européen, de maturité T et de strike K. Le prix (à la date 0) du put est donné par

https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-valuation-C-01.png

Rappelons d’ailleurs que la formule de Black et Scholes donne un prix pour un tel produit financier (un put européen),

https://perso.univ-rennes1.fr/arthur.charpentier/latex/BSC01.png

https://perso.univ-rennes1.fr/arthur.charpentier/latex/BSC02.png

et

https://perso.univ-rennes1.fr/arthur.charpentier/latex/BSC03.png

Voilà rapidement comment on valorise un produit financier en marché complet, en calculant la valeur d’un portefeuille de réplication.

 

  • Valorisation de produits d’assurance

Pour valoriser les produits d’assurance, classiquement, on calcule des primes pures, ou des valeurs actuelles probables. Par exemple en assurance-vie (car c’est essentiellement de ça dont on parle non ?), considérons simplement une garantie plancher, où l’on s’engage à verser 1€ au décès d’un souscripteur. La valeur actuelle probable (que l’on calcule sous la  probabilité historique) est alors

https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-valuation-C-02.png

que l’on pourrra réécrire

https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-valuation-C-03.png

(où j’adopte une convention de temps continu pour les taux pour garder une analogie avec le cas financier, mais ça n’est pas vraiment important ici). Bref, par rapport à la valorisation financière on valorise toujours des flux futurs (payoffs) actualisés, mais en prenant l’espérance sous la probabilité historique. Et c’est généralement comme ça que l’on valorise un produit d’assurance: les risques ne sont pas échangés sur des marchés liquides, on n’a donc pas vraimentt de marché complet, et les assureurs espèrent simplement que la mutualisation (via la loi des grands nombres) leur permettra d’avoir égalité entre les flux entrants et les flux sortants.
Remarque: j’avais fait l’an dernier un rapide billet sur la transformée d’Esscher (ici) qui présentait une méthode conciliant les deux approches. Mais je suis toujours à la recherche du papier original (datant de 1932) pour proposer une nouvelle discussion sur ce sujet…

  • Valorisation des garanties planchers dans les contrats d’assurance vie, une approche heuristique

Là où les choses se compliquent c’est quand on mélange des garanties d’assurance (un décès) avec des garanties financières (sur des taux ou bien un actif financier pour des contats en unité de compte). Considérons pour cela un contrat décès, où l’assureur s’engage à payer, au moment du décès à la date https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-tau.png le maximum entre la somme investie en actif risqué à la signature du contrat, https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-valuation-07.png et le niveau atteint par l’actif lors du décès, https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-valuation-06.png. Autrement dit, l’assureur s’engage à verser https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-valuation-05.png. De manière plus général, disons qu’il existe une valeur plancher en dessous de laquelle on ne descendra pas, que l’on notera K.
Notons que l’on peut réécrire le payoff

https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-valuation-04.png

Le terme de gauche est la valeur de l’actif, et le second terme peut être interprété comme la valeur de la garantie plancher. Par la suite, on ne s’intéressera qu’à la valeur de cette garantie. Le flux versé (pour cette garantie) à la date t, pour un individu d’âge x à la souscription est alors tout simplement

https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-valuation-C-04.png

(les dates étant supposées discrètes ici). Si l’on cherche à calculer la valeur actuelle probable de l’ensemble des flux, on voit que deux espaces de probabilités intervienne: l’espace des probabilités actuarielles (sous lesquelles on couvre le risque de décès en supposant les risques mutualisables, i.e. l’univers historique https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-valuation-C-05.png) et l’espace des probabilités financières (sous lesquelles on couvre le risque financiers, i.e. l’univers risque neutre https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-valuation-C-06.png). Autrement dit, la valeur actuelle probable s’écrit

https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-valuation-C-07.png

que l’on pourrrait réécrire

https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-valuation-C-08.png

où le terme à droite est le prix d’un call de strike K de maturité t.
Notons que la commission de contrôle s’était penché sur le sujet il y a quelques années, via le mémoire d’actuariat de Sylvain Merlus et Olivier Pequeux (ici), et c’est ce genre de formule qui avait été proposé, à savoir de faire somme somme de prix de put. Et c’est cette approche que  la réglementation a imposé de valoriser cette garantie.
Pour aller plus loin, il avait aussi été précaunisé d’utiliser la formule de Black et Scholes pour valoriser ces puts, ce qui a donné la formule “standard

https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-valuation-03.png

https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-valuation-02.png

et

https://perso.univ-rennes1.fr/arthur.charpentier/latex/market-valuation-01.png

Voilà comment on valorise une telle garantie. Je n’ai parlé de la couverture, car on voit qu’elle ne serait pas parfaite, loin de là.
Comme toutes les formules standards, elle a l’avantage de la simplicité.
Bon formellement, si on veut faire ça proprement, ça risque d’êetre plus compliqué, et il faudrait parler des mesures risques neutre forward, mais on fera ça une autre fois. Sinon tout ça est repris dans une multitude de documents que l’on peut trouver ici, , voire ici ou encore . Et pour aller plus loin, il y a eu récemment pas mal de livres qui abordent ces sujets….

En complément, je pourrais aussi mettre des vieux transparents d’un cours donné il y a plusieurs années maintenant,

Modélisation des tables de mortalité

En cours d’actuariat, nous commencerons à parler des tables de mortalité prospectives.

Côté théorie, je peux renvoyer vers quelques slides de Pierre Devolder, ici. Je renvoie ici pour un court papier de vulgarisation sur le sujet. Sur l’ajustement de lois de mortalité prospectives, je renvoie au très beau mémoire en ligne ici (ou  pour des choses très proches en anglais). Pour l’approche statique, je renvoie également ici pour une méthode simple d’ajustement de loi de Makeham. Pour les aspects dynamique, il y a des documents techniques ici ou .

Sinon pour les données que nous utiliserons en TD, elles sont ici (pour l’exposition, i.e. le nombre de vivants) et  (pour les décès). Il faut retravailler un peu les données, à cause d’une valeur “110+”, mais rien de plus simple,

> DEATH=read.table("https://perso.univ-rennes1.fr/arthur.charpentier/R/FR-Deaths_1x1.txt",header=TRUE)
> EXPOSURE=read.table("https://perso.univ-rennes1.fr/arthur.charpentier/R/FR-Exposures_1x1.txt",header=TRUE)
> DEATH$Age=as.numeric(as.character(DEATH$Age))
Message d'avis :
NAs introduits lors de la conversion automatique 
> DEATH$Age[is.na(DEATH$Age)]=110
> EXPOSURE$Age=as.numeric(as.character(EXPOSURE$Age))
Message d'avis :
NAs introduits lors de la conversion automatique 
> EXPOSURE$Age[is.na(EXPOSURE$Age)]=110

Coté code, je renvoie à un précédant billet sur le sujet, ici, et  pour le package (non officiel) de démographie développé par Rob Hyndman.

Enfin, à la fin du cours, on parlera un peu de modèles de durée. J’essayerai de taper quelques choses dans les jours à venir, mais côté TD, on utilisera ce que fait Frédéric Planchet (qui connaît fort bien le sujet), et qui se trouve en ligne ici. Les données sont en ligne ici, avec un peu de code R pour retraiter les donnéesici et , et un peu de code pour des modèles classiques, comme du modèle de Cox (ici) ou du Kaplan Meier ().

Les codes R utilisés en TD sont les suivants (avec en plus les graphiques obtenus),

> D=DEATH[DEATH$Year==an,]
> E=EXPOSURE[EXPOSURE$Year==an,]
> MU = D[,3:5]/E[,3:5]
> plot(0:110,log(MU[,1]),type="l",col="red")
> lines(0:110,log(MU[,2]),col="blue")

Avec ci-dessus les log des taux de décès instantannés pour l’année 1986.

> PH=matrix(NA,111,111) # ligne X, colonne H
> PF=PH
> for(x in 0:110){
+ PH[x+1,1:(111-x)]=exp(-cumsum(MU[(x+1):111,2]))
+ PF[x+1,1:(111-x)]=exp(-cumsum(MU[(x+1):111,1]))
+ }
> x=0
> plot(1:111,PH[x+1,],ylim=c(0,1),type="l",col="blue")
> lines(1:111,PF[x+1,],col="red")

pour la fonction de survie à la naissance (ci-dessus) ou à 50 ans (ci-dessous),

> somme=function(x){sum(x,na.rm=TRUE)}
> EH=EF=rep(NA,111)
> EH=apply(PH,1,somme)
> EF=apply(PF,1,somme)
> plot(0:110,EH,type="l",col="blue")
> lines(0:110,EF,col="red")

Ce qui correspond aux espérances de vie résiduelles.

> for(k in 1:3){
+ Z=MU[,k]
+ Z[is.nan(Z)==TRUE]=NA
+ MU[,k]=Z
+ }
> AGE=unique(DEATH$Age)
> ANNEE=unique(DEATH$Year)
> MUF=matrix(MU[,1],length(AGE),length(ANNEE))
> MUH=matrix(MU[,2],length(AGE),length(ANNEE))
> persp(AGE[1:100],ANNEE,log(MUH[1:100,]),
+ theta=-30,col="light green",shade=TRUE)
> persp(AGE[1:100],ANNEE,log(MUH[1:100,107:1]),
+ theta=-30,col="light green",shade=TRUE)
> image(AGE,ANNEE,log(MUH))

On a retourné le graph en mettant les années à l’envers, mais on peut aussi regarder la surface pour les femmes, ainsi que l’allure des courbes de niveau

> persp(AGE[1:100],ANNEE,log(MUF[1:100,107:1]),
+ theta=-30,col="light green",shade=TRUE)
> image(AGE,ANNEE,log(MUF))

Mais nous n’avons ici qu’une lecture statique des tables, et nous ne suivons pas vraiment un individu: en même temps que la personne vieilli, le temps passe…. Considérons une personne née en 1900, et regardons (ex post) l’évolution de sa probabilité de décès (en regardant non pas la probabilité de décéder à 18 d’une personne de 18 ans en 1900, mais en tenant compte du fait que la personne a fêté ses 18 ans en 1918).

> MUFP=MUHP=rep(NA,111)    # MU PROSPECTIF
> for(i in 1:111){
+     MUHP[i]=MUH[i,colonne+i]
+     MUFP[i]=MUF[i,colonne+i]
+ }
> plot(0:110,log(MUHP),col="purple",type="l")
> D=DEATH[DEATH$Year==annee.naissance,]
> E=EXPOSURE[EXPOSURE$Year==annee.naissance,]
> MU = D[,3:5]/E[,3:5]
> lines(0:110,log(MU[,2]),col="blue")

On peut aussi l’ajustement par un modèle à la Lee-Carter (via la fonction fournie par LifeMetrics)

> source("https://perso.univ-rennes1.fr/arthur.charpentier/fitModels.r")
> DEATH=DEATH[DEATH$Age<90,] # virer age>90
> EXPOSURE=EXPOSURE[EXPOSURE$Age<90,] # virer age>90
> XV=unique(DEATH$Age)
> YV=unique(DEATH$Year)
> ETF=t(matrix(EXPOSURE[,3],length(XV),length(YV)))
> DTF=t(matrix(DEATH[,3],length(XV),length(YV)))
> ETH=t(matrix(EXPOSURE[,4],length(XV),length(YV)))
> DTH=t(matrix(DEATH[,4],length(XV),length(YV)))
> WA=matrix(1,length(YV),length(XV))
> LCF = fit701(xv=XV,yv=YV,etx=ETF,dtx=DTF,wa=WA)
-1166100 -> -263940.2  ->-204175.2 
-188375.9 -> -185822.4  ->-184403.5 
-183223.3 -> -183069.8  ->-182350.3 
-181701.1 -> -181600.6  ->-181145.5 
-180729.6 -> -180671.5  ->-180371.8 
-180101.9 -> -180068.6  ->-179871.2 
> LCH = fit701(xv=XV,yv=YV,etx=ETH,dtx=DTH,wa=WA)
-2564403 -> -1028776  ->-953629.9 
-881334.2 -> -868279.3  ->-834413.6 
-804946.8 -> -802817.2  ->-786343 
-772775 -> -772007.5  ->-764501.7 
-758653.8 -> -758399  ->-755212 
-752931 -> -752851  ->-751620.2 
-750811.1 -> -750786.7  ->-750351.4 
> plot(LCF$x,LCF$beta1,type="l",col="red")
> lines(LCH$x,LCH$beta1,col="blue")

> plot(LCF$x,LCF$beta2,type="l",col="red")
> lines(LCH$x,LCH$beta2,col="blue")

> plot(LCF$y,LCF$kappa2,type="l",col="red")
> lines(LCH$y,LCH$kappa2,col="blue")

> MUF      = DTF/ETF
> MUH      = DTH/ETH
> MUFMODEL = MUF; MUHMODEL = MUH
> for(i in 1:nrow(MUF)){
+ for(j in 1:ncol(MUF)){
+ MUFMODEL[i,j] = exp(LCF$beta1[j]
+ +LCF$beta2[j]*LCF$kappa2[i])
+ MUHMODEL[i,j] = exp(LCH$beta1[j]
+ +LCH$beta2[j]*LCH$kappa2[i])
+ }}

soit, pour avoir un dessin,

> persp(LCH$y,LCH$x,log(MUH),
+ theta=45,col="light green",shade=TRUE)
> persp(LCH$y,LCH$x,log(MUHMODEL),
+ theta=45,col="light blue",shade=TRUE)

On peut ainsi comparer les taux bruts (en vert) avec les taux estimés (en bleu)

mais aussi regarder  se trouvent les erreurs le plus grandes,

Peut on faire l’économie du formalisme quand on parle d’extrêmes ?

Tous les blogs économiques saluent la parution en poche du joli petit livre de Daniel Zajdenweber, Economie des Extrêmes. En particulier, beaucoup de monde salue ce livre qui explique simplement des choses complexes…. Par exemple Alexandre dès 2001 “Passé le premier chapitre, un peu ardu, et qui nécessite du lecteur des connaissances de base en statistique et probabilités (notion de lois de probabilité, d’espérance, de variance…) qui décrit en termes littéraires les caractéristiques de ces lois, l’auteur applique ces résultats à un grand nombre de phénomènes concrets, et en tire les conséquences“. Mais peut-on parler d’économie des extrêmes sans être technique ?

Histoire que mon message ne soit pas déformé, je trouve passionnant ce petit livre introductif à la problématique des risques extrêmes (qui est un de mes dadas depuis quelques années) mais j’espère qu’il servira d’encouragement à une lecture d’ouvrages plus détaillés sur le sujet. Car la vulgarisation a des limites que l’on atteint vite quand on parle de sujets aussi complexes.

L’exemple que j’ai le plus étudié est celui des sinistres de perte d’exploitation (longuement évoqué par Daniel Zajdenweber dans son livre). Il y a quelques années j’avais utilisé cette partie du livre comme base pour faire un sujet d’examen pour le cours de “réassurance et grands risques” que je donnais alors à l’ENSAE1. Et malheureusement, mes compétences littéraires sont très limitées, donc je vais faire des maths. Dans le livre, la figure suivante est présentée,

qui correspond effectivement à la fonction tracée dès 1925 par Karl Gustav Hagstroem (j’avais souligné (ici) ses travaux précurseurs où l’intérêt de la loi de Pareto pour modéliser les très grands riques apparaissait pour la première fois). C’est en effet assez naturel: si on a une loi de Pareto, i.e.

https://perso.univ-rennes1.fr/arthur.charpentier/latex/z01.png

alors on pourrait écrire, en passant au logarithme

https://perso.univ-rennes1.fr/arthur.charpentier/latex/z02.png

Si on représente la version empirique, c’est à dire le nuage de points

https://perso.univ-rennes1.fr/arthur.charpentier/latex/z03.png

alors pour une loi de Pareto, les points devraient être alignés suivant une droite, et la pente doit correspondre au paramètre de la fonction puissance. C’est visiblement l’idée exploitée ici.

Autrement dit, les pointillés ne sont un intervalle de confiance, mais juste un outils graphique pour se demander si la pente vaut 1, ou pas. Daniel Zajdenweber affirme que la pente doit ici être -1.

Le fait que la valeur soit unitaire ou pas a en effet un impact très important en terme d’assurabilité du risque de perte d’exploitation. Rappelons que pour une variable positive (et c’est le cas ici). Et si on a une telle loi de Pareto (de puissance unitaire), alors la prime pure d’un traité de réassurance, couvrant entre m et M s’écrit

https://perso.univ-rennes1.fr/arthur.charpentier/latex/z06.png

soit

https://perso.univ-rennes1.fr/arthur.charpentier/latex/z07.png

ce qui correspond aux calculs de Daniel Zajdenweber… Mais encore une fois “l’absence d’espérance mathématique de la distribution des sinistres” est une conclusion très forte sur laquelle on peut essayer de revenir.

https://perso.univ-rennes1.fr/arthur.charpentier/latex/z04.png

aussi ici

https://perso.univ-rennes1.fr/arthur.charpentier/latex/z05.png

autrement dit, l’espérance est finie si la pente est strictement plus grande (en valeur absolue) que 1. Si la pente est inférieure (ou égale) à 1, le risque n’est pas assurable ! Ce qui est une conclusion très très forte pour les assureurs.
J’ai donc demandé à la FFSA la base de données utilisée ici, et pour éviter des problèmes d’inflation des coûts de sinistres entre 1992 et 2000. Si je prends tous les sinsitres, on obtient l’ajustement de Pareto suivant

soit une pente (en valeur absolue de 1.47). Mais encore une fois, l’ajustement de Pareto se fait sur les grands sinistres. Hill a proposé un estimateur très populaire pour estimer ce coefficient, où on ne prend en compte que les k observations les plus grandes, et on regarde l’estimation de la pente du graphique de Pareto pour ces quelques valeurs. On représente alors l’estimation en fonction du nombre de grands sinistres, ou du seuil définissant les graphs sinistres. Numériquement, en posant

https://perso.univ-rennes1.fr/arthur.charpentier/latex/z08.png

on peut écrire comme estimateur de la pente

https://perso.univ-rennes1.fr/arthur.charpentier/latex/z09.png

soit, en simplifiant le numérateur,

https://perso.univ-rennes1.fr/arthur.charpentier/latex/z10.png

tel que l’a construit Hill en 1975. Graphiquement, on a ici

Bref, la question est ce savoir si on atteint la valeur 1 pour les grands sinistres. Graphiquement, on a malgré tout envie de rejeter cette hypothèse.
Une solution peut être de faire un test statistique, basé sur de ratio de vraisemblance, comme le suggèrent Reiss & Thomas (2001) or Coles (2001). En fait, on peut même utiliser d’autres estimateurs que celui de Hill, comme celui obtenu en faisant un ajustement de loi GPD (Pareto généralisée
) sur la loi des Excès, ou une loi GEV sur des maximas par blocs (Generalized Extreme Value). On introduit alors la statistique de test suivante

https://perso.univ-rennes1.fr/arthur.charpentier/latex/z12.png

et on regarde les p-value (ainsi que la correction de Bartlett à droite),

On peut aussi, plus simplement, estimer plusieurs coefficients de pentes pour des seuils différents, et regarder la borne supérieure de l’intervalle de confiance,

Bref, même si avec un des ajustements de loi GPD on hésite à retenir une pente unitaire, la plupart des tests rejettent cette hypothèse, et donc le risque de perte d’exploitation semble être assurable, d’espérance mathématique finie. Bref, les dessins c’est très bien pour faire passer une idée, mais ne retenir que ça pour en tirer des conclusions aussi fortes me laisser sceptique….

Incertitude(s) dans les triangles de paiements

Un petit billet aujourd’hui sur les méthodes de provisionnements, et sur les calculs d’incertitude à 1 an, ou à l’ultime, de montants de provisions pour sinistres restant à payer…. Ce sujet avait été évoqué ici (surtout dans les commentaires à l’époque).

  • Mack et incertitude à long terme 

Le mse dans le cadre du modèle de Mack est proposé dans la library(ChainLadder), via la fonction MackChainLadder(). Sur le triangle utilisé dans les exemples et formations en provisionnement, on obtient les mse estimés suivants

> source("https://perso.univ-rennes1.fr/arthur.charpentier/bases.R")
> library(ChainLadder)
> MackChainLadder(PAID)
MackChainLadder(Triangle = PAID)

  Latest Dev.To.Date Ultimate    IBNR Mack.S.E CV(IBNR)
1  4,456       1.000    4,456     0.0    0.000      NaN
2  4,730       0.995    4,752    22.4    0.639   0.0285
3  5,420       0.993    5,456    35.8    2.503   0.0699
4  6,020       0.989    6,086    66.1    5.046   0.0764
5  6,794       0.978    6,947   153.1   31.332   0.2047
6  5,217       0.708    7,367 2,149.7   68.449   0.0318

              Totals
Latest:    32,637.00
Ultimate:  35,063.99
IBNR:       2,426.99
Mack S.E.:     79.30
CV(IBNR):       0.03

Mais il s’agit fondamentalement d’une formule sur l’estimation ultime de la charge. Pour une année de survenance donnée, la formule dérivée par Thomas Mack ressemble à ça

(honteusement emprunté dans des slides de Mario Wüthrich, par paresse, ici), où pour rappel,

  • Muerz & Wüthrich et incertitude à court terme

Pour faire simple, il existe une formule fermée, un peu complexe, et une simplification, reposant sur l’idée que, si on note

alors, Michael Merz et Mario Wüthrich ont montré que

Le code sur le triangle mentionné ici donne le résultat suivant (ce code est parfois appelé dans divers documents MackMerzWuthrich.R)

> source("https://perso.univ-rennes1.fr/arthur.charpentier/merz-wuthrich-triangle.R")
> MSEP_Mack_MW(PAID,0)
   MSEP Mack MSEP observable approchée MSEP observable exacte
1  0.0000000                  0.000000               0.000000
2  0.6393379                  1.424131               1.315292
3  2.5025153                  2.543508               2.543508
4  5.0459004                  4.476698               4.476698
5 31.3319292                 30.915407              30.915407
6 68.4489667                 60.832875              60.832898
7 79.2954414                 72.574735              72.572700

La dernière ligne (notée ici 7 correspond au total – toutes années de survenance confondues). Ici, on regarde l’erreur à un an. Et Mario Wüthrich concluait – lors d’une jolie école d’été à Lyon il y a quelques années – en notant le point suivant

Sur le triangle étudié ici, on était à 91%. Mais supposer que l’on soit inférieur à 1 est un point qu’il me semble intéressant de discuter…

  • incertitude plus grande à court terme ?

Dans un billet datant un peu (ici), on laissait entendre qu’il était sûrement possible que l’incertitude à court terme excède l’incertitude à long terme.
Arnaud Lacoume a trouvé un vrai triangle dans lequel l’incertitude de court terme excède l’incertitude à long terme (détaillé ici dans une feuille xls). Pour information, il s’agit du risque GAV, garanties accidents de la vie, bref, rien de très exotique…

> Triangle = read.table("https://perso.univ-rennes1.fr/arthur.charpentier/GAV-triangle.csv",sep=";")

Si on utilise le code gentiment offert par Arnaud, on note qu’effectivement, l’incertitude à 1 an est grande, très grande…

> MSEP_Mack_MW(Triangle,0)
    MSEP Mack MSEP observable approchée MSEP observable exacte
1        0.00                       0.0                    0.0
2    12459.74                  129692.2               152605.9
3   209431.14                  214136.5               214419.6
4   258003.38                  198072.3               198773.0
5  3055297.40                 3048489.5              3065525.1
6 58429393.29                57056117.3             67375794.0
7 58669646.13                57301552.4             67586106.6

On note ici que le MSE exact à un an (car Merz et Wüthrich proposent aussi une formule approchée, basée sur une approximation qui contraint l’estimation à un an à être inférieure à l’estimation à horizon ultime) est – si on regarde toutes années de survenance confondues – 15% plus élevé que le MSE à ultime. Avec un résultat du même ordre pour la dernière année (ce n’est donc pas un problème d’agrégation). On peut jouer à simuler des triangles, via un processus de déroulé de paiements.

Regardons par exemple le code suivant pour simuler des triangles. L’idée est ici que l’on n’a pas d’incertude sur la valeur à l’ultime,

> M=matrix(NA,10,10)
> for(i in 1:10){
+ X=rexp(10)
+ Y=cumsum(X)/sum(X)
+ M[i,1:(11-i)]=Y[1:(11-i)]*100}
> Triangle = M

On tire des variables positives, car on veut un modèle à incréments positifs. Par exemple, une des simulations a donné le triangle suivant

> Triangle
      [,1] [,2] [,3] [,4] [,5] [,6] [,7] [,8] [,9] [,10]
 [1,] 35.4 36.3 38.2 47.6 50.7 66.1 71.3 72.1 98.7   100
 [2,] 15.5 20.2 22.3 23.7 63.9 68.2 75.9 79.3 95.1    NA
 [3,]  4.1 30.3 32.3 33.3 49.4 53.1 62.8 77.1   NA    NA
 [4,]  1.1 19.4 38.9 41.9 48.0 66.5 68.3   NA   NA    NA
 [5,]  0.3  7.2 10.8 17.8 28.1 30.1   NA   NA   NA    NA
 [6,] 15.6 32.2 46.0 50.6 68.1   NA   NA   NA   NA    NA
 [7,]  4.0  9.9 10.6 15.8   NA   NA   NA   NA   NA    NA
 [8,]  6.1 24.0 29.4   NA   NA   NA   NA   NA   NA    NA
 [9,]  0.0  0.5   NA   NA   NA   NA   NA   NA   NA    NA
[10,]  2.3   NA   NA   NA   NA   NA   NA   NA   NA    NA

En faisant le calcul sur les mesures d’incertitude à la Merz & Würthrich, on obtient une très grande incertitude à court terme,

> MSEP_Mack_MW(Triangle,0)
   MSEP Mack MSEP observable approchée MSEP observable exacte
1   0.000000                  0.000000               0.000000
2   6.820602                 12.776729             118.050196
3  13.377536                 13.119557              13.126625
4  17.415678                 13.638880              13.651131
5  11.797752                  5.574916               5.580840
6  26.767233                 18.407916              18.429069
7  27.184076                 23.725348              23.773581
8  44.100071                 16.410374              16.444281
9   7.227724                  3.450343               3.460982
10 54.012270                 49.443633              49.593284
11 99.184167                 85.042393             145.008282

Ici je suis 50% plus élevé ! Bon, mais c’est beau de montrer pourquoi une idée partagée par beaucoup de monde n’est pas toujours vraie…mais il reste à comprendre plus en détails le pourquoi … à suivre donc !

En attendant, je mets ici le code utilisé pour quantifier cette incertitude à un

### MERZ & WUTHRICH

MWChainLadder <- function(Trian, bool = 0, warning = FALSE){    

##  variables
I <- nrow(Trian)
J <- ncol(Trian)
diag <- matrix(NA,I,1)
diag_inv <- matrix(NA,I,1)
S_I <- matrix(NA,1,J)
S_II <- matrix(NA,1,J)
Phi <- matrix(NA,I,1)
Psi <- matrix(NA,I,1) 
Delta <- matrix(NA,I,1)
Gamma <- matrix(NA,I,1)
Lambda <- matrix(NA,I,1)
Upsilon <- matrix(NA,I,1)
cov_obs <- matrix(0,I,J)
cov_reel <- matrix(0,I,J)
msep_obs <- matrix(NA,1,J+1)
msep_reel <- matrix(NA,1,J+1)
Delta[1] <- Phi[1] <- Psi[1] <- Upsilon[1] <- Lambda[1] <- 0

Mack<-MackChainLadder(Trian)
Mack$sigma[J-1] <- sqrt(min(Mack$sigma[J-2]^4/
Mack$sigma[J-3]^2,min(Mack$sigma[J-2]^2,Mack$sigma[J-3]^2)))
#plot(Mack)

for (i in 1:I){
    diag[i] = Trian[i,I-i+1]
    diag_inv[i] = Trian[I-i+1,i]
}

for (j in 1:J){
    S_I[j] <- sum(Trian[1:(I-j),j])
    S_II[j] <- sum(Trian[1:(I-j+1),j])
}
S_I[I] <- 0

Delta[2] <- Mack$sigma[I-1]^2/(S_I[I-1]*(Mack$f[I-1])^2)
Phi[2] <- 0
Psi[2] <- Mack$sigma[I-1]^2/(diag[2]*(Mack$f[I-1])^2)
Upsilon[2] <- Mack$sigma[I-1]^2/(S_II[I-1]*(Mack$f[I-1])^2)
Lambda[2] <- diag[2]*Mack$sigma[I-1]^2/((Mack$f[I-1]^2)*S_II[I-1]*S_I[I-1])
for (i in 3:I){
    Delta[i] <- Mack$sigma[I-i+1]^2/(S_I[I-i+1]*(Mack$f[I-i+1])^2) + sum(
            (diag_inv[(I-i+2):(J-1)]/S_II[(I-i+2):(J-1)])^2*
         Mack$sigma[(I-i+2):(J-1)]^2/(S_I[(I-i+2):(J-1)]*(Mack$f[(I-i+2):(J-1)])^2))        
    Phi[i] <- sum((diag_inv[(I-i+2):(J-1)]/S_II[(I-i+2):(J-1)])^2*
           Mack$sigma[(I-i+2):(J-1)]^2/(diag_inv[(I-i+2):(J-1)]*(Mack$f[(I-i+2):(J-1)])^2))    
    Psi[i] <- Mack$sigma[I-i+1]^2/(diag[i]*(Mack$f[I-i+1])^2)
    Upsilon[i] <- Phi[i] + Mack$sigma[I-i+1]^2/(S_II[I-i+1]*(Mack$f[I-i+1])^2)

Lambda[i] <- diag[i]*Mack$sigma[I-i+1]^2/((Mack$f[I-i+1]^2)*S_II[I-i+1]*S_I[I-i+1]) + sum(
            (diag_inv[(I-i+2):(J-1)]/S_II[(I-i+2):(J-1)])^2*
           Mack$sigma[(I-i+2):(J-1)]^2/(S_I[(I-i+2):(J-1)]*(Mack$f[(I-i+2):(J-1)])^2))
}
Gamma = Phi + Psi

#MSEP par acc. year
for (i in 1:I){
    msep_obs[i] = Mack$FullTriangle[i,J]^2 * (Gamma[i] + Delta[i])
    msep_reel[i] = Mack$FullTriangle[i,J]^2 * (Phi[i] + Delta[i])
}

#Covariance
for (i in 2:(I-1)){ 
        for (k in (i+1):I){
        cov_obs[i,k] <- Mack$FullTriangle[i,J]*Mack$FullTriangle[k,J]*(Upsilon[i] + Lambda[i])
        cov_reel[i,k] <- Mack$FullTriangle[i,J]*Mack$FullTriangle[k,J]*(Phi[i] + Lambda[i])    
    }
}

#MSEP aggregated
msep_obs[I+1] = sum(msep_obs[2:I]) + 2*sum(cov_obs)
msep_reel[I+1] = sum(msep_reel[2:I]) + 2*sum(cov_reel)

##################################################################################################################
## EXACTES
##################################################################################################################
#DÈclaration des variables
facteur <- matrix(1,I,1)
Phi_exact <- matrix(NA,I,1)
Psi_exact <- matrix(NA,I,1)
Gamma_exact <- matrix(NA,I,1)
Upsilon_exact <- matrix(NA,I,1)
cov_obs_exact <- matrix(0,I,J)
cov_reel_exact <- matrix(0,I,J)
msep_obs_exact <- matrix(NA,1,J+1)
msep_reel_exact <- matrix(NA,1,J+1)
Phi_exact[1] <- Upsilon_exact[1]<- Gamma_exact[1] <- 0

# i=2
Phi_exact[2] <- Psi_exact[2] <- 0
Gamma_exact[2] <- Mack$sigma[I-1]^2/(diag[I-1]*(Mack$f[I-1])^2)
Upsilon_exact[2] <- Mack$sigma[I-1]^2/(S_II[I-1]*(Mack$f[I-1])^2)

# loop for prior
for (i in 3:I){
    facteur[i] <-prod(1 + diag_inv[(I-i+2):(J-1)]*(Mack$sigma[(I-i+2):(J-1)]^2)/((S_II[(I-i+2):(J-1)]*
            Mack$f[(I-i+2):(J-1)])^2))
    Phi_exact[i] <- (1 + Mack$sigma[I-i+1]^2/(diag_inv[I-i+1]*(Mack$f[I-i+1])^2))* (facteur[i]-1)    
    Psi_exact[i] <- (1 + Mack$sigma[I-i+1]^2/(S_II[i]*(Mack$f[I-i+1])^2)) * Phi[i] /
                (1 + Mack$sigma[I-i+1]^2/(diag_inv[I-i+1]*(Mack$f[I-i+1])^2))
    Upsilon_exact[i] <- (1 + Mack$sigma[I-i+1]^2/(S_II[I-i+1]*(Mack$f[I-i+1])^2))* facteur[i]-1    
    Gamma_exact[i] <- (1 + Mack$sigma[I-i+1]^2/(diag_inv[I-i+1]*(Mack$f[I-i+1])^2))* facteur[i]-1    
}

#MSEP per accident year
for (i in 1:I){
    msep_obs_exact[i] = Mack$FullTriangle[i,J]^2 * (Gamma_exact[i] + Delta[i])
    msep_reel_exact[i] = Mack$FullTriangle[i,J]^2 * (Phi_exact[i] + Delta[i])
}

#Covariance
for (i in 2:(I-1)){ 
    for (k in (i+1):I){
        cov_obs_exact[i,k] <- Mack$FullTriangle[i,J]*Mack$FullTriangle[k,J]*(Upsilon_exact[i] + Lambda[i])
        cov_reel_exact[i,k] <- Mack$FullTriangle[i,J]*Mack$FullTriangle[k,J]*(Psi_exact[i] + Lambda[i])    
    }
}

#MSEP ag.
msep_obs_exact[I+1] <- sum(msep_obs_exact[2:I]) + 2*sum(cov_obs_exact)
msep_reel_exact[I+1] <- sum(msep_reel_exact[2:I]) + 2*sum(cov_reel_exact)

##################################################################################################################
##Outputs Fonction 
##################################################################################################################
msep_Mack <- array(0,c(1,I+1))
msep_Mack[1:I] <- Mack$Mack.S.E[,I]
msep_Mack[I+1] <- Mack$Total.Mack.S.E
Vari <- array(0,c(1,I+1))
for (i in 1:I){
    Vari[i] <- Mack$FullTriangle[i,J]^2 * Psi[i]
}
Vari[I+1] <- sum(Vari[1:I])

# complete vision
if(bool==1){    
    result <- cbind(t(msep_Mack), t(sqrt(Vari)), t(sqrt(msep_obs)), t(sqrt(msep_obs_exact)), t(sqrt(msep_reel)), t(sqrt(msep_reel_exact)))
    result <- as.data.frame(result)
    l <- list("Variance CDR", "MSEP Mack","MSEP obs. approx.","MSEP obs. exact", "MSEP real approx.","MSEP real exact")
    names(result) <- l
    return(result)
    }

#si l'on souhaite uniquement l'aspect solvabilité
if(bool==0){
    result <- cbind(t(msep_Mack), t(sqrt(msep_obs)), t(sqrt(msep_obs_exact)))
    result <- as.data.frame(result)
    l <- list("MSEP Mack","MSEP obs. approx","MSEP obs. exact")
    names(result) <- l
    return(result)
    }
}

Les méchanismes bonus-malus

L’outil de base pour l’étude des systèmes bonus-malus est l’utilisation des chaînes de Markov, dès lors que l’on travaille sur des classes de bonus. Pour plus de compléments sur les chaînes de Markov, je renvoie à des vieilles notes de cours (inachevée, et sûrement remplis de coquillettes) ici.

  • Les chaînes de Markov

L’outil de base pour l’étude des systèmes bonus-malus est l’utilisation des chaînes de Markov, dès lors que l’on travaille sur des classes de bonus. Pour plus de compléments sur les chaînes de Markov, je renvoie à des vieilles notes de cours (inachevée, et sûrement remplis de coquillettes) ici.
Pour reprendre les grandes idées, on considère une suite de variables aléatoires à valeurs dans https://perso.univ-rennes1.fr/arthur.charpentier/latex/calE.png, https://perso.univ-rennes1.fr/arthur.charpentier/latex/suiteXn.png. Cette suite est sera une chaîne de Markov si pour tout https://perso.univ-rennes1.fr/arthur.charpentier/latex/ngeq1.png, et pour tout

https://perso.univ-rennes1.fr/arthur.charpentier/latex/suitei.png

telle que

https://perso.univ-rennes1.fr/arthur.charpentier/latex/markov-proba1.png

on ait

https://perso.univ-rennes1.fr/arthur.charpentier/latex/proab1.png
= https://perso.univ-rennes1.fr/arthur.charpentier/latex/proab2.png

On parlera d’absence de mémoire (au sens où seule la dernière valeure intervient dans le calcul de la loi conditionnelle, et pas l’intégralité de la trajectoire). Si https://perso.univ-rennes1.fr/arthur.charpentier/latex/suiteXn.png est une chaîne de Markov de distribution initiale https://perso.univ-rennes1.fr/arthur.charpentier/latex/lambdaV.png et de matrice de transition P, pour tout https://perso.univ-rennes1.fr/arthur.charpentier/latex/nkentier.png,

https://perso.univ-rennes1.fr/arthur.charpentier/latex/proba-markov-statioonnaire.png
https://perso.univ-rennes1.fr/arthur.charpentier/latex/pr135.png

Pour étudier l’évolution sur le long terme, on peut utiliser l’équation de Chapman-Kolmogorov: si https://perso.univ-rennes1.fr/arthur.charpentier/latex/suiteXn.png  est une chaîne de Markov de distribution initiale  https://perso.univ-rennes1.fr/arthur.charpentier/latex/lambdaV.png et de matrice de transition P,  pour tout
https://perso.univ-rennes1.fr/arthur.charpentier/latex/nkentier.png, pour tout https://perso.univ-rennes1.fr/arthur.charpentier/latex/zvze54v6.png,

https://perso.univ-rennes1.fr/arthur.charpentier/latex/chapman-kol.png

Un corrolaire important est que si

https://perso.univ-rennes1.fr/arthur.charpentier/latex/ezjpjio.png

alors

https://perso.univ-rennes1.fr/arthur.charpentier/latex/alcnaoli.png

Bref, on peut facilement étudier l’évolution à long terme de la chaîne de Markov en faisant un peu de calcul matriciel, i.e. des calculs de puissances de matrices.
En fait, formellement, on peut soit étudier des matrices, soit des graphes. Par exemple,

Si on veut étudier l’évolution sur deux périodes, il suffit de calculer le carré de la matrice, i.e.

  • Les systèmes bonus-malus par classes

Les tableaux ci-dessous sont tirés des annexes de la bible de Jean Lemaire, Bonus-Malus systems in automobile insurance. Le plus simple est probablement celui de Hong Kong, qui est Markovien (à l’ordre 1),

https://blogperso.univ-rennes1.fr/arthur.charpentier/public/perso2/.BM-hong-kong-B-9_m.jpg

L’exemple ci-dessous est celui qui a existé en Allemagne dans les années 80,

https://blogperso.univ-rennes1.fr/arthur.charpentier/public/perso2/.BM-germany-B-7_m.jpg

puis celui qui existait en Finlande,

https://blogperso.univ-rennes1.fr/arthur.charpentier/public/perso2/.BM-finnish-B-6_m.jpg

et enfin le système italien du début des années 90,

https://blogperso.univ-rennes1.fr/arthur.charpentier/public/perso2/.BM-italy-B-11_m.jpg

Pour tous ces sytèmes, une fois spécifiée la loi du nombre de sinistres par an, il est possible de construire la matrice de transition des classes de bonus sur un an, ou plus.
En particulier, l’étude de la loi invariante de ce genre de sytème (que je mentionnais en introduction) permet d’étudier la viabilité du système à plus ou moins long terme, et la vitesse de la convergence vers l’état stationnnaire.

  • Le système bonus-malus français

En France, le mécanisme de bonus-malus peut paraître un peu plus compliqué, et il fonctionne de la manière suivante.

Le calcul du bonus-malus s’ef
fectue en fonction des sinistres survenus (ou non) pendant la période de douze mois consécutifs précédant de deux mois l’échéance annuelle du contrat d’assurance.
Seuls les accidents responsables donnent lieu à malus. C’est à dire qu’en cas de foce majeur, ou si l’accident est entièrement imputable à un tiers, alors le conducteur n’aura pas de malus. Notons aussi que Le malus est rattaché au contrat d’assurance. Autrement dit, si quelqu’un cause un accident avec un véhicule assuré par un tiers, il se supportera pas le malus.
La prime de référence est la prime qui sera utilisée pour ensuite actualiser le montant de prime en tenant compte du coefficient de majoration ou minoration. Il s’agit de la prime technique calculée par l’assureur en fonction de charactéristiques a priori, comme la puissance du véhicule véhicule, la zone géographique de circulation, l’usage socioprofessionnel ou encore le kilométrage parcouru….
A partir de cette prime de référence, on calcule une prime qui tient compte de l’expérience passée, via un coefficient multiplicatif.
La valeur de référence est de 1 pour le taux de majoration-minoration. Ce taux subie une réduction – corrrespondant à un bonus – pour une année sans accident, d’un niveau de 5 %. Toutefois, la valeur du coefficient est arrêté à la deuxième décimale et arrondi par défaut. Au final, le taux prend donc un nombre fini de valeurs. Pour résumer schématiquement

1 année sans accident 100 x 0.95 = 0.95 de Bonus
2 années sans accident 0.95 x 0.95 = 0.90 de Bonus
3 années sans accident 0.90 x 0.95 = 0.85 de Bonus
4 années sans accident 0.85 x 0.95 = 0.80 de Bonus
5 années sans accident 0.80 x 0.95 = 0.76 de Bonus
6 années sans accident 0.76 x 0.95 = 0.72 de Bonus
7 années sans accident 0.72 x 0.95 = 0.68 de Bonus
8 années sans accident 0.68 x 0.95 = 0.64 de Bonus
9 années sans accident 0.64 x 0.95 = 0.60 de Bonus
10 années sans accident 0.60 x 0.95 = 0.57 de Bonus
11 années sans accident 0.57 x 0.95 = 0.54 de Bonus
12 années sans accident 0.54 x 0.95 = 0.51 de Bonus
13 années sans accident 0.51 x 0.95 = 0.48 plafonné à 0.50

La réduction maximale est donc de 0.50, soit 50% de réduction sur la prime de référence.
En revanche, chaque accident responsable donne lieu à une majoration – on palera de malus – de 25 % du coefficient par accident. Toutefois, cette majoration est réduite de moitié (12,5 %) si la responsabilité n’est que partiellement engagée.Là encore, le coefficient obtenu est arrêté à la deuxième décimale et arrondi par défaut. Aussi, avec deux accidents responsables dans l’année, le coefficient sera multiplié par un coefficient de 1,56. Afin d’éviter l’absence d’assurance, on notera qu’en aucun cas, le coefficient de réduction-majoration ne peut excéder 3,50.
Pour les jeunes conducteurs, l’article A. 335-9-1 du code des assurances prévoit que l’assureur a la possibilité de majorer la prime de référence (dans des limites fixées par la loi là encore).
Notons enfin qu’il existe des avantages pour les bons conducteurs: si un conducteur est au bonus maximal depuis au moins trois ans, il conserve son bonus après le premier accident responsable. Notons que c’est le genre de disposition qui figure dans la loi et que certains assureurs utilisent à des fins commerciales (comme le montre l’affiche ci-dessus ou ci-dessous, qui a été analysée sur le blog de Gizmo (ici)).

Il existe aussi une règle dite de descente rapide qui peut atténuer les effets des malus: après deux années consécutives sans accident, votre coefficient ne peut être supérieur à 1.
Ce qui rend le système intéressant, c’est qu’il est obligatoire, en France. Aussi, en cas de changement d’assureur, le coefficient de bonus-malus sera recalculé à partir d’un relevé d’information: l’assureur est tenu de fournir ce document à la résiliation du contrat (contenant des informations sur le conducteur et le véhicule, ainsi que le nombre des sinistres survenus au cours des cinq années précédant l’établissement du relevé d’information, ainsi que leur nature et la part de responsabilité retenue)
Davantage d’information se trouve ici, (ce document propose d’ailleurs une petite simulation à la fin sur le coût d’un accident) ou encore (pour une approche plus actuarielle), voire (pour une approche plus économique).

De l’interprétation d’un effet nonlinéaire en régression

Je vais poursuivre un peu le TD de statistique de l’actuariat 2, où nous avions parlé tarification et régression Poissonnienne.

  • Analyser un effet nonlinéaire

Dans un premier temps, nous avions regardé un modèle simple, où l’on voulait voir si l’âge du conducteur expliquait la fréquence de sinistres. Bref, on peut utiliser une régression log-Poisson, qui donne la prédiction suivante, en fonction de l’âge,

Mais l’âge est une variable un peu particulière car c’est une fausse variable continue: seule les valeurs entières apparaissent. Comme on chercher à calculer

https://perso.univ-rennes1.fr/arthur.charpentier/latex/latex-glm-nb-1.png

notons qu’il est possible de considérer l’estimateur empirique naturel de cette grandeur, i.e.

https://perso.univ-rennes1.fr/arthur.charpentier/latex/latex-glm-nb-2.png

https://perso.univ-rennes1.fr/arthur.charpentier/latex/latex-glm-nb-3.png

Graphiquement, on obtient la courbe suivante,

Le côté ératique de la série à droite et à à gauche vient de la faible représentation dans la base, comme en atteste l’histogramme de l’âge du conducteur,

Je peux en profiter pour revenir à un commentaire posté ici sur les modèles GAM. En effet, au lieu de regarder la fréquence en fonction de l’âge, on peut s’intéresser à la la différence par rapport à la fréquence globale.

> m=mean(base$nbsin)
> Y=predict(REG,newdata=data.frame(age=X,expo=E),type = "response")
> lines(X,Y/m-1,col="red")

On obtient alors le graphique suivant (en superposant à l’ajustement obtenu par GAM)

Notons que si l’on ajuste un modèle log-Poisson uniquement sur la tranche d’âge [50-70], on est très proche de ce que propose le modèle GAM, d’où l’interprétation de l’ajustement GAM comme un ajustement local.

 

  • Changer de modèle (et les comparer)

On a choisi un modèle log-Poisson, mais on se doute que d’autres modèles seraient possible. En particulier dans notre base, les assurés ont eu 0 ou 1 sinistre l’année passée. Autrement dit, des lois binomiales peuvent être ajustée,

> REG0 = glm(nbsin~age,base,
+  family=poisson(link="log"),
+  offset=expo)
> REG1a = glm(nbsin~age,base,
+  family=binomial(link = "logit"),
+  offset=expo)
> REG1b = glm(nbsin~age,base,
+  family=binomial(link = "probit"),
+  offset=expo)
> REG1c = glm(nbsin~age,base, 
+  family=gaussian(link = "identity"),
+  offset=expo)

Le critère AIC est souvent invoqué pour juger de l’adéquation d’un modèle. Il vaut ici

> AIC(REG0)
[1] 36470.7
> AIC(REG1a)
[1] 36044.45
> AIC(REG1b)
[1] 37304.02
> AIC(REG1c)
[1] 32549.83

Mais on peut également regarder le critère BIC,

> AIC(REG0, k = log(nrow(base)))
[1] 36488.81
> AIC(REG1a, k = log(nrow(base)))
[1] 36062.57
> AIC(REG1b, k = log(nrow(base)))
[1] 37322.14
> AIC(REG1c, k = log(nrow(base)))
[1] 32577

Autrement dit, en utilisant le règle “the smaller, the better” on retiendrait le dernier modèle qui est un modèle Gaussien, manifestement assez mauvais: sur le graphique ci-dessous, la régression log-Poisson est en mauve, presque confondu avec la régression logistique1, en bleu, la régression probit en vert et le modèle linéaire (Gaussien) en rouge,

Bref, ce dernier modèle n’est pas convainquant… Par contre, notons que si l’on compare le modèle log-Poisson et logit, ce dernier conduit à tarifer moins cher les âges extrêmes,

Autrement dit, au lieu de comparer des modèles actuariels sur la base d’outils statistiques, et uniquement de ces outils, il peut être instructif d’essayer de comprendre qui verra sa prime baisser avec tel ou tel modèle… Mais je reviendrais plus longuement, une autre fois, sur les histoires de choix de modèles…..

Crédibilité et calcul bayésien

On a commencé en cours la tarification a posteriori, en commençant par la tarification bayésienne, et en présentant ensuite l’approche d’Hans Bühlmann comme un cas particulier en projetant dans un sous espace de combinaisons linéaires (un peu comme l’économétrie linéaire était présenté en cours comme un cas particulier de calcul d’espérance conditionnelle).La classification a priori sera poursuivie un TD, sur des bases de données, afin d’orienter davantage le cours sur la pratique des GLM.
L’utilisation des techniques bayésiennes en actuariat ne s’est pas faite sans peine. Pour ceux qui souhaitent s’en convaincre, il suffit de relire les comptes rendus de la conférence à Trieste en 1962. Dans le discours de Bruno de Finetti (publié dans l’ASTIN Bulletin), on retrouve un point historique, présentant les nouveaux travaux d’Hans Bühlmann,

Et Bruno de Finetti de conclure en rappelant les apports d’Arthur Bailey en 1950,

  • La formule de Bayes, ou apprendre par l’expérience

La formule de Bayes ne dit pas grand chose quand on le regarde sous la forme suivante

https://perso.univ-rennes1.fr/arthur.charpentier/latex/latex-bayes-01.png

A la rigueur on peut aussi la réécrire en utilisant deux fois cette formule

https://perso.univ-rennes1.fr/arthur.charpentier/latex/latex-bayes-02.png

Bon, là encore, on ne voit peut être toujours pas ce qu’on peut en faire…. L’idée géniale est de faire la lecture suivante

https://perso.univ-rennes1.fr/arthur.charpentier/latex/latex-bayes-03.gif

Autrement dit la loi a posteriori du paramètre, conditionnellement aux observations (en bleu) est proportionnelle au produit de  la loi a priori du paramètre (en vert) et de la vraisemblance des observations (en rouge, qui sera un simple produit de lois si on suppose l’indépendance conditionnelle au paramètre). Le paramètre en noir est juste un facteur de normalisation, histoire que la loi soit effectivement une mesure de probabilité. D’où parfois l’écriture

https://perso.univ-rennes1.fr/arthur.charpentier/latex/latex-bayes-4.gif

Bref, on a réussit à inverser les deux termes dans la probability (on retrouve le terme “inverse probability” dans la littérature anglo-saxonne).
En cours, on avait vu que si on supposait le nombre annuel de sinistres Poisson, conditionnellement au paramètre, avec une indépendance (toujours conditionnellement à ce paramètre), et que l’on suppose que la loi a priori du paramètre est une loi Gamma, alors la loi a posteriori se calcule facilement et c’est une loi Gamma.
Considérons la petite simulation suivante: on suppose que le nombre de sinistres suit une loi de Poisson, et on cherche à estimer le paramètre de cette loi. On commence par une loi a priori Gamma, centrée sur 5, et on regarde évoluer la loi a posteriori au fur et à mesure que l’on observe les sinistres arriver (ou les années passer en l’occurence).

On peut également considerer l’évolution de l’espérance a posteriori du paramètre, ainsi qu’un intervalle de confiance (puisqu’on connaît la loi)

La “magie” des méthodes bayésiennes est que peu importe la loi a priori, on va converger vers la vraie valeur. Par exemple sur l’animation ci-dessous, on part d’une loi a priori de moyenne 2,

avec là aussi l’évolution de la moyenne a posteriori avec un intervalle de confiance,

La morale de l’histoire est que comme la loi a priori n’influence pas (j’exagère un peu car elle va influencer la vitesse de la convergence, ce qui n’est pas rien en pratique), autant prendre quelque chose de simple pour les calculs: la loi conjuguée de la loi exponentielle. La loi conjuguée d’une loi de la famille exponentielle de la forme

s’écrit

où et sont des constantes, et où  est une constante de normalisation (afin que cette fonction soit effectivement une mesure de probabilité).

  • L’hérésie de la statistique bayésienne

La principale hérésie des statisticiens bayésien a été de proposer une lecture des probabilités d’un point de vue subjectif. En particulier Jimmie Savage et Bruno de Finetti partait de l’hypothèse que chacun avait à sa disposition son propre espace de probabilité.C’est ce qui lui permet de se construire un a priori. Et au fur et à mesure que le temps passe, on apprend, et on peut dériver une loi a posteriori. Par exemple  je pars d’un a priori comme quoi un enfant sur deux qui nait est un garçon, en France. Je peux ensuite récupérer des données, et affiner ma connaissance. C’est ce qu’avait fait Laplace il y a fort fort longtemps….
L’hérésie est ici un
problème philosophique lié au concept même de probabilité. Mais on peut noter qu’à l’époque dans les années 60, les débats étaient assez virulants. Et pour Jimmie Savage et Bruno de Finetti, il n’est pas nécessaire d’avoir le formalisme mathématique d’Andreï Kolmogorov pour appréhender la notion de probabilité (c’est à dire la mesure d’ensembles dans un espace abstrait d’évènements, voir ici ou ). On notera que pour John Keynes (qui avait écrit un “traité de probabilités“) il existait une probabilité subjective, mais davantage liée à une information donnée par l’éducation. C’est un problème culturel plus que personnel. Mais finalement, il a surtout fallu attendre les travaux de Daniel Kahneman et Amos Tversky dans les années 70 pour que des psychologues se penchent sur la construction individuelle des probabilités.

  • La crédibilité au sens de Bühlmann

L’idée d’Hans Bühlmann est de proposer une formule fermée qui approche l’espérance conditionnelle, tout simplement en utilisant celle qui avait été obtenue dans le cas bayésien. En tarification, on s’intéresse à des primes pures, et il est finalement assez naturel de chercher à calculer

Hans Bühlmann a suggéré d’approcher cette grandeur par

avec classiquement

est appelé facteur de crédibilité, avec , et où respectivement,


Bref, on a un facteur de crédibilité qui dépend de la variance intra et de la variance inter, pour reprendre des termes connus en statistique. L’énorme avantage est que cette formule peut très facilement s’utiliser sur des données. En particulier, considérons un tableau de la forme suivante

où on regarde si un contrat a – ou n’a pas – eu de sinistre chacune des années. Formellement, on dispose d’un tableau d’observations,

(on notera que le nombre d’années d’observations pour chacun des contrats peut être différent). On peut alors obtenir les estimateurs empiriques des moments nécéssaires, i.e. on pose

correspondant à la moyenne empirique pour un contrat,

et

Le facteur de crédibilité pour le contrat i est alors

où tout naturellement

Formation ERM de l’Institut des Actuaires

Intervention jeudi après midi dans le cadre de la formation ERM (Enterprise Risk Management) de l’Institut des Actuaires, sur le thème “everything you wanted to ask about tail dependence in risk management and you’re afraid to ask“.
Comme Anne-Laure Fougères est déjà intervenue sur les extrêmes, et que Stéphane Loisel ont présenté les copules, le cours sera une discussion autour des diverses recommandations du Groupe Consultatif, des documents du CEIOPS, et de divers documents qui circulent sur le calcul du SCR (Solvency Capital Requierment) de l’agrégation des risques. Ou “everything you wanted to ask about tail dependence in risk management and you’re afraid to ask“… Les slides sont en ligne ici, et ne contienent que les éléments illustratifs du cours.

Les généralités sont généralement fausses

Suite à un commentaire sur un vieux billet (ici), concomitante avec une demande d’intervention pour la formation ERM proposée par l’Institut des Actuaires, je me suis replongé dans les calculs de SCR (solvency capital requierement). Plus je lis la documentation technique, plus je suis surpris par l’importance de l’hypothèse de normalité ! J’ai repris mon Loss Model, qui est le livre de base pour l’examen d’actuariat non-vie dans les pays anglo-saxons, et j’ai vu des dizaines de lois bizarres… mais rien sur la loi normale…

  • de l’hypothèse de normalité dans la formule d’agrégation

Comme je le disais dans un vieux billet (ici), on ne peut montrer la validité de la formule de calcul du SCR que dans le cas Gaussien (ou plus précisément dans le cas elliptique, voire ici). En effet, la preuve est assez simple. Je renvois vers un papier tout chaud de Laurent Devineau et Stéphane Loisel sur le sujet, ici.

  • où l’on continue à entrevoir l’hypothèse de normalité

En fait l’hypothèse de normalité est vraiment partout dans les calculs de SCR. Mais bien sûr, sans le dire. Par exemple dans un papier de la SCOR,

https://blogperso.univ-rennes1.fr/arthur.charpentier/public/perso2/.doc14-1_m.jpg

Comme je le disais l’autre jour (ici), le problème des documents “professionels” (par opposition à “académiques“), c’est que comme on veut que le CEO le lise, on évite de parler d’hypothèses, ou pire, de mettre un peu de maths…. Bref, on apprend que “generally” (terme d’une grande rigueur scientifique, que je recaserai un jour dans mes démonstrations), “the expected shortfall[…] at the 99 % level corresponds quite closely to the[…] value-at-risk at a 99.6% level “. Mais on retrouve la même phrase dans un rapport de l’edhec

https://blogperso.univ-rennes1.fr/arthur.charpentier/public/perso2/.SCR-edhec_m.jpg

Là aussi, aucune source… Damned ! En creusant un peu, Actuaris m’a mis sur la piste (ici), car on apprend que c’est l'”autorité suisse” qui a énoncé ce résultat que je qualifierais de surprenant1.

https://blogperso.univ-rennes1.fr/arthur.charpentier/public/perso2/.actuaris-var_m.jpg

Effectivement, on peut retrouver le document originel,

https://blogperso.univ-rennes1.fr/arthur.charpentier/public/perso2/.doc13-1_m.jpg

J’avoue avoir essayé – sans grand succès – de lire l’annexe D. Je n’ai pas tout compris, mais j’y ai lu un “assuming that the distribution functions follows a normal law”,

https://blogperso.univ-rennes1.fr/arthur.charpentier/public/perso2/.doc13-2_m.jpg

Donc si j’ai bien tout compris, “generally”, “the distribution functions follows a normal law” ! Bon, je crois que je vais pouvoir envoyer un mail à Stuart, Harry et Gordon pour leur expliquer qu’ils peuvent virer des pages à leur livre ! Ou peut-être Jean-Luc Besson (pour lequel j’ai une immense estime !) qui présente beaucoup de lois utiles en assurance dans l’ouvrage qu’il a coécrit avec Christian Partrat…
Mais peut-être a-t-il raison…? On peut se dire qu’en actuariat, on avait l’habitude de présenter l’approximation “normal power“, et peut être que ce résultat est vrai pour d’autres lois….
Regardons vite-fait sur quelques lois….. Par exemple, pour la loi lognormale, si on regarde en fonction de l’écart-type de la loi normale sous-jacente,

Effectivement, on colle au quantile à 99.6%. L’ordonnée est en échelle logarithmique pour que le dessin soit présentable, avec en pointillé la VaR (formule exacte, à 99.6% en rouge, et à 99.8% en bleu), et les ronds sont des estimations de la TVaR à 99%. Par contre, par exemple pour la loi lognormale, si on regarde en fonction de l’écart-type de la loi normale sous-jacente, plus la volatilité est grande, plus l’approximation est mauvaise.

Bon, essayons une loi Gamma (classique en assurance, car c’est une loi de la famille exponentielle, et souvent comparée à la loi lognormale).

Visiblement, là aussi l’approximation semble marcher….  mais rappelons que la loi Gamma est une loi à queue fine… Allons-y franchement, comparons avec le cas d’une loi de Pareto…

cette fois, l’abscisse se lit dans l’autre sens: à gauche, c’est le cas le plus risqué (avec une variance infinie entre 1 et 2). On notera que pour les risques les plus importants, là aussi on sous-estime la TVaR en considérant la VaR (à un seuil plus important). On retrouve d’ailleurs exactement le même type de graphique avec une loi de Student,

Bref, on me dira que je pinaille, n’empêche qu’à mon avis, les compangies sous-estiment “généralement” la TVaR (à supposer que ça soit la mesure de risque qui les intéresse) en prenant une hypothèse de normalité…

1 en tant qu’enseignant, ce genre de conclusion me gêne toujours, car je me fatigue à faire des cours sur la mesure de risque, j’essaye de présenter des notions aussi proprement que possible, en expliquant que les hypothèses c’est important, qu’il faut définir différentes notions de mesures de risques (entre les convexes, les cohérentes, les dynamiques, etc), et là tout d’un coup, on nous apprend que “generally” c’est kif kif….

Grippe, cravate et (ré)assureurs

Entendu sur France Info, dans la revue de presse de Jean Christophe Martin, ce matin, ici, “parmi ces conseils, une bonne nouvelle pour les allergiques à la cravate : en plus des recommandations classiques, se laver les mains et se dire bonjour sans bisous, le Parisien signale que le groupe Scor, 5ème société mondiale de réassurance, demande aussi à ses salariés d’éviter de porter une cravate. Une consigne qui a eu un gros succès dans les bureaux de la Défense, au siège du groupe. Et oui, contrairement aux chemises, on ne fait pas souvent nettoyer sa cravate, qui peut se transformer en boîte à postillons et en bouillon de culture… Le nouveau mot d’ordre dans les bureaux : à bas la cravate, vive la grippe A...”. Chez AXA, on apprend que c’est la bise, qui est interdite (ici).