De lo lineal a lo generalizado: introducción a los modelos GLM

Tasta ahora, hemos trabajado con el Modelo Lineal Clásico (LM), una herramienta fundamental en estadística. En un modelo lineal, asumimos que nuestra variable de respuesta (o dependiente) es continua y que los errores o residuales siguen una distribución normal.

Pensemos en un ejemplo clásico en comportamiento alimentario: queremos predecir el Índice de Masa Corporal (IMC) de un paciente basándonos en su ingesta calórica diaria y sus niveles de ansiedad. La estructura del modelo es:

\[Y_i = \beta_0 + \beta_1 X_{1i} + \beta_2 X_{2i} + \epsilon_i\]

Donde \(Y_i\) (el IMC) es continuo, y \(\epsilon_i\) representa el error aleatorio normalmente distribuido. Cuando nuestros datos cumplen estos supuestos, la función lm() en R es perfecta.

¿Cuál es el problema con los Modelos Lineales Clásicos?

En la investigación del comportamiento alimentario, la realidad clínica a menudo no se ajusta a una “distribución normal” ni a variables continuas. ¿Qué pasa cuando queremos investigar preguntas como las siguientes?

  1. Variables binarias (Sí/No): Queremos predecir la presencia o ausencia de un patrón de alimentación disfuncional (ej. 1 = Sí presenta el patrón, 0 = No lo presenta) a partir de ciertos rasgos apetitivos.

  2. Variables de conteo: Queremos modelar el número de episodios de atracón que un paciente reporta en una semana (ej. 0, 1, 2, 5 episodios).

Si intentamos trazar una línea recta clásica a través de datos que solo son ceros y unos (presencia/ausencia), la línea predecirá valores absurdos (como probabilidades negativas o mayores al 100%).

¿Qué son los Modelos Lineales Generalizados (GLM)?

Como su nombre lo indica, un GLM es una ampliación o generalización del modelo lineal clásico. Los GLM nos permiten utilizar el mismo marco lógico de la regresión lineal, pero flexibilizando las reglas para poder trabajar con variables de respuesta que no son normales, como conteos, proporciones o categorías binarias.

Para entender cómo logran esto, debemos ver los tres componentes fundamentales de cualquier GLM:

1. El Componente Aleatorio (Distribución de probabilidad)

En lugar de estar atados exclusivamente a la distribución Normal (Campana de Gauss), los GLM nos permiten elegir una distribución de la “familia exponencial” que se ajuste a la naturaleza de nuestra variable.

  • ¿Es presencia/ausencia de un trastorno? Usamos la familia Binomial.

  • ¿Es el conteo de episodios de purga? Usamos la familia Poisson.

  • ¿Es el IMC continuo? Usamos la familia Gaussiana (Normal).

2. El Componente Sistemático (Predictor Lineal)

Esta es la parte matemática que ya conoces. Es la combinación de nuestras variables independientes (predictores) multiplicadas por sus coeficientes. Se suele representar con la letra griega eta (\(\eta\)): \[\eta_i = \beta_0 + \beta_1 X_{1i} + \dots + \beta_k X_{ki}\]

En Resumen: Por qué los GLM son “más amplios”

La belleza conceptual de los GLM es que los modelos LM clásicos son simplemente un caso especial dentro del universo de los GLM.

Si a un GLM le indicas que su distribución es Gaussiana (Normal) y que su función de enlace es la Identidad (es decir, que no transforme nada), el resultado matemático es exactamente idéntico a un modelo lineal estándar.

Por lo tanto, aprender a usar la función glm() en R (en lugar de solo lm()) te abre la puerta a modelar casi cualquier tipo de fenómeno clínico, conductual o nutricional, ajustando únicamente los parámetros de family y link.

Ejemplos

1. Variables del dataset

Nombre del dataset:

behavior_food_glm

