8.5 Regresión con variables categóricas

La función lm() admite también variables categóricas (factores), lo que equivaldría a modelos de análisis de la varianza o de la covarianza.

Como ejemplo, en el resto del tema emplearemos los datos de empleados:

load("datos/empleados.RData")
datos <- with(empleados, data.frame(lnsal = log(salario), lnsalini = log(salini), catlab, sexo))

Al incluir variables categóricas la función lm() genera las variables indicadoras (variables dummy) que sean necesarias. Por ejemplo, la función model.matrix() construye la denominada matriz de diseño \(X\) de un modelo lineal: \[\mathbf{Y}=X\mathbf{\beta}+\mathbf{\varepsilon}\] En el caso de una variable categórica, por defecto se toma la primera categoría como referencia y se generan variables indicadoras del resto de categorías:

X <- model.matrix(lnsal ~ catlab, datos)
head(X)
##   (Intercept) catlabSeguridad catlabDirectivo
## 1           1               0               1
## 2           1               0               0
## 3           1               0               0
## 4           1               0               0
## 5           1               0               0
## 6           1               0               0

En el correspondiente ajuste (análisis de la varianza de un factor):

modelo <- lm(lnsal ~ catlab, datos)
summary(modelo)
## 
## Call:
## lm(formula = lnsal ~ catlab, data = datos)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.58352 -0.15983 -0.01012  0.13277  1.08725 
## 
## Coefficients:
##                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)     10.20254    0.01280 797.245  < 2e-16 ***
## catlabSeguridad  0.13492    0.04864   2.774  0.00576 ** 
## catlabDirectivo  0.82709    0.02952  28.017  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2438 on 471 degrees of freedom
## Multiple R-squared:  0.625,  Adjusted R-squared:  0.6234 
## F-statistic: 392.6 on 2 and 471 DF,  p-value: < 2.2e-16

el nivel de referencia no tiene asociado un coeficiente (su efecto se corresponde con (Intercept)). Los coeficientes del resto de niveles miden el cambio que se produce en la media al cambiar desde la categoría de referencia (diferencias de efectos respecto al nivel de referencia).

Para contrastar el efecto de los factores, es preferible emplear la función anova:

modelo <- lm(lnsal ~ catlab + sexo, datos)
anova(modelo)
## Analysis of Variance Table
## 
## Response: lnsal
##            Df Sum Sq Mean Sq F value    Pr(>F)    
## catlab      2 46.674 23.3372  489.59 < 2.2e-16 ***
## sexo        1  5.596  5.5965  117.41 < 2.2e-16 ***
## Residuals 470 22.404  0.0477                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Notas:

  • Para centrarse en las efectos de los factores, se puede emplear la función aov (analysis of variance; ver también model.tables() y TukeyHSD()). Esta función llama internamente a lm() (utilizando la misma parametrización).

  • Para utilizar distintas parametrizaciones de los efectos se puede emplear el argumento contrasts = c("contr.treatment", "contr.poly") (ver help(contrasts)).