Estoy ejecutando una regresión logística ordinal en R y me encuentro con problemas cuando incluyo variables ficticias. Mi modelo funciona bien con mi primer conjunto de predictores. A continuación, quiero añadir variables ficticias para cada uno de los años representados en mi conjunto de datos.
He creado las variables ficticias con car:recode
de esta manera (una declaración como ésta para cada uno de los 11 años)
fsd$admityear2000 <- recode(fsd$ApplicationYear ,"2000=1;else=0")
El modelo lrm se especifica como sigue
library(Design)
ddist<- datadist(fsd)
options(datadist='ddist')
m4 <- lrm(Outcome ~ relGPA + mcAvgGPA + Interview_Z + WorkHistory_years + GMAT + UGI_Gourman + admityear1999 + admityear2000 + admityear2001 + admityear2002 + admityear2003 + admityear2004 + admityear2005 + admityear2006 + admityear2007 + admityear2008 + admityear2009, data=fsd)
(lo siento por todas las demás variables aleatorias, pero no quiero introducir confusión cambiando mi código)
Me sale el error
singular information matrix in lrm.fit (rank= 22 ). Offending variable(s):
admityear2009 admityear2000 admityear1999
Error in lrm(Outcome ~ relGPA + mcAvgGPA + Interview_Z + WorkHistory_years + :
Unable to fit model using “lrm.fit”
Entiendo que incluir todas las opciones de una variable ficticia sobredefine el modelo, pero me da el error tanto si incluyo los 11 años como si sólo incluyo 10.
He encontrado una sugerencia aquí para establecer el parámetro de penalización de lrm
a un pequeño valor positivo . Si se fija en 1 o 5, el error cambia de forma que sólo nombra una de las variables como infractora. El error no desaparece ni siquiera con penalty=100
.
Soy bastante nuevo en R, pero me encanta la libertad hasta ahora. ¡Gracias por cualquier ayuda!
Respuestas y lecciones
- Los factores son increíbles y no puedo creer que no me haya dado cuenta antes. Hombre eso limpia mi código mucho. Gracias.
- Mi DV, 'Outcome' es efectivamente ordinal y después de hacerlo un factor(), también lo hice ordenado().
- El comando str() también es impresionante y así es como se ven mis datos ahora (con algunas de las variables no relevantes omitidas)
de salida:
str(fsd)
Outcome : Ord.factor w/ 3 levels "0"<"1"<"2"
relGPA : num
mcAvgGPA : num
admitschool : Factor w/ 4 levels "1","2","3","4"
appyear : Factor w/ 11 levels "1999","2000",..
- tanto lrm() como polr() ahora se ejecutan con éxito, y ambos tratan con appyear dejando caer algunos valores del factor. lrm() deja caer 1999, 2000 y 2001 mientras que polr() sólo deja caer 1999 y 2000. lrm() no da advertencias mientras que polr() dice "el diseño parece ser deficiente en el rango, por lo que deja caer algunos coefs". Esto es una mejora, pero sigo sin entender por qué hay que eliminar más de un valor. xtabs muestra que no hay una separación completa, ¿verdad?
de salida:
xtabs(~fsd$appyear + fsd$Outcome)
fsd$Outcome
fsd$appyear 0 1 2
1999 1207 123 418
2000 1833 246 510
2001 1805 294 553
2002 1167 177 598
2003 4070 158 1076
2004 2803 106 1138
2005 3749 513 2141
2006 4429 519 2028
2007 6134 670 1947
2008 7446 662 1994
2009 4411 86 1118