Variable Tipo Descripción
edad continua edad del participante
sexo binaria 0 = mujer, 1 = hombre
estres escala nivel de estrés (1–10)
sueno continua horas de sueño
actividad continua minutos de actividad física por día
impulsividad escala impulsividad alimentaria (1–10)
IMC continua índice de masa corporal
obesidad binaria 1 = IMC ≥ 30
snacks_semana conteo número de snacks consumidos por semana
episodios_emocionales conteo episodios de comer emocional
tiempo_redes continua positiva minutos diarios en redes sociales
set.seed(123)

n <- 200

edad <- round(rnorm(n, 35, 10))
sexo <- rbinom(n, 1, 0.5)

estres <- round(runif(n, 1, 10),1)

sueno <- round(rnorm(n, 7, 1.2),1)

actividad <- round(rnorm(n, 45, 20))

impulsividad <- round(runif(n, 1, 10),1)

IMC <- round(22 + 
             0.4*estres + 
             0.3*impulsividad - 
             0.2*actividad/10 +
             rnorm(n,0,2),1)

obesidad <- ifelse(IMC >= 30,1,0)

snacks_semana <- rpois(n, 
                       lambda = exp(1 + 
                       0.08*estres + 
                       0.06*impulsividad))

episodios_emocionales <- rpois(n,
                               lambda = exp(0.5 +
                               0.12*estres +
                               0.08*impulsividad))

tiempo_redes <- rgamma(n,
                       shape = 5,
                       scale = 20 + 2*estres)

datos <- data.frame(
edad,
sexo,
estres,
sueno,
actividad,
impulsividad,
IMC,
obesidad,
snacks_semana,
episodios_emocionales,
tiempo_redes
)

head(datos)
##   edad sexo estres sueno actividad impulsividad  IMC obesidad snacks_semana
## 1   29    1    3.1   6.9        66          6.8 21.9        0             7
## 2   33    0    7.2   5.6        44          2.1 23.0        0             7
## 3   51    1    3.0   6.2        44          3.3 23.9        0             6
## 4   36    1    3.9   7.0        15          8.4 29.1        0             2
## 5   36    0    2.6   7.8        61          8.2 26.4        0             3
## 6   52    0    8.2   5.0        41          1.4 23.6        0             6
##   episodios_emocionales tiempo_redes
## 1                     6    217.88300
## 2                     5     89.09039
## 3                     1    159.86682
## 4                     1     84.79594
## 5                     7     40.85726
## 6                     2    115.55849

2. Modelos que puedes estimar en clase

A continuación, exploraremos cómo la naturaleza de nuestra variable de respuesta dicta el tipo de modelo que debemos utilizar.

1. Regresión Lineal Clásica (LM)

Uso: Variable de respuesta continua y teóricamente normal.

Ejemplo en clase: Queremos evaluar cómo el estrés, la actividad física y la edad explican las variaciones en el Índice de Masa Corporal (IMC). Asumimos que el IMC es una métrica continua que puede tomar cualquier valor dentro de un rango razonable.

modelo_lm <- lm(IMC ~ estres + actividad + edad, data = datos)

summary(modelo_lm)
## 
## Call:
## lm(formula = IMC ~ estres + actividad + edad, data = datos)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -5.9207 -1.4433 -0.1031  1.2254  5.3009 
## 
## Coefficients:
##              Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 23.978851   0.810792  29.575  < 2e-16 ***
## estres       0.268310   0.060392   4.443 1.48e-05 ***
## actividad   -0.011303   0.007718  -1.464    0.145    
## edad         0.005123   0.016763   0.306    0.760    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.225 on 196 degrees of freedom
## Multiple R-squared:  0.1048, Adjusted R-squared:  0.09114 
## F-statistic: 7.652 on 3 and 196 DF,  p-value: 7.309e-05

Guía de Interpretación:

Los coeficientes nos hablan de cambios absolutos:

  • Un coeficiente positivo en estres indica que, a mayor estrés, mayor es el IMC esperado (manteniendo lo demás constante).

  • Un coeficiente negativo en actividad indica el efecto protector: a mayor actividad, menor IMC.

