Regresión Logística

Al estudiar problemas de conteos, la llamada inflación de ceros ocurre frecuentemente: los datos contienen más ceros de los que ocurrirían según las distribuciones habituales (Poisson, binomial negativa). Un modelo con inflación de ceros es una mezcla (mixtura) de un modelo de conteos y una distribución de Dirac (en cero).

Las técnicas habituales para resolverlos involucran (explícita o implícitamente) una estructura jerárquica de modelos: primero, uno (similar a una logística), separa las observaciones que corresponderían a la Dirac del resto. Un segundo modelo de conteos trata de ajustar el segundo.

En general, dos variables interaccionan cuando el efecto de una cambia al modificarse el nivel de la otra. Un caso particular (aunque notable) de interacción es el habitual en los modelos lineales, generalizados o no. En ellos, al introducir en el modelo términos del tipo x1 * x2, estamos indicando que el coeficiente de la segunda variable, $x_2$, es $\alpha + \beta x_1$. El efecto de un incremento de una unidad de $x_2$ depende entonces de $x_1$.

Impartí un curso sobre máquinas de vector soporte (SVMs en lo que sigue) en Lima el pasado mes de agosto.

Las SVMs (o más propiamente, los clasificadores de margen máximo) son exóticos dentro del repertorio del científico de datos. Lo que buscan es un hiperplano que maximiza el margen entre tirios o troyanos,

con o sin penalización para los puntos que insisten en permanecer en la región del espacio que no les corresponde. El modelo se ajusta resolviendo un problema de minimización inhabitual: uno de los llamados programas cuadráticos convexos. (Del que no nos tenemos que preocupar habitualmente porque delegamos la resolución en el software).

Una pregunta reciente en r-help-es se refería a la comparación en R de las proporciones en tres grupos. Obviando algunas pequeñas complicaciones en el problema, la respuesta canónica podría ser esta:

total <- c(56, 49,51)
positivos <- c(14, 10, 17)
prop.test(tmp$positivos, tmp$positivos + tmp$negativos)

# 3-sample test for equality of proportions without continuity correction
#
# data:  tmp$positivos out of tmp$positivos + tmp$negativos
# X-squared = 2.2289, df = 2, p-value = 0.3281
# alternative hypothesis: two.sided
# sample estimates:
#   prop 1    prop 2    prop 3
# 0.2500000 0.2040816 0.3333333

Los grupos no parecen ser desiguales.

En ocasiones, el conjunto de datos sobre el que se ajusta una regresión logística está desequilibrado con respecto a la población subyacente. Por ejemplo, puede suceder que la tasa de casos positivos en los datos sea del 20% mientras que en la población general es del 5%.

Esto puede suceder por varios motivos. El sobremuestreo es uno de ellos: se sobremuestrea cuando se toman, por ejemplo, todos los casos positivos y solo un subconjunto de los negativos.

Rescato y reconvierto un comentario de mi buen amigo José Luis Cañadas en una entrada mía reciente en la de hoy.

Sugiere José Luis el uso del paquete effects de R para estudiar el efecto de (que el caso concreto de interés, aunque hay otros) las variables de un modelo logístico.

Nos copia el código

library(effects)
mod.cowles <- glm(volunteer ~ sex + neuroticism*extraversion,
    data = Cowles, family = binomial)
eff.cowles <- allEffects(mod.cowles,
    xlevels = list(extraversion = seq(0, 24, 6)),
    given.values = c(sexmale = 0.5))
plot(eff.cowles, type = "response")

que genera

Los coeficientes de la regresión logística tienen una interpretación recta en términos de odds ratio. Que es un concepto sobre el que puede que alguien tenga algún tipo de intuición. Pero yo no.

¿Cómo podemos interpretar, aunque sea de manera rápida y grosera, los coeficientes? En términos de la variación de la probabilidad cuando la variable correspondiente cambia de valor (p.e., en una unidad). El problema es que la probabilidad depende del valor del resto de las variables: la relación no es lineal. No obstante, esa intuición es posible (en algunos casos: véase la nota).

Preguntaba el otro día Emilio Torres esto en R-help-es. Resumo la pregunta. Se trata de una simulación de unos datos y su ajuste mediante una regresión logística para ver si los coeficientes obtenidos son o no los esperados (teóricamente y por construcción).

El código de Emilio (cuyos resultados no podemos reproducir porque no nos ha contado qué similla usa) es

logisticsimulation <- function(n){
  dat <- data.frame(x1=sample(0:1, n,replace=TRUE),
                    x2=sample(0:1, n,replace=TRUE))
  odds <- exp(-1 - 4 * dat$x1 + 7*dat$x2 - 1 *dat$x1* dat$x2 )
  pr <- odds/(1+odds)
  res <- replicate(100, {
    dat$y <- rbinom(n,1,pr)
    coef(glm(y ~ x1*x2, data = dat, family = binomial()))
  })
  t(res)
}

res <- logisticsimulation(100)
apply(res,2,median)
## (Intercept)          x1          x2       x1:x2
## -1.0986123 -18.4674562  20.4823593  -0.0512933

Efectivamente, los coeficientes están lejos de los esperados, i.e., -1, -4, 7 y 1.

Creo que todos sabéis la historia de las admisiones de la Universidad de Berkeley y la paradoja de Simpson. Con palabras, muchas palabras, está contado, por ejemplo, aquí. Y si buscáis ubc admissions simpson en Google la encontraréis también en modo --verbose en muchos más sitios.

En R puede resumirse en

library(reshape2)
library(plyr)

data(UCBAdmissions)

raw <- as.data.frame(UCBAdmissions)

dat <- dcast(raw, Gender + Dept ~ <a href="http://inside-r.org/packages/cran/AdMit">Admit)

mod.0 <- glm(cbind(Admitted, Rejected) ~ Gender, data = dat, family = binomial)
mod.1 <- glm(cbind(Admitted, Rejected) ~ Gender + Dept, data = dat, family = binomial)

Echad un vistazo a los coeficientes de Gender en ambos modelos y veréis.

Esta entrada viene a cuento de una pregunta en r-help-es con, por referencia, este contexto:

Tengo un dataset con 4505 observaciones en el que la variable dependiente son presencias (n=97 y clasificadas como 1) y ausencias (n=4408 y clasificadas como 0).

Y la cuestión tiene que ver con la conveniencia de utilizar una muestra equilibrada o no de los datos al ajustar una regresión logística y si procede o no utilizar pesos.