Reg log - prédictions en validation croisée incohérentes
Bonjour,
_____
J'ai déjà posté ce sujet sur le forum CIRAD, sans succès jusque là. J'espère que vous pourrez m'aider. Si j'ai une réponse sur un des deux sujets je ferai part de la réponse dans tous les cas.
Merci par avance pour votre aide.
_____
Je suis en train de travailler sur une régression logistique et j'ai détecté une anomalie lorsque je fais une validation croisée. Je poste donc ici dans l'espoir que vous m'aidiez à cerner le problème.
J'utilise 100 observations de deux variables générées aléatoirement selon des lois normales de paramètres différents.
Par la suite je crée un modèle qui a l'air assez cohérent (voir la courbe ROC à la fin).
Code:
1
2
3
4
5
6
7
8
9
| #### Génération des données ####
z = c(rep(0,50),rep(1,50))
x = c(rnorm(50,0,4),rnorm(50,2,4))
y = c(rnorm(50,0,1),rnorm(50,2,1))
plot(x,y,col=z+1)
#### Création du modèle et prédiction par GLM ###
res.glm = glm(formula=z~x+y,family=binomial)
predict.res.glm = predict(res.glm, type="response") |
J'effectue maintenant la validation croisée sur quatre groupes. Le principe étant de prédire les données de chaque groupe à partir des trois autres.
Code:
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
| #### Validation croisée ####
# Définition de 4 groupes. Trois d'entres eux construisent le modèle, on prédit sur le dernier.
g = sample(rep(1:4,25))
# Initialisation des prédictions
vect.proba = rep(0,100)
for (k in 1:4){
# Echantillon d'apprentissage sur 3 groupes
appr = as.data.frame(cbind(x[g!=k], y[g!=k]))
z.appr = z[g!=k]
res.glm.kfold = glm(formula=z.appr~appr[,1]+appr[,2],family=binomial)
# Echantillon de vérification
verif = as.data.frame(cbind(x[g==k], y[g==k]))
vect.proba[g==k] = predict(res.glm.kfold, newdata = verif, type="response") #https://stat.ethz.ch/R-manual/R-patched/library/stats/html/predict.lm.html
} |
Alors ici, plusieurs choses :
--> J'obtiens des warnings (8 exactement) au niveau de la prédiction dans la boucle.
Citation:
1: 'newdata' had 25 rows but variables found have 75 rows
2: In vect.proba[g == k] = predict(res.glm.kfold, newdata = verif, :
number of items to replace is not a multiple of replacement length
J'avoue ne pas comprendre car la définition de newdata est, selon le manuel d'utilisation :
Citation:
optionally, a data frame in which to look for variables with which to predict. If omitted, the fitted linear predictors are used.
Je fournis bien le même nombre de variables, mais en effet il n'y a pas le même nombre d'observations ce qui est tout à fait normal. J'ai vu plusieurs scripts sur le net qui utilisent ce paramètre pour faire une validation croisée (ici par exemple).
D'où la question : Est-ce que c'est la bonne chose à faire ?
--> Ensuite : mes prédictions par validation croisée sont foireuses. Je vous indique un extrait des prédictions entre la 40e et la 60e ligne (à la 50e, on passe de 0 à 1)
Code:
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
| Nb Z Predic_GLM Predic_ValCroisee
40 0 0.319408438 0.145276517
41 0 0.205839858 0.119031341
42 0 0.012653655 0.029117270
43 0 0.181420164 0.749168854
44 0 0.020543180 0.914641430
45 0 0.442173491 0.645873783
46 0 0.157567167 0.131166378
47 0 0.045872366 0.102201020
48 0 0.085971756 0.060728705
49 0 0.489664348 0.021215917
50 0 0.013036409 0.036993982
51 1 0.984022606 0.007688740
52 1 0.918224041 0.039215022
53 1 0.944663266 0.201269337
54 1 0.297424130 0.701688379
55 1 0.803391279 0.897638029
56 1 0.946804556 0.210754072
57 1 0.260940236 0.097831911
58 1 0.996113114 0.130593783
59 1 0.919221343 0.071754892
60 1 0.508894219 0.071091520 |
Vous pouvez voir que dans certains cas les prédictions ne concordent pas du tout. J'ai testé les différents modèles en modifiant les valeurs de k à la main pour tester. Les paramètres ne sont pas excessivement différents, et du coup je ne comprends pas pourquoi il y a un tel écart.
Code:
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
| > res.glm
Coefficients:
(Intercept) x y
-2.0162 0.1303 2.1463
> glm(formula=z.appr.test.1~appr.test.1[,1]+appr.test.1[,2],family=binomial)
Coefficients:
(Intercept) appr.test.1[, 1] appr.test.1[, 2]
-2.121 0.224 2.437
> glm(formula=z.appr.test.2~appr.test.2[,1]+appr.test.2[,2],family=binomial)
Coefficients:
(Intercept) appr.test.2[, 1] appr.test.2[, 2]
-2.3885 0.1578 2.3182
> glm(formula=z.appr.test.3~appr.test.3[,1]+appr.test.3[,2],family=binomial)
Coefficients:
(Intercept) appr.test.3[, 1] appr.test.3[, 2]
-2.02831 0.09401 2.26203
> glm(formula=z.appr.test.4~appr.test.4[,1]+appr.test.4[,2],family=binomial)
Coefficients:
(Intercept) appr.test.4[, 1] appr.test.4[, 2]
-1.73949 0.09701 1.78679 |
Je finis en indiquant les courbes ROC, où vous pouvez bien voir que les prédictions obtenues par ma validation croisée sont très mauvaises.
http://img15.hostingpics.net/pics/245025Rplot.png
Le code pour info
Code:
1
2
3
4
5
6
7
8
9
10
| library("ROCR")
par(mfrow=c(1,2))
roc.pred.glm <- prediction(predict.res.glm,z)
roc.perf.glm <- performance(roc.pred.glm,"tpr","fpr")
plot(roc.perf.glm,colorize=T,lwd=3,xlab = "Taux de faux positifs",ylab = "Taux de vrais positifs",main = "Courbe ROC - Prédiction par GLM")
roc.pred.croisee <- prediction(vect.proba,z)
roc.perf.croisee <- performance(roc.pred.croisee,"tpr","fpr")
plot(roc.perf.croisee,colorize=T,lwd=3,xlab = "Taux de faux positifs",ylab = "Taux de vrais positifs",main = "Courbe ROC - Validation croisée") |
Est-ce que vous avez une idée du problème pour faire ma validation croisée de manière cohérente ?
J'espère que le pavé ne posera pas de souci, mais vu que je m'arrache les cheveux depuis quelques heures sur cela, je tenais à être le plus précis possible.
Merci par avance pour votre aide.
Cordialement,
Emmanuel Rambeau