2. Regresión Logística (GLM Binomial)

Uso: Variable de respuesta binaria o dicotómica (0 o 1).

Ejemplo en clase: Ahora no nos interesa el IMC exacto, sino clasificar si un paciente tiene obesidad clínica (1 = Sí, 0 = No). Un modelo lineal clásico podría predecir probabilidades sin sentido (ej. \(-15\%\) o \(120\%\)). Al usar family = binomial, la función de enlace logit comprime las predicciones para que vivan estrictamente entre \(0\) y \(1\) (probabilidad).

modelo_logit <- glm(obesidad ~ estres + impulsividad + actividad,
                    family = binomial,
                    data = datos)

summary(modelo_logit)
## 
## Call:
## glm(formula = obesidad ~ estres + impulsividad + actividad, family = binomial, 
##     data = datos)
## 
## Coefficients:
##                Estimate Std. Error z value Pr(>|z|)
## (Intercept)  -3.0514402  2.0065835  -1.521    0.128
## estres        0.0003926  0.1984367   0.002    0.998
## impulsividad -0.0144320  0.1968880  -0.073    0.942
## actividad    -0.0186707  0.0246006  -0.759    0.448
## 
## (Dispersion parameter for binomial family taken to be 1)
## 
##     Null deviance: 39.216  on 199  degrees of freedom
## Residual deviance: 38.624  on 196  degrees of freedom
## AIC: 46.624
## 
## Number of Fisher Scoring iterations: 7

Interpretación mediante Odds Ratios (OR):

En regresión logística, los coeficientes crudos (los log-odds) son difíciles de interpretar clínicamente. Por ello, los exponenciamos (calculando \(e^{\beta}\)) para obtener los Odds Ratios (Razones de Momios).

  • Un OR \(> 1\) indica que la variable aumenta el riesgo de presentar obesidad.

  • Un OR \(< 1\) indica un factor protector.

# Cálculo de los Odds Ratios
exp(coef(modelo_logit))
##  (Intercept)       estres impulsividad    actividad 
##   0.04729077   1.00039271   0.98567163   0.98150251

3. Regresión Poisson (GLM Poisson)

Uso: Variables de conteo discreto (números enteros \(\ge 0\)). No puedes tener eventos negativos ni “fracciones de evento” (ej. no puedes consumir -3 snacks o 2.5 snacks enteros en ocasiones separadas).

Ejemplo en clase: Queremos predecir el número de snacks ultraprocesados consumidos a la semana basándonos en el estrés y la impulsividad.

modelo_poisson <- glm(snacks_semana ~ estres + impulsividad,
                      family = poisson,
                      data = datos)

summary(modelo_poisson)
## 
## Call:
## glm(formula = snacks_semana ~ estres + impulsividad, family = poisson, 
##     data = datos)
## 
## Coefficients:
##              Estimate Std. Error z value Pr(>|z|)    
## (Intercept)   0.98781    0.10369   9.527  < 2e-16 ***
## estres        0.09290    0.01113   8.347  < 2e-16 ***
## impulsividad  0.05145    0.01071   4.806 1.54e-06 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for poisson family taken to be 1)
## 
##     Null deviance: 305.73  on 199  degrees of freedom
## Residual deviance: 223.33  on 197  degrees of freedom
## AIC: 940.01
## 
## Number of Fisher Scoring iterations: 4

Guía de Interpretación:

El modelo Poisson asume que la media y la varianza de los datos son iguales. Los coeficientes exponenciados representan razones de tasas (Incident Rate Ratios). Si el coeficiente exponencial del estrés es \(1.08\), significa que por cada punto extra de estrés, el número de snacks consumidos aumenta en un \(8\%\).

4. Regresión Gamma (GLM Gamma)

Uso: Variables continuas y estrictamente positivas, que a menudo presentan asimetría positiva (una “cola larga” hacia la derecha).

