6.3 Qualité de prédiction

Pour rester dans le compromis entre le biais et la variance, et éviter le surajustement des données par un modèle, une autre approche intéressante est de comparer la qualité de prédiction des modèles lorsqu’on les applique sur de nouvelles données.

Pour reprendre l’exemple des données utilisées depuis le début de ces notes, la base de données dbfictif.Rda a été divisée en deux parties: la partie TRAIN et la partie TEST. En estimant les paramètres en n’utilisant seulement que la partie TRAIN, il nous reste assez de données pour comparer la qualité de prédiction des modèles sur la partie TEST.

6.3.1 Score de prédiction

L’idée du score de prédiction est d’obtenir une valeur numérique pour évaluer la qualité de la prédiction d’un modèle sur de nouvelles données. Par convention, on supposera que l’objectif est de minimiser le score de prédiction. Plus spécifiquement, on supposera que pour évaluer le modèle \(P\), nous prédirons une valeur \(p\) que nous comparerons avec une observation \(x\). Par la suite, une pénalité \(s(p,x)\) sera calculée pour déterminer l’erreur de prédiction.

Pour le modèle \(P\), on obtiendra la score de prédiction du modèle, \(S(P)\), en prenant la moyenne du score de prédiction du modèle pour toutes les observation \(x\) d’un jeu de données (qui n’a jamais été utilisé pour estimer n’importe quels paramètres ou propriétés du modèle \(P\)):

\[S(P) = \frac{1}{n} \sum_{i=1}^n s(p_i, x_i)\]

Théoriquement, une série de propriétés relativement simples doivent être respectées pour utiliser un score de prédiction. Nous ne les couvrirons pas mais pour les détails, on peut consulter:

Czado, C., Gneiting, T., & Held, L. (2009). Predictive model assessment for count data. Biometrics, 65(4), 1254-1261.

L’analyse du nombre de réclamation est particulière: les observations sont des valeurs entières, et la moyenne observée est souvent proche de 0. Ainsi, il faut être prudent dans le choix du score de prédiction.

En supposant que \(\lambda_p\) est la prédiction de l’observation, on peut lister quelques scores pertinents pour les données de comptage pour un modèle \(P\) et une observation \(x\):

Score logarithmique:

\[logs(p, x) = -\log(\Pr(N=x))\]

Score quadratique:

\[quad(p, x) = -2 \Pr(N=x) + \left(\sum_{j=0}^{\infty} \Pr(N=j)^2\right)^2 \] 3) Score d’erreur-carrée

\[sq.err(p,x) = (x - \lambda_p)^2\]

Score sphérique

\[sph(p,x) = - \frac{\Pr(N=x)}{\sum_{j=0}^{\infty} \Pr(N=j)^2}\]

Score Déviance Poisson

\[dev.poi(p, x) = 2\left( x \ln\left(\frac{x}{\lambda_p}\right) + (x - \lambda_p) \right)\]

Le score du modèle \(P\) pour la partie TEST de la base de données se calcule ainsi simplement comme la somme des scores individuels:

\[S(P) = \sum_{i=1}^n S(p_i, x_i)\]

6.3.1.1 Applications en R

Les scores de prédictions sont relativement faciles à calculer.

db.test$pred <- predict(Poisson1, newdata=db.test, type="response") 

logs <- sum(dpois(db.test$nb.sin, db.test$pred, log=TRUE))
sq.err <- sum( (db.test$nb.sin - db.test$pred)^2)
c(logs, sq.err)

## [1] -61111.99  24965.97

Afin de ne pas favoriser un modèle qui surajuste les données, nous avons considéré l’analyse hors-échantillon sur la portion TEST de la base de données. Même si un modèle semble générer des scores de prédiction meilleurs que les autres modèles, il faut considérer la possibilité qu’il ne s’agit que d’une seule analyse fait un seul autre échantillon.

Ainsi, il peut aussi être pertinent de comparer la qualité de prédiction des modèles sur plusieurs jeux de données. Sachant que nous n’avons pas une infinité de jeux de données avec lesquels on peut faire des analyses, la statégie est de diviser la base de données en \(K\) morceaux. Ensuite, pour \(j=1, 2, \ldots, K\), nous effectuerons les tâches suivantes:

On exclut la partie \(j\) du jeu de données;
On estime les paramètres sur les \(K-1\) parties du jeu de données;
On calcule le score de prédiction sur la partie \(j\) du jeu de données.

On somme ainsi \(K\) scores de prédiction, et on peut ainsi voir quel modèle génère le plus petit score total.

En choisissant un certain nombre de folds, il est assez simple de développer une routine pour mettre en place la stratégie d’estimation précédente.

set.seed(334)
NbCV <- 5
folds <- cut(seq(1,nrow(db.fictif)),breaks=NbCV,labels=FALSE)
all.res <- data.frame()
for(i in 1:NbCV){
  testIndexes <- which(folds==i, arr.ind=TRUE)
  test.cv <- db.fictif[testIndexes,]
  train.cv <- db.fictif[-testIndexes,]
  Poisson.cv  <- glm(score.glm, family=poisson(link=log), data=train.cv)
  
  test.cv$pred <- predict(Poisson.cv, newdata=test.cv, type="response") 
  logs <- sum(dpois(test.cv$nb.sin, test.cv$pred, log=TRUE))
  sq.err <- sum( (test.cv$nb.sin - test.cv$pred)^2)
  res <- c(i, logs, sq.err)
  all.res <- rbind(all.res, res)
}
colnames(all.res) <- c('fold', 'logs', 'sq.err')
c(sum(all.res$logs), sum(all.res$sq.err))

## [1] -203475.89   82895.27