Un billet rapide suite à une question que j’ai reçu par mail, au sujet des triangles de liquidation, "pourquoi tout le monde fait du ODP alors que le résultat est le même que Poisson normale ?". Le plus simple est de regarder effectivement, sur un triangle, ce qui se passe si l’on prend des modèles intégrant de la surdispersion, au lieu d’un modèle de Poisson plus classique. Pour commencer, ajustons un modèle de Poisson sur le triangle des incréments,
>  library(ChainLadder)
>  an <- 10; ligne = rep(1:an, each=an); colonne = rep(1:an, an)
>  passe = (ligne + colonne - 1)<=an; n = sum(passe)
>  PAID=GenIns; INC=PAID
>  INC[,2:an]=PAID[,2:an]-PAID[,1:(an-1)]
>  Y = as.vector(INC)
>  lig = as.factor(ligne)
>  col = as.factor(colonne)
>  base = data.frame(Y,col,lig)
>  reg1=glm(Y~col+lig,data=base,family="poisson")
>  sum(exp(predict(reg1,newdata=base))[passe!=TRUE])
[1] 18680856
Si on utilise un modèle de Poisson surdispersé, la prédiction sera la même,
>  reg2=glm(Y~col+lig,data=base,family="quasipoisson")
>  sum(exp(predict(reg2,newdata=base))[passe!=TRUE])
[1] 18680856
En effet, seule la variance change, i.e. l’incertitude sera plus grande si l’on utilise un modèle de Poisson surdispersion. Notons que les tests permettent de valider l’hypothèse de surdispersion,
> dispersiontest(reg1)
        Overdispersion test
data:  reg1
z = 4.3942, p-value = 5.558e-06
alternative hypothesis: true dispersion is greater than 1

Une autre idée peut être de faire une régression binomiale négative (qui permet de prendre en compte la surdispersion),
>  library(MASS)
>  reg3=glm.nb(Y~col+lig,data=base)
> summary(reg3)
(Dispersion parameter for Negative Binomial(13.8349) family taken to be 1)
              Theta:  13.83
          Std. Err.:  2.61
 2 x log-likelihood:  -1460.766

On peut alors calculer à nouveau le montant de provision total,
>  sum(exp(predict(reg3,newdata=base))[passe!=TRUE])
[1] 18085795
qui est légèrement différent du cas Poissonnien.

On peut d’ailleurs aller un peu plus loin, et regarder les mse de prédiction, pour chacun des modèles (que j’avais évoqué ici). Mais auparavant, pour des raisons techniques, il est plus simple de se débarrasser des valeurs nulles dans le triangle,
>  an <- 10; ligne = rep(1:an, each=an); colonne = rep(1:an, an)
>  passe = (ligne + colonne - 1)<=an; np = sum(passe)
>  futur = (ligne + colonne - 1)> an; nf = sum(passe)
>  base$Y2=base$Y; base$Y2[is.na(Y)]=.001
>  reg1b=glm(Y2~lig+col, family=poisson,data=base)
>  reg2b=glm(Y2~lig+col, family=quasipoisson,data=base)
>  reg3b=glm.nb(Y2~lig+col, data=base)

La fonction suivante calcule l’estimateur du montant total de réserve, ainsi que la prediction error (que l’on exprime également en fonction du montant de provisions)
> predCL=function(reg=reg1,regb=reg1b){
+  p = 2*6-1;
+  phi.P = sum(residuals(reg,"pearson")^2)/(np-p)
+  Sig = vcov(reg)
+  X = model.matrix(regb)
+  Cov.eta = X%*%Sig%*%t(X)
+  mu.hat = exp(predict(reg,newdata=data.frame(lig,col)))*futur
+  pe2 = phi.P * sum(mu.hat) + t(mu.hat) %*% Cov.eta %*% mu.hat
+  cat("Total reserve =", sum(mu.hat), "prediction error =", sqrt(pe2),sqrt(pe2)/sum(mu.hat),"\n")
+ }

Avec nos trois modèles, Poisson, ODP et binomiale négative, on obtient,
> predCL(reg1,reg1b)
Total reserve = 18680856 prediction error = 896876.9 0.04801048
> predCL(reg2,reg2b)
Total reserve = 18680856 prediction error = 4736425 0.2535443
> predCL(reg3,reg3b)
Total reserve = 18085795 prediction error = 2058134 0.1137984
Moralité, oui la surdispersion a un impact, ainsi que la méthode que l’on retient pour ma modéliser....