Ejemplo en clase: El tiempo invertido en redes sociales tiene un límite inferior estricto (\(0\) minutos) pero no tiene un límite superior fijo. La mayoría de las personas lo usa moderadamente, pero unos pocos pacientes reportan un uso extremo. La distribución Normal no ajusta bien aquí, pero la distribución Gamma sí.

modelo_gamma <- glm(tiempo_redes ~ estres + sueno,
                    family = Gamma(link="log"),
                    data = datos)

summary(modelo_gamma)
## 
## Call:
## glm(formula = tiempo_redes ~ estres + sueno, family = Gamma(link = "log"), 
##     data = datos)
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  4.67569    0.21271  21.982  < 2e-16 ***
## estres       0.05514    0.01240   4.448 1.45e-05 ***
## sueno        0.01080    0.02806   0.385    0.701    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for Gamma family taken to be 0.2106307)
## 
##     Null deviance: 46.125  on 199  degrees of freedom
## Residual deviance: 42.062  on 197  degrees of freedom
## AIC: 2249.5
## 
## Number of Fisher Scoring iterations: 5

5. Regresión Binomial Negativa (Para Sobredispersión)

Uso: Variables de conteo donde ocurre sobredispersión (la varianza es mucho mayor que la media).

Ejemplo en clase: Los episodios de comer emocional son conteos. Podríamos intentar usar Poisson. Sin embargo, en psicología clínica es común que gran parte de la muestra reporte \(0\) o \(1\) episodio, mientras que un subgrupo clínico reporta números muy altos (ej. \(15\) o \(20\) episodios). Esto “rompe” el supuesto de la regresión de Poisson (donde media \(=\) varianza). El modelo Binomial Negativo introduce un parámetro extra que relaja este supuesto y corrige el error.

Nota: Para esto usamos la función glm.nb de la paquetería MASS.

library(MASS)

modelo_nb <- glm.nb(episodios_emocionales ~ estres + impulsividad,
                    data = datos)
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
## Warning in theta.ml(Y, mu, sum(w), w, limit = control$maxit, trace =
## control$trace > : iteration limit reached
summary(modelo_nb)
## 
## Call:
## glm.nb(formula = episodios_emocionales ~ estres + impulsividad, 
##     data = datos, init.theta = 74789.86771, link = log)
## 
## Coefficients:
##              Estimate Std. Error z value Pr(>|z|)    
## (Intercept)   0.45362    0.11675   3.885 0.000102 ***
## estres        0.11677    0.01226   9.523  < 2e-16 ***
## impulsividad  0.08399    0.01186   7.080 1.44e-12 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## (Dispersion parameter for Negative Binomial(74789.87) family taken to be 1)
## 
##     Null deviance: 299.55  on 199  degrees of freedom
## Residual deviance: 175.31  on 197  degrees of freedom
## AIC: 852.52
## 
## Number of Fisher Scoring iterations: 1
## 
## 
##               Theta:  74790 
##           Std. Err.:  869484 
## Warning while fitting theta: iteration limit reached 
## 
##  2 x log-likelihood:  -844.521

6. Visualización de Modelos Generalizados

El paquete ggplot2 es extremadamente útil para graficar modelos GLM, ya que nos permite visualizar la curva de predicción no lineal de vuelta en la escala original de los datos. Aquí graficamos la predicción de conteo (Poisson) del número de snacks en función del estrés.

library(ggplot2)
## Warning: package 'ggplot2' was built under R version 4.4.3
ggplot(datos, aes(x = estres, y = snacks_semana)) +
  geom_point(alpha = 0.5, color = "darkorange") +
  geom_smooth(method = "glm", 
              method.args = list(family = "poisson"), 
              color = "darkblue") +
  theme_minimal() +
  labs(title = "Predicción de Consumo de Snacks mediante Modelo Poisson",
       x = "Nivel de Estrés Reportado (1-10)",
       y = "Snacks Ultraprocesados por Semana")
## `geom_smooth()` using formula = 'y ~ x'