Après un petit clin d’œil à Michel Rabagliati qui a gagné le prix du public hier à Angoulême.... Je vais poursuivre le billet commencé la semaine dernière (ici), histoire de montrer que ce n’est pas si compliqué que ça de calculer des PM avec R. Compte tenu des doubles indices dans toutes les fonctions, nous allons créer des matrices afin de stocker les différentes grandeurs.
> L=read.table("http://perso.univ-rennes1.fr/arthur.charpentier/TD8890.csv",
+ sep=";",header=TRUE)
> Lx=L$Lx
> m=length(Lx)
> p=matrix(0,m,m); d=p
> for(i in 1:m){
+ p[1:(m-i+1),i]=Lx[1+i:m]/Lx[i]
+ d[1:(m-i+1),i]=(Lx[i:(m)]-Lx[i:(m)+1])/Lx[i]}
> diag(d[m:1,])=0
> diag(p[m:1,])=0
> q=1-p
aussi, p[10,40] correspondra à http://perso.univ-rennes1.fr/arthur.charpentier/latex/PM80.png (attention, cf commentaire plus bas). Je stocke ainsi la probabilité de décéder au bout de k années, ou bien exactement dansk années. On peut alors calculer les différentes grandeurs,
> a=E=A=matrix(0,m,m)
> r=3/100
> for(i in 1:(m-1)){
+ a[,i]=cumsum(1/(1+r)^(0:(m-1))*c(1,p[1:(m-1),i]))
+ }
> for(i in 1:(m-1)){
+ #for(j in 1:(m-i)){
+ #A[j,i]=sum(1/(1+r)^(1:j)*p[1:j,i]*q[1,i-1+1:j])
+ A[,i]=cumsum(1/(1+r)^(1:m)*d[,i])
+ }
> for(i in 1:m){
+ E[,i]=(1/(1+r)^(1:m))*p[,i]
+ }

A partir de ces constantes (ou de ces matrices de constantes), on peut déjà calculer la prime annuelle du contrat décès que nous avions considéré,
> x=50; n=30
> prime=A[n,x]/a[n,x]
> sum(prime/(1+r)^(0:(n-1))*c(1,p[1:(n-1),x]))
[1] 0.3116454
> sum(1/(1+r)^(1:n)*d[1:n,x])
[1] 0.3116454

De plus, on peut calculer les provisions mathématiques à l’aide des différentes méthodes (et vérifier qu’elles donnent bien le même résultat).
> ### méthode rétrospective
> VR=(prime*a[1:n,x]-A[1:n,x])/E[1:n,x]
> plot(0:n,c(0,VR),col="red")
> ### méthode prospective
> VP=diag(A[n-(0:(n-1)),x+(0:(n-1))])-prime*diag(a[n-(0:(n-1)),x+(0:(n-1))])
> points(0:n,c(VP,0),pch=4,col="blue")
> # méthode itérative
> VI=0
> for(k in 1:n){
+ VI=c(VI,(VI[k]+prime-A[1,x+k-1])/E[1,x+k-1])
+ }
> points(0:n,VI,pch=5,col="green")
Pour un contrat signé par une personne de 50 ans, pour une couverture de 30 ans, au taux d’actualisation de 3%, on obtient les provisions mathématiques suivantes,

(les trois vecteurs coïncident, presque par miracle). Si la personne avait 45 ans, les provisions auraient l’allure suivante (en bleu),

Enfin, si on regarde la sensibilité au taux d’actualisation, un taux de 5% donnerait les provisions suivantes (en vert),

Ben voilà. La prochaine étape ça sera de distribuer les devoirs maison....