8 Regresión dummy
Jorge Méndez González (UAAAN)
2022-10-25
Introducción
La regresión lineal no se limita estrictamente a los datos cuantitativos. Un modelo de regresión lineal también puede tener datos cualitativos. En este video lo explico todo.
Ocasionalmente se debe introducir variables categóricas (o cualitativas) con dos o más categorias. Por ejemplo, género (hombre, mujer), estado civil, partido politico. Estas se pueden representar como variables dummy o indicadoras. Estas variables toman dos (DOS) valores usualmente, cero (0) y uno (1). Los dos valores significan que la observación pertenece a una de dos categorias. Las variables dummy o indicadoras sirven para identificar categorias o clase a las que pertenecen las observaciones.
\(\begin{aligned} D = \begin{cases} 1 & \text{if characteristic is present} \\ 0 & \text{if characteristic is not present} \end{cases} \end{aligned}\)
- Las variables ficticias se pueden utilizar para capturar cambios en la intercepta, la pendiente o ambas de un modelo de regresion. esto es regresion dummy.
La informacion
Librerias
library(car) # extracts model results
library(tidyverse)
library(knitr)
library(ggplot2)
library(performance)
library(DT)
library(effects)
library(equatiomatic) # takes a model and returns to 'LaTeX'
library(vembedr)
Datos
“Salarios para profesores”
El salario académico de nueve meses de 2008-09 para profesores asistentes, profesores asociados y profesores en una universidad en los EE. UU. Los datos se recopilaron como parte del esfuerzo contÍnuo de la administración de la universidad para monitorear las diferencias salariales entre los profesores masculinos y femeninos. Un set de datos con 397 observaciones y 6 variables.
# setwd("c:/CURSO REG JMG")
head(Salaries) # los datos## rank discipline yrs.since.phd yrs.service sex salary
## 1 Prof B 19 18 Male 139750
## 2 Prof B 20 16 Male 173200
## 3 AsstProf B 4 3 Male 79750
## 4 Prof B 45 39 Male 115000
## 5 Prof B 40 41 Male 141500
## 6 AssocProf B 6 6 Male 97000head(Datos<-Salaries) # cambiar el nombre al set## rank discipline yrs.since.phd yrs.service sex salary
## 1 Prof B 19 18 Male 139750
## 2 Prof B 20 16 Male 173200
## 3 AsstProf B 4 3 Male 79750
## 4 Prof B 45 39 Male 115000
## 5 Prof B 40 41 Male 141500
## 6 AssocProf B 6 6 Male 97000# head(Datos[with(Datos, order(yrs.service)), ]) # ordenar los datosDatos$sex_D <- factor(Datos$sex, labels=c("0", "1")) # convertir la varable sexo a factor 0 y 1
attach(Datos)
view(Datos)AQUI HAY VARIOS SET DE DATOS (https://r-data.pmagunia.com/dash)
La variable dummy
datatable(Datos, extensions = "Buttons",
options = list(dom = "Bfrtip",
buttons = c("copy", "csv", "excel", "pdf")))# Exportar los datos a csv para explicar
write.csv(Datos, "Datos con la var dummy.csv")Opcion 1: Una sola regresion
La plot
color<-as.factor(sex)
plot(yrs.service, salary, col=color, cex=1.5, pch=16)
La regresion sin considerar genero
Fit0 <- lm(salary ~ yrs.service) # sin considear el genero
summary(Fit0) # los coeficicientes de regresion##
## Call:
## lm(formula = salary ~ yrs.service)
##
## Residuals:
## Min 1Q Median 3Q Max
## -81933 -20511 -3776 16417 101947
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 99974.7 2416.6 41.37 < 2e-16 ***
## yrs.service 779.6 110.4 7.06 7.53e-12 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 28580 on 395 degrees of freedom
## Multiple R-squared: 0.1121, Adjusted R-squared: 0.1098
## F-statistic: 49.85 on 1 and 395 DF, p-value: 7.529e-12extract_eq(Fit0) # El modelo final\[ \operatorname{salary} = \alpha + \beta_{1}(\operatorname{yrs.service}) + \epsilon \]
extract_eq(Fit0, wrap = TRUE, use_coefs = TRUE) # modelo final y coefficientes\[ \begin{aligned} \operatorname{\widehat{salary}} &= 99974.65 + 779.57(\operatorname{yrs.service}) \end{aligned} \]
Opcion 2: Regresion Dummy (genero)
Aqui haremos regresion unicamente con la variable dummy (genero), cualitativa o indicadora.
Metodo 1
Fit1 <- lm(salary ~ sex) # opcion 1: con la variable original
summary(Fit1)##
## Call:
## lm(formula = salary ~ sex)
##
## Residuals:
## Min 1Q Median 3Q Max
## -57290 -23502 -6828 19710 116455
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 101002 4809 21.001 < 2e-16 ***
## sexMale 14088 5065 2.782 0.00567 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 30030 on 395 degrees of freedom
## Multiple R-squared: 0.01921, Adjusted R-squared: 0.01673
## F-statistic: 7.738 on 1 and 395 DF, p-value: 0.005667plot(predictorEffects(Fit1))
extract_eq(Fit1) # extraer el modelo \[ \operatorname{salary} = \alpha + \beta_{1}(\operatorname{sex}_{\operatorname{Male}}) + \epsilon \]
extract_eq(Fit1, wrap = TRUE, use_coefs = TRUE) # modelo final y coefficientes\[ \begin{aligned} \operatorname{\widehat{salary}} &= 101002.41 + 14088.01(\operatorname{sex}_{\operatorname{Male}}) \end{aligned} \]
Metodo 2
Fit1_a <- lm(salary ~ sex_D) # opcion 2: con la variable dummy
summary(Fit1_a)##
## Call:
## lm(formula = salary ~ sex_D)
##
## Residuals:
## Min 1Q Median 3Q Max
## -57290 -23502 -6828 19710 116455
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 101002 4809 21.001 < 2e-16 ***
## sex_D1 14088 5065 2.782 0.00567 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 30030 on 395 degrees of freedom
## Multiple R-squared: 0.01921, Adjusted R-squared: 0.01673
## F-statistic: 7.738 on 1 and 395 DF, p-value: 0.005667Esto equivale a…
Male <- Datos %>% filter(sex == "Male") %>% pull(salary)
Female <- Datos %>% filter(sex == "Female") %>% pull(salary)
t.test(
x = Male,
y = Female,
alternative = "two.sided",
mu = 0,
var.equal = TRUE,
conf.level = 0.95
)##
## Two Sample t-test
##
## data: Male and Female
## t = 2.7817, df = 395, p-value = 0.005667
## alternative hypothesis: true difference in means is not equal to 0
## 95 percent confidence interval:
## 4131.107 24044.910
## sample estimates:
## mean of x mean of y
## 115090.4 101002.4ggplot(Datos, aes(x = sex, y = salary, colour = sex)) +
geom_boxplot(outlier.shape = NA) +
geom_jitter()
Opcion 3: Intercept dummy (reg. misma pendiente)
Quizás el uso más común de las variables ficticias es modificar el parámetro de la intercepta \(\beta_{0}\)del modelo de regresión. Agregar una variable indicadora D al modelo, junto con un nuevo parámetro δ nos da.
\(\begin{aligned} y = \beta_{0} + \delta D + \beta_{1}x + e \end{aligned}\)
\(\begin{aligned} \mathbf{E}(\text{salary}) = \begin{cases} (\beta_{0} + \delta) + \beta_{1}\text{yrs.service}, & \text{Male}=1 \\ \beta_{0} + \beta_{1}\text{yrs.service}, & \text{Female} = 0 \end{cases} \end{aligned}\)
Metodo 1
Fit2_a <- lm(salary ~ factor(sex) + yrs.service ) # con la variable original
summary(Fit2_a)##
## Call:
## lm(formula = salary ~ factor(sex) + yrs.service)
##
## Residuals:
## Min 1Q Median 3Q Max
## -81757 -20614 -3376 16779 101707
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 92356.9 4740.2 19.484 < 2e-16 ***
## factor(sex)Male 9071.8 4861.6 1.866 0.0628 .
## yrs.service 747.6 111.4 6.711 6.74e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 28490 on 394 degrees of freedom
## Multiple R-squared: 0.1198, Adjusted R-squared: 0.1154
## F-statistic: 26.82 on 2 and 394 DF, p-value: 1.201e-11Metodo 2
Fit2 <- lm(salary ~ sex_D + yrs.service ) # con la variable dummy
summary(Fit2)##
## Call:
## lm(formula = salary ~ sex_D + yrs.service)
##
## Residuals:
## Min 1Q Median 3Q Max
## -81757 -20614 -3376 16779 101707
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 92356.9 4740.2 19.484 < 2e-16 ***
## sex_D1 9071.8 4861.6 1.866 0.0628 .
## yrs.service 747.6 111.4 6.711 6.74e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 28490 on 394 degrees of freedom
## Multiple R-squared: 0.1198, Adjusted R-squared: 0.1154
## F-statistic: 26.82 on 2 and 394 DF, p-value: 1.201e-11extract_eq(Fit2) # extraeer el modelo \[ \operatorname{salary} = \alpha + \beta_{1}(\operatorname{sex\_D}_{\operatorname{1}}) + \beta_{2}(\operatorname{yrs.service}) + \epsilon \]
extract_eq(Fit2, wrap = TRUE, use_coefs = TRUE) # modelo final y coefficientes\[ \begin{aligned} \operatorname{\widehat{salary}} &= 92356.95 + 9071.8(\operatorname{sex\_D}_{\operatorname{1}}) + 747.61(\operatorname{yrs.service}) \end{aligned} \]
color<-as.factor(sex)
plot(yrs.service, salary, col=color, cex=1.5, pch=16)
points(yrs.service, fitted.values(Fit2), lwd=3, col="blue")
Opcion 4: Slope dummy (reg. dif. intercepta y pendiente)
También se pueden introducir variables ficticias para afectar la pendiente de un modelo de regresión. En lugar de suponer que el efecto del sexo provoca un cambio en el intercepto en nuestro modelo anterior, supongamos que el cambio está en la pendiente de la relación. Para lograr esto, se incluye una variable explicativa adicional en el modelo que sea igual al producto de una variable ficticia y una variable continua. Podemos especificar una forma genérica de tal modelo como:
\[ \begin{equation} y = \beta_{0} + \beta_{1}x + \gamma(x \times D) + e \end{equation} \]
La nueva variable (x×D) es el producto de x y la variable ficticia D. Esto es lo que se denomina una variable de interacción, ya que captura el efecto de interacción de la variable continua x y la variable ficticia D. Alternativamente, también se denomina como una variable indicadora de pendiente, o una variable ficticia de pendiente, porque permite un cambio en la pendiente de la relación.
\[ \begin{equation} salary = \beta_{0} + \beta_{1}yrs.service + \gamma(yrs.service \times D) + e \end{equation} \]
\[ \begin{aligned} \mathbf{E}(\text{salary }) = \beta_{0} + \beta_{1}\text{yrs.service} + \gamma(\text{yrs.service} \times \text{sex}) = \begin{cases} \beta_{0} + (\beta_{1} + \gamma)\text{yrs.service}, & \text{Male} = 1\\ \beta_{0} + \beta_{1}\text{yrs.service}, & \text{Female} = 0 \end{cases} \end{aligned} \]
Metodo1
Fit3 <- lm(salary ~ sex_D + yrs.service + sex_D*yrs.service) # dummy
summary(Fit3)##
## Call:
## lm(formula = salary ~ sex_D + yrs.service + sex_D * yrs.service)
##
## Residuals:
## Min 1Q Median 3Q Max
## -80381 -20258 -3727 16353 102536
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 82068.5 7568.7 10.843 < 2e-16 ***
## sex_D1 20128.6 7991.1 2.519 0.01217 *
## yrs.service 1637.3 523.0 3.130 0.00188 **
## sex_D1:yrs.service -931.7 535.2 -1.741 0.08251 .
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 28420 on 393 degrees of freedom
## Multiple R-squared: 0.1266, Adjusted R-squared: 0.1199
## F-statistic: 18.98 on 3 and 393 DF, p-value: 1.622e-11plot(yrs.service, salary, col=color, cex=1.5, pch=16)
points(yrs.service, fitted.values(Fit3), lwd=3, col="#7FFF00") # Agregar los estimados a la plot
extract_eq(Fit3) # extraeer el modelo \[ \operatorname{salary} = \alpha + \beta_{1}(\operatorname{sex\_D}_{\operatorname{1}}) + \beta_{2}(\operatorname{yrs.service}) + \beta_{3}(\operatorname{sex\_D}_{\operatorname{1}} \times \operatorname{yrs.service}) + \epsilon \]
extract_eq(Fit3, wrap = TRUE, use_coefs = TRUE) # modelo final y coefficientes\[ \begin{aligned} \operatorname{\widehat{salary}} &= 82068.51 + 20128.63(\operatorname{sex\_D}_{\operatorname{1}}) + 1637.3(\operatorname{yrs.service}) - 931.74(\operatorname{sex\_D}_{\operatorname{1}} \times \operatorname{yrs.service}) \end{aligned} \]