0.1 Introducción

Este es un documento de trabajo, contiene ejemplos para aplicar modelos de regresión lineal con los datos del Wooldridge 2010-Introducción a la Econometría (4ta edución), texto de referencia del curso de Econometría 1. El mismo fue elaborado en Markdown utilizando R studio (software que vemos en el curso). Aquí se encuentran detalles para utilizar R Markdown http://rmarkdown.rstudio.com.

El R permite instalar la librería Wooldridge de la siguiente manera.

if(!require(‘wooldridge’)) { install.packages(‘wooldridge’) library(‘wooldridge’) }

Esto es equivalente a lo siguiente: install.packages(“wooldridge”) #la instalo si es la primera vez que la llamo. library(wooldridge) #la llamo una vez que ya la tengo instalada.

En esta librería encontrará todas las bases de datos que contiene el libro. Se deberá utilizar el mismo nombre con el que son referenciadas en el texto. A modo de ejemplo, trabajaremos con la base wage1 que contiene datos de la Encuesta de población de 1976, recopilados por Henry Farber cuando él y Wooldridge fueron colegas en el MIT en 1988.

0.2 Modelo de regresión lineal simple (MRLS)

0.2.1 Cargar y manipular la base

Para cargar los datos se escribe data(“wage1”) # carga la base salarios base <-wage1 # la asigno al objeto base con el que trabajaremos

Acceso a la base y solicito los nombres de las variables:

library(wooldridge)
base <-data("wage1")
base <-wage1
#View(base) 
names(base)
##  [1] "wage"     "educ"     "exper"    "tenure"   "nonwhite" "female"  
##  [7] "married"  "numdep"   "smsa"     "northcen" "south"    "west"    
## [13] "construc" "ndurman"  "trcommpu" "trade"    "services" "profserv"
## [19] "profocc"  "clerocc"  "servocc"  "lwage"    "expersq"  "tenursq"

Las variables que utilizaremos son las siguientes:

  • wage: salario promedio por hora
  • educ: años de educación
  • exper: años de experiencia potencial
  • ternure: años con el empleador actual (antiguedad)
  • nonwhite: =1 si no la persona no es blanca
  • female: =1 si la persona es mujer
  • married: =1 si la persona es casada

Cambiamos el nombre de la variable que está en la posición 4

names(base)[4]<-"antiguedad"

Seleccionamos las variables que nos interesan por la posición de las columnas en el data frame. Indicamos que queremos las variables que están desde la posición 1 hasta la 7 y las que están en la posición 22, 23 y 24

base1<-base[,c(1: 7, 22, 23, 24)]

Miramos las observaciones de las primeras 5 personas de la base:

head(base1, n=5)
##   wage educ exper antiguedad nonwhite female married    lwage expersq tenursq
## 1 3.10   11     2          0        0      1       0 1.131402       4       0
## 2 3.24   12    22          2        0      1       1 1.175573     484       4
## 3 3.00   11     2          0        0      0       0 1.098612       4       0
## 4 6.00    8    44         28        0      0       1 1.791759    1936     784
## 5 5.30   12     7          2        0      0       1 1.667707      49       4

También se puede seleccionar utilizando el nombre de las variables:

datos1<-base[,c("wage","educ","exper", "antiguedad" )] 
head(datos1, n=5)
##   wage educ exper antiguedad
## 1 3.10   11     2          0
## 2 3.24   12    22          2
## 3 3.00   11     2          0
## 4 6.00    8    44         28
## 5 5.30   12     7          2

Otra forma mas sencilla es utilizar la función select del paquete dplyr, se usa el comando llamado pipe %>%, el mismo enlaza las órdenes de manera secuencial.

Si es la primera vez que se usa recordar instalarlo previamente: #install.packages(“dplyr”)

library(dplyr)
## 
## Attaching package: 'dplyr'
## The following objects are masked from 'package:stats':
## 
##     filter, lag
## The following objects are masked from 'package:base':
## 
##     intersect, setdiff, setequal, union
datos2<-base %>% 
        select("wage", "educ", "exper", "antiguedad")
head(datos2, n=5)
##   wage educ exper antiguedad
## 1 3.10   11     2          0
## 2 3.24   12    22          2
## 3 3.00   11     2          0
## 4 6.00    8    44         28
## 5 5.30   12     7          2

Se puede combinar con un filtro por filas con la función filter del paquete dplyr. A modo de ejemplo, seleccionamos solo a las personas casadas.

datos3<-base %>% 
  select("wage", "educ", "exper", "antiguedad", "married")%>% 
  filter(married==1)

#View(datos3)

0.2.2 Investigar las correlaciones

  plot(datos2)

Calculo la matriz de correlaciones de las variables que figuran en datos

  cor(datos2)
##                 wage        educ      exper  antiguedad
## wage       1.0000000  0.40590333  0.1129034  0.34688957
## educ       0.4059033  1.00000000 -0.2995418 -0.05617257
## exper      0.1129034 -0.29954184  1.0000000  0.49929145
## antiguedad 0.3468896 -0.05617257  0.4992914  1.00000000

Existe una correlación positiva entre el salario y la educación (0.4059)

Creamos el diagrama de dispersión entre salario y educación. Utilizaremos la librería ggplot2.

 #install.packages("ggplot2") 
library(ggplot2)
ggplot(datos2, aes(x=educ, y=wage)) + 
  geom_point() + theme_light() +
  ggtitle("Relación entre salario y educación")

0.2.3 Estimar un modelo de regresión lineal simple (MRLS)

Estimamos un MRLS con el método mínimos cuadrados ordinarios (en adelante MCO) que explique los salarios en función de los años de educación de las personas.

\[\begin{equation} wage_i = \beta _0+\beta _1 educ_i +u_i \\ \end{equation}\]

mod1 <- lm(wage ~ educ, data=datos2)
summary(mod1) #para imprimir la salida
## 
## Call:
## lm(formula = wage ~ educ, data = datos2)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -5.3396 -2.1501 -0.9674  1.1921 16.6085 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.90485    0.68497  -1.321    0.187    
## educ         0.54136    0.05325  10.167   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3.378 on 524 degrees of freedom
## Multiple R-squared:  0.1648, Adjusted R-squared:  0.1632 
## F-statistic: 103.4 on 1 and 524 DF,  p-value: < 2.2e-16
anova(mod1) #para imprimir el análisis de varianzas
## Analysis of Variance Table
## 
## Response: wage
##            Df Sum Sq Mean Sq F value    Pr(>F)    
## educ        1 1179.7 1179.73  103.36 < 2.2e-16 ***
## Residuals 524 5980.7   11.41                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

0.2.4 Objeto lm y sus componentes

Los resultados fueron asignados al objeto mod1 de clase lm. Para consultar la clase:

class(mod1)
## [1] "lm"

Dentro del objeto se encuentran 12 elementos que se pueden utilizar llamándolos por los siguientes nombres:

names(mod1)
##  [1] "coefficients"  "residuals"     "effects"       "rank"         
##  [5] "fitted.values" "assign"        "qr"            "df.residual"  
##  [9] "xlevels"       "call"          "terms"         "model"

Se puede acceder a ellos como a los objetos del data frame, utilizando el operador $ entre el objeto y el elemento. Por ejemplo, para extraer los coeficientes estimados se escribe lo siguiente:

mod1$coefficients
## (Intercept)        educ 
##  -0.9048516   0.5413593

Para consultar la estimación de un coeficiente de regresión se utilizan paréntesis rectos.

ahat <-coef(summary(mod1))[1,1]
bhat <- coef(summary(mod1))[2,1]

ahat
## [1] -0.9048516
bhat
## [1] 0.5413593

Miramos los componentes del modelo y los asignamos a nuevos objetos

names(mod1)
##  [1] "coefficients"  "residuals"     "effects"       "rank"         
##  [5] "fitted.values" "assign"        "qr"            "df.residual"  
##  [9] "xlevels"       "call"          "terms"         "model"
coeficientes<- mod1$coefficients
ygorro<- mod1$fitted.values
resid<-mod1$residuals
datos2$predicciones <- predict(mod1) 
head(datos2, 5)
##   wage educ exper antiguedad predicciones
## 1 3.10   11     2          0     5.050100
## 2 3.24   12    22          2     5.591459
## 3 3.00   11     2          0     5.050100
## 4 6.00    8    44         28     3.426022
## 5 5.30   12     7          2     5.591459

0.2.5 Diagrama de dispersión y recta de regresión

Creamos el gráfico de dispersión entre el salario y la educación e incluímos la estimación de la recta de regresión.

ggplot(datos2, aes(x=educ, y=wage)) + 
  geom_point() +
  geom_smooth(method='lm', formula=y~x, se=FALSE, col='dodgerblue1') +
      theme_light() +
   ggtitle("Relación entre salario y educación")

Es posible obtener un gráfico de dispersión interactivo. Posicionándose encima de cada observación se ven los valores de (x, y) para cada uno de los individuos. Necesita la función plotly.

library(plotly)
## 
## Attaching package: 'plotly'
## The following object is masked from 'package:ggplot2':
## 
##     last_plot
## The following object is masked from 'package:stats':
## 
##     filter
## The following object is masked from 'package:graphics':
## 
##     layout
plot_ly(data=datos2, x=~educ, y=~wage)
## No trace type specified:
##   Based on info supplied, a 'scatter' trace seems appropriate.
##   Read more about this trace type -> https://plotly.com/r/reference/#scatter
## No scatter mode specifed:
##   Setting the mode to markers
##   Read more about this attribute -> https://plotly.com/r/reference/#scatter-mode

0.2.6 Trabajar con los residuos del modelo

Calculamos los residuos del modelo:

datos2$residuos <- datos2$wage-datos2$predicciones
head(datos2, 5)
##   wage educ exper antiguedad predicciones   residuos
## 1 3.10   11     2          0     5.050100 -1.9501003
## 2 3.24   12    22          2     5.591459 -2.3514594
## 3 3.00   11     2          0     5.050100 -2.0501002
## 4 6.00    8    44         28     3.426022  2.5739776
## 5 5.30   12     7          2     5.591459 -0.2914593

Comparamos los residuos calculados con los que nos dio el modelo. \[\widehat{u}_{i}= y_{i}-(\widehat{\beta}_0 + \widehat{\beta}_1 x_{i}) = y_{i} - \widehat{y}_i \]

datos2$resid<- resid
datos2$verif<-datos2$residuos-datos2$resid

Agrego al gráfico de dispersión los valores estimados de \(y\) (\(\widehat{y}_i\)) en rojo y muestro los residuos (\(\widehat{u}_{i}\))

ggplot(datos2, aes(x=educ, y=wage)) +
  geom_smooth(method="lm", se=FALSE, color="lightgrey") +
  geom_segment(aes(xend=educ, yend=predicciones), col='red', lty='dashed') +
  geom_point() +
  geom_point(aes(y=predicciones), col='red') +
  theme_light()
## `geom_smooth()` using formula 'y ~ x'

Inspección gráfica de los resididuos

plot(datos2$residuos)

Graficamos las funciones de densidad estimadas

library(e1071)  # for skewness function
par(mfrow=c(1, 2))  # divide graph area in 2 columns

plot(density(datos2$wage), main="Density Plot: salario", ylab="Frequency", sub=paste("Skewness:", round(e1071::skewness(datos2$wage), 2)))  # density para 'salario'

polygon(density(datos2$wage), col="red")

plot(density(datos2$residuos), main="Density Plot: residuos", ylab="Frequency", sub=paste("Skewness:", round(e1071::skewness(datos2$residuos), 2)))  # density para los 'residuos'

polygon(density(datos2$residuos), col="red")

0.2.7 Predicción e intervalos de confianza

Genero la predicción puntual y el intervalo de confianza para una educación de 10 años, prediction devuelve el valor para la predicción puntual

Calculo el salario esperado para una persona con la educación promedio (\(y/x=\overline{x}\))

mean(datos2$educ)
## [1] 12.56274
sal_pred<-0.90485+0.054136*12.56
sal_pred
## [1] 1.584798

Calculo el intervalo de confianza para una persona con la educación promedio (\(x=12.56\))

nuevo <- data.frame(educ=12.56)
future_y<-predict(object=mod1, newdata=nuevo, interval="prediction", level=0.95)

Generamos el intervalo de confianza para \(E(y/x)\). Debemos cambiar el argumento de interval= “confidence”.

future_esp_y <- predict(object=mod1,newdata = nuevo, interval = "confidence", level = 0.95)
future_esp_y<-as.data.frame(future_esp_y)

IC_inf_esp_y<-future_esp_y$lwr
IC_sup_esp_y<-future_esp_y$upr

Agregamos los extremos de los intervalos de confianza a la base con la que estamos trabajando:

nuevos_datos <- cbind(datos2, future_y, IC_inf_esp_y, IC_sup_esp_y)

Generamos el gráfico con los intervalos de confianza para la predicción puntual y para el valor esperado.

IC_y<-ggplot(nuevos_datos, aes(x=educ, y=wage))+
  geom_point() +
  geom_line(aes(y=lwr), color="red", linetype="dashed") +
  geom_line(aes(y=upr), color="red", linetype="dashed") +
  geom_smooth(method=lm, formula=y~x, se=TRUE, level=0.95, col='blue', fill='pink2') +
  theme_light() + ggtitle("Precicción de y al 95%")


library(ggplot2)
IC_esp_y_95<-ggplot(nuevos_datos, aes(x=educ, y=wage))+
  geom_point() +
  geom_line(aes(y=IC_inf_esp_y), color="blue", linetype="dashed") +
  geom_line(aes(y=IC_inf_esp_y), color="blue", linetype="dashed") +
  geom_smooth(method=lm, formula=y~x, se=TRUE, level=0.95, col='blue', fill='pink2') +
  theme_light() + ggtitle("predicción de E(y/x) al 95%")

Imprimo los gráficos juntos para poder compararlos mejor. ¿Cuál de los dos tiene mayor amplitud?

#nstall.packages("gridExtra")
library(gridExtra)
## 
## Attaching package: 'gridExtra'
## The following object is masked from 'package:dplyr':
## 
##     combine
grid.arrange(IC_esp_y_95, IC_y, ncol=2,nrow=1)

Veamos el impacto que tiene el nivel de confianza en la amplitud de los intervalos.

library(ggplot2)

#install.packages("jtools")
#install.packages("ggstance")
#install.packages("broom.mixed")

library(jtools)
library(ggstance)
## 
## Attaching package: 'ggstance'
## The following objects are masked from 'package:ggplot2':
## 
##     geom_errorbarh, GeomErrorbarh
library(broom.mixed)
## Registered S3 methods overwritten by 'broom':
##   method            from  
##   tidy.glht         jtools
##   tidy.summary.glht jtools
a=plot_summs (mod1 , escala  =  TRUE , plot.distributions  =  TRUE , inner_ci_level  =  .90)
b=plot_summs (mod1 , escala  =  TRUE , plot.distributions  =  TRUE , inner_ci_level  =  .7) 


library(gridExtra)
grid.arrange(a, b, ncol=1,nrow=2)

0.3 Modelos de regresión lineal múltiple (MRLM)

Estimaremos modelos de regresión lineal múltiple utilizando la misma base de datos de Wooldridge, wage1 que contiene datos de la Encuesta de población de 1976.

0.3.1 Acceder y conocer la base

Para cargar los datos se escribe data(“wage1”) # carga la base salarios base <-wage1 # la asigno al objeto base con el que trabaremos

Accedemos a la base, la asignamos al objeto base y solicitamos los nombres de las variables:

library(wooldridge)
base <-data("wage1")
base <-wage1
#View(base) 
names(base)
##  [1] "wage"     "educ"     "exper"    "tenure"   "nonwhite" "female"  
##  [7] "married"  "numdep"   "smsa"     "northcen" "south"    "west"    
## [13] "construc" "ndurman"  "trcommpu" "trade"    "services" "profserv"
## [19] "profocc"  "clerocc"  "servocc"  "lwage"    "expersq"  "tenursq"

Las variables que utilizaremos son las siguientes:

  • wage: salario promedio por hora
  • educ: años de educación
  • exper: años de experiencia potencial
  • ternure: años con el empleador actual (antiguedad)
  • nonwhite: =1 si no la persona no es blanca
  • female: =1 si la persona es mujer
  • married: =1 si la persona es casada
  • lwage: log(salario)
  • expersq: exper ^2, el cuadrado de los años de experiencia potencial
  • tenursq: ternure^2, el cuadrado de la antiguedad en el trabajo actual

Cambiamos el nombre de la variable que está en la posición 4

names(base)[4]<-"antiguedad"

Seleccionamos las variables que nos interesan por la posición de las columnas en el data frame. Indicamos que queremos las variables que están desde la posición 1 hasta la 7 y las que están en la posición 22, 23 y 24. Las asignamos al objeto datos.

datos<-base[,c(1: 7, 22, 23, 24)]

Para mirar la dimensión de la base resultante solitamos la dimensión de la misma.

dim(datos)
## [1] 526  10

Podríamos solicitar el número de columnas (personas en este caso) y de columnas por separado (variables) con las siguientes funciones:

nrow(datos)
## [1] 526
ncol(datos)
## [1] 10

Miramos las observaciones de las primeras 5 personas de la base:

head(datos, n=5)
##   wage educ exper antiguedad nonwhite female married    lwage expersq tenursq
## 1 3.10   11     2          0        0      1       0 1.131402       4       0
## 2 3.24   12    22          2        0      1       1 1.175573     484       4
## 3 3.00   11     2          0        0      0       0 1.098612       4       0
## 4 6.00    8    44         28        0      0       1 1.791759    1936     784
## 5 5.30   12     7          2        0      0       1 1.667707      49       4

Solicitamos las estadísticas básicas de las variables de la base:

summary(datos)
##       wage             educ           exper         antiguedad    
##  Min.   : 0.530   Min.   : 0.00   Min.   : 1.00   Min.   : 0.000  
##  1st Qu.: 3.330   1st Qu.:12.00   1st Qu.: 5.00   1st Qu.: 0.000  
##  Median : 4.650   Median :12.00   Median :13.50   Median : 2.000  
##  Mean   : 5.896   Mean   :12.56   Mean   :17.02   Mean   : 5.105  
##  3rd Qu.: 6.880   3rd Qu.:14.00   3rd Qu.:26.00   3rd Qu.: 7.000  
##  Max.   :24.980   Max.   :18.00   Max.   :51.00   Max.   :44.000  
##     nonwhite          female          married           lwage        
##  Min.   :0.0000   Min.   :0.0000   Min.   :0.0000   Min.   :-0.6349  
##  1st Qu.:0.0000   1st Qu.:0.0000   1st Qu.:0.0000   1st Qu.: 1.2030  
##  Median :0.0000   Median :0.0000   Median :1.0000   Median : 1.5369  
##  Mean   :0.1027   Mean   :0.4791   Mean   :0.6084   Mean   : 1.6233  
##  3rd Qu.:0.0000   3rd Qu.:1.0000   3rd Qu.:1.0000   3rd Qu.: 1.9286  
##  Max.   :1.0000   Max.   :1.0000   Max.   :1.0000   Max.   : 3.2181  
##     expersq          tenursq       
##  Min.   :   1.0   Min.   :   0.00  
##  1st Qu.:  25.0   1st Qu.:   0.00  
##  Median : 182.5   Median :   4.00  
##  Mean   : 473.4   Mean   :  78.15  
##  3rd Qu.: 676.0   3rd Qu.:  49.00  
##  Max.   :2601.0   Max.   :1936.00

0.3.2 Correlaciones entre variables

Calculamos las correlaciones de los pares de variables, indicando que queremos trabajar con 3 decimales:

round(cor(datos, method = "pearson"), 3)
##              wage   educ  exper antiguedad nonwhite female married  lwage
## wage        1.000  0.406  0.113      0.347   -0.039 -0.340   0.229  0.937
## educ        0.406  1.000 -0.300     -0.056   -0.085 -0.085   0.069  0.431
## exper       0.113 -0.300  1.000      0.499    0.014 -0.042   0.317  0.111
## antiguedad  0.347 -0.056  0.499      1.000    0.012 -0.198   0.240  0.326
## nonwhite   -0.039 -0.085  0.014      0.012    1.000 -0.011  -0.062 -0.039
## female     -0.340 -0.085 -0.042     -0.198   -0.011  1.000  -0.166 -0.374
## married     0.229  0.069  0.317      0.240   -0.062 -0.166   1.000  0.271
## lwage       0.937  0.431  0.111      0.326   -0.039 -0.374   0.271  1.000
## expersq     0.030 -0.331  0.961      0.459    0.009 -0.028   0.217  0.023
## tenursq     0.267 -0.069  0.423      0.922   -0.007 -0.176   0.167  0.236
##            expersq tenursq
## wage         0.030   0.267
## educ        -0.331  -0.069
## exper        0.961   0.423
## antiguedad   0.459   0.922
## nonwhite     0.009  -0.007
## female      -0.028  -0.176
## married      0.217   0.167
## lwage        0.023   0.236
## expersq      1.000   0.414
## tenursq      0.414   1.000

Para que reconozca el nombre de las variables sin necesidad de citar el nombre del data frame con el signo de $

attach(datos)

Exploramos, mediante un gráfico en tres dimensiones, la relación entre el salario, los años de educación y la antiguedad en el cargo que tienen las personas. Se puede cambiar la perspectiva del gráfico con el puntero del mousse.

#install.packages("plotly")
library(plotly)
library(ggplot2)
plot_ly(x=educ, y=antiguedad, z=wage, type="scatter3d", color=wage) %>% 
  layout(scene = list(xaxis = list(title = 'educación (en años)'),
                      yaxis = list(title = 'antiguedad (en años)'),
                      zaxis = list(title = 'Salario (en USD/h)')))
## No scatter3d mode specifed:
##   Setting the mode to markers
##   Read more about this attribute -> https://plotly.com/r/reference/#scatter-mode
library(scatterplot3d)
attach(datos)
## The following objects are masked from datos (pos = 4):
## 
##     antiguedad, educ, exper, expersq, female, lwage, married, nonwhite,
##     tenursq, wage
graf=scatterplot3d(x=educ, y=antiguedad, z=wage, pch=16, cex.lab=1,
              highlight.3d=TRUE, type="h", xlab='Años de educación',
              ylab='Antiguedad (años)', zlab='Salario (USD/h)')

0.3.3 Modelo para el salario

Estimo un primer modelo en el que se explica el nivel del salario en función de los años de educación y de experiencia. El modelo a estimar es el siguiente: \[\begin{equation} wage_i = \beta _0+\beta _1 educ_i +\beta _2 antiguedad_i +u_i \\ \end{equation}\]

mod1 <- lm(wage ~ educ+ antiguedad, data=datos)
summary(mod1)
## 
## Call:
## lm(formula = wage ~ educ + antiguedad, data = datos)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -8.1438 -1.7288 -0.6372  1.2575 14.7482 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -2.22162    0.64015   -3.47 0.000563 ***
## educ         0.56914    0.04881   11.66  < 2e-16 ***
## antiguedad   0.18958    0.01871   10.13  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 3.092 on 523 degrees of freedom
## Multiple R-squared:  0.3019, Adjusted R-squared:  0.2992 
## F-statistic: 113.1 on 2 and 523 DF,  p-value: < 2.2e-16

Ajustamos el plano de regresión estimado.

graf=scatterplot3d(x=educ, y=antiguedad, z=wage, pch=16, cex.lab=1,
              highlight.3d=TRUE, type="h", xlab='Años de educación',
              ylab='Experiencia (años)', zlab='Salario (USD/h)')

graf$plane3d(mod1, lty.box="solid", col='mediumblue')

0.3.4 Análisis de los residuos

Calculamos la suma de cuadrados de los residuos

SCR <- deviance(mod1) #Suma del cuadrado de los residuos
gl <- df.residual(mod1) # Grados de libertad del modelo
sigma2_gorro <- (SCR/gl)
#otra manera es utilizar la función sigma que devuelve el error estándar de los residuos del modelo. 
sigma2_gorro <- (sigma(mod1))^2 
SCR
## [1] 4998.965
gl
## [1] 523
sigma2_gorro
## [1] 9.55825

También podemos utilizar la función anova que nos devuelve:

anova(mod1)
## Analysis of Variance Table
## 
## Response: wage
##             Df Sum Sq Mean Sq F value    Pr(>F)    
## educ         1 1179.7 1179.73  123.43 < 2.2e-16 ***
## antiguedad   1  981.7  981.72  102.71 < 2.2e-16 ***
## Residuals  523 4999.0    9.56                      
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Inspeccionamos la normalidad de los residuos mediante análisis gráficos

library(e1071)  # for skewness function

plot(density(mod1$residuals), main="Density Plot: residuos", ylab="Frequency", sub=paste("Skewness:", round(e1071::skewness(mod1$residuals), 2)))  # density para los 'residuos'

polygon(density(mod1$residuals), col="red")

Solicitamos los siguientes gráficos para inspeccionar la normalidad de los errores a través de la distribución empírica de los residuos (que es lo que observamos).
Permite comparar los cuantiles de los residuos con los cuantiles de una distribución normal. Si el proceso generador de datos fuera normal, los puntos estarían cerca de la recta.

qqnorm(mod1$residuals)
qqline(mod1$residuals)

Realizamos el contraste de normalidad:

shapiro.test(mod1$residuals)
## 
##  Shapiro-Wilk normality test
## 
## data:  mod1$residuals
## W = 0.88987, p-value < 2.2e-16

El resultado del p-valor indica que debemos rechazar la hipótesis nula de distribución normal de los errores. Esto confirma lo sugerido por el análisis gráfico precedente.

Graficamos los residuos contra cada uno de los regresores

library(ggplot2)
library(gridExtra)
plot1 <- ggplot(data = datos, aes(educ, mod1$residuals)) +
  geom_point() + geom_smooth(color = "firebrick") + geom_hline(yintercept = 0) +
  theme_bw()
plot2 <- ggplot(data = datos, aes(antiguedad, mod1$residuals)) +
  geom_point() + geom_smooth(color = "firebrick") + geom_hline(yintercept = 0) +
  theme_bw()
grid.arrange(plot1, plot2)
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'

El análisis gráfico indica problemas de heteroscedasticidad y de correlación entre los residuos y el nivel de los regresores.

Una manera resumida de verlo es graficando los residuos contra los valores ajustados de \(y\), que es una combinación lineal de los regresores que estamos considerando.

ggplot(data = datos, aes(mod1$fitted.values, mod1$residuals)) +
  geom_point() +
  geom_smooth(color = "firebrick", se = FALSE) +
  geom_hline(yintercept = 0) +
  theme_bw()
## `geom_smooth()` using method = 'loess' and formula 'y ~ x'

Se observan los mismos problemas.

###¿El género explica los salarios?

Nos preguntamos si el sexo biológico explica el nivel de salarios de las personas. Para responder esta pregunta comenzaremos calculando una tabla donde veremos algunas medidas resumenes del salario en función del sexo biológico:

library(tidyverse)
## -- Attaching packages --------------------------------------- tidyverse 1.3.1 --
## v tibble  3.1.6     v purrr   0.3.4
## v tidyr   1.2.0     v stringr 1.4.0
## v readr   2.1.2     v forcats 0.5.1
## -- Conflicts ------------------------------------------ tidyverse_conflicts() --
## x gridExtra::combine()       masks dplyr::combine()
## x plotly::filter()           masks dplyr::filter(), stats::filter()
## x ggstance::geom_errorbarh() masks ggplot2::geom_errorbarh()
## x dplyr::lag()               masks stats::lag()
tabla_salarios <- datos %>%
  group_by(female) %>%
  summarise("cantidad"=n(), "Mediana" = median(wage),  "Q1"=quantile(wage, 0.25),  "Q3"=quantile(wage, 0.75),"Promedio"= mean(wage), "Desvío"=sd(wage), "Mínimo"=min(wage),"Máximo"=max(wage))

print(tabla_salarios)
## # A tibble: 2 x 9
##   female cantidad Mediana    Q1    Q3 Promedio Desvío Mínimo Máximo
##    <int>    <int>   <dbl> <dbl> <dbl>    <dbl>  <dbl>  <dbl>  <dbl>
## 1      0      274    6     4.14  8.77     7.10   4.16  1.5     25.0
## 2      1      252    3.75  3     5.51     4.59   2.53  0.530   21.6

Observamos estadísticos de posición mayores para los hombres que para las mujeres.

Asignamos etiquetas a la variable female. La recodificación es en el mismo orden en el que está codificada la variable: 0 es hombre y 1 es mujer.

datos$género<- factor(datos$female, 
                         labels = c("hombre", "mujer"))

Observamos el diagrama de dispersión de los salarios por el género.

boxplot(wage ~ género, data = datos)

Estimamos un modelo de regresión múltiple agregando la variable female.

\[\begin{equation} wage_i = \beta _0+\beta _1 educ_i +\beta _2 antiguedad_i+ \beta _3 antiguedad_i +u_i \\ \end{equation}\]

mod2 <- lm(wage ~ educ + antiguedad+ female, data = datos)
summary(mod2)
## 
## Call:
## lm(formula = wage ~ educ + antiguedad + female, data = datos)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -7.5184 -1.8074 -0.4477  1.0270 14.1229 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) -0.84503    0.64774  -1.305    0.193    
## educ         0.53799    0.04709  11.425  < 2e-16 ***
## antiguedad   0.16441    0.01835   8.962  < 2e-16 ***
## female      -1.78839    0.26559  -6.734 4.38e-11 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 2.968 on 522 degrees of freedom
## Multiple R-squared:  0.3577, Adjusted R-squared:  0.354 
## F-statistic: 96.88 on 3 and 522 DF,  p-value: < 2.2e-16

Utilizamos el comando summ de la librería jtools, que proporciona varias opciones para dar formato a los resúmenes de regresión

library(jtools)
summ(mod2, scale = TRUE, vifs = TRUE, part.corr = TRUE, confint = TRUE, pvals = TRUE)
Observations 526
Dependent variable wage
Type OLS linear regression
F(3,522) 96.88
R² 0.36
Adj. R² 0.35
Est. 2.5% 97.5% t val. p VIF partial.r part.r
(Intercept) 6.75 6.40 7.11 37.21 0.00 NA NA NA
educ 1.49 1.23 1.75 11.43 0.00 1.01 0.45 0.40
antiguedad 1.19 0.93 1.45 8.96 0.00 1.05 0.37 0.31
female -1.79 -2.31 -1.27 -6.73 0.00 1.05 -0.28 -0.24
Standard errors: OLS; Continuous predictors are mean-centered and scaled by 1 s.d.

Para concatenar múltiples modelos en una sola tabla. https://jtools.jacob-long.com/reference/export_summs.html

#install.packages("huxtable")
library(huxtable)
## 
## Attaching package: 'huxtable'
## The following object is masked from 'package:ggplot2':
## 
##     theme_grey
## The following object is masked from 'package:dplyr':
## 
##     add_rownames
library(jtools)

mod0 <- lm(wage ~ educ, data = datos)
mod1 <- lm(wage ~ educ + antiguedad, data = datos)
mod2 <- lm(wage ~ educ + antiguedad + female, data = datos)
export_summs(mod0, mod1, mod2)
Model 1Model 2Model 3
(Intercept)-0.90    -2.22 ***-0.85    
(0.68)   (0.64)   (0.65)   
educ0.54 ***0.57 ***0.54 ***
(0.05)   (0.05)   (0.05)   
antiguedad       0.19 ***0.16 ***
       (0.02)   (0.02)   
female              -1.79 ***
              (0.27)   
N526       526       526       
R20.16    0.30    0.36    
*** p < 0.001; ** p < 0.01; * p < 0.05.

0.4 Variables ficticias o dummies

0.4.1 Introducción

Continuaremos trabajando con la base de salarios del texto de Wooldridge. El objeto de este capítulo es aprender a crear variables dummies y a utilizarlas en un MRLM.

Las variables que utilizaremos son las siguientes:

  • wage: salario promedio por hora
  • educ: años de educación
  • exper: años de experiencia potencial
  • ternure: años con el empleador actual (antiguedad)
  • nonwhite: =1 si no la persona no es blanca
  • female: =1 si la persona es mujer
  • married: =1 si la persona es casada
  • numdep: =1 número de dependientes que tiene la persona
  • lwage: log(salario)

0.4.2 Modelo para el salario

Continuaremos trabajando con la base de salarios de los temas anteriores (MRLS y MRLM)

library(wooldridge)

data("wage1")
names(wage1)
##  [1] "wage"     "educ"     "exper"    "tenure"   "nonwhite" "female"  
##  [7] "married"  "numdep"   "smsa"     "northcen" "south"    "west"    
## [13] "construc" "ndurman"  "trcommpu" "trade"    "services" "profserv"
## [19] "profocc"  "clerocc"  "servocc"  "lwage"    "expersq"  "tenursq"

\[\begin{equation} los(wage_i) = \beta _0+\beta _1 educ_i +\beta _2 exper_i+ \beta _3 antiguedad_i \beta _4 nro-dependiendientes_i+ +u_i \\ \end{equation}\]

wageModel <- lm(lwage ~ educ + exper + tenure +numdep , data = wage1)
summary(wageModel)
## 
## Call:
## lm(formula = lwage ~ educ + exper + tenure + numdep, data = wage1)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -2.04536 -0.29037 -0.03488  0.29086  1.42843 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 0.257392   0.112794   2.282   0.0229 *  
## educ        0.093192   0.007566  12.318  < 2e-16 ***
## exper       0.004259   0.001738   2.450   0.0146 *  
## tenure      0.022009   0.003097   7.107 3.92e-12 ***
## numdep      0.009872   0.015763   0.626   0.5314    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4411 on 521 degrees of freedom
## Multiple R-squared:  0.3165, Adjusted R-squared:  0.3113 
## F-statistic: 60.32 on 4 and 521 DF,  p-value: < 2.2e-16
anova(wageModel)
DfSum SqMean SqF valuePr(>F)
127.6   27.6   142    4.67e-29
19.42  9.42  48.4  1.03e-11
19.89  9.89  50.8  3.39e-12
10.07630.07630.3920.531   
521101     0.195            

0.4.3 Crear variables dummies

Para crear dummies se puede utilizar el paquete fastDummies. Creamos una dummie para cada valor de la variable numdep. Crea una dummie por cada categoría, en este caso son 7.

#install.packages("fastDummies")
library(fastDummies)

datos1<-dummy_cols(wage1, select_columns="numdep")
names(datos1)
##  [1] "wage"     "educ"     "exper"    "tenure"   "nonwhite" "female"  
##  [7] "married"  "numdep"   "smsa"     "northcen" "south"    "west"    
## [13] "construc" "ndurman"  "trcommpu" "trade"    "services" "profserv"
## [19] "profocc"  "clerocc"  "servocc"  "lwage"    "expersq"  "tenursq" 
## [25] "numdep_0" "numdep_1" "numdep_2" "numdep_3" "numdep_4" "numdep_5"
## [31] "numdep_6"

Para crear las dummies para mas de una variable las listo con c() donde listo el nombre de las variables separados por coma.

#dummy_cols(datos, select_columns=c("numdep", "exper"))

0.4.4 Crear interacciones entre los regresores

Como ejemplo se van crear las interacciones entre estado civil y género.

Se crean las dummies variables

  • mujer casada: =1 si la persona es mujer casada
  • hombre casado: =1 si la persona es hombre casado
  • mujer soltera: =1 si la persona es mujer soltera
  • hombre soltero=1 si la persona es hombre soltero

Para ello es necesario crea antes la variable hombre y la variable soltero que no están en la base como tal. Utilizaremos la librería car:

# install.packages("car")
library(car)
## Loading required package: carData
## 
## Attaching package: 'car'
## The following object is masked from 'package:purrr':
## 
##     some
## The following object is masked from 'package:dplyr':
## 
##     recode
datos1$male <- recode(datos1$female, "1=0; 0=1")
datos1$single <- recode(datos1$married, "1=0; 0=1")

names(datos1)
##  [1] "wage"     "educ"     "exper"    "tenure"   "nonwhite" "female"  
##  [7] "married"  "numdep"   "smsa"     "northcen" "south"    "west"    
## [13] "construc" "ndurman"  "trcommpu" "trade"    "services" "profserv"
## [19] "profocc"  "clerocc"  "servocc"  "lwage"    "expersq"  "tenursq" 
## [25] "numdep_0" "numdep_1" "numdep_2" "numdep_3" "numdep_4" "numdep_5"
## [31] "numdep_6" "male"     "single"
attach(datos1)
## The following objects are masked from datos (pos = 14):
## 
##     educ, exper, expersq, female, lwage, married, nonwhite, tenursq,
##     wage
## The following objects are masked from datos (pos = 16):
## 
##     educ, exper, expersq, female, lwage, married, nonwhite, tenursq,
##     wage
marrmale<-married*male
marrfem<-married*female
singfem<-single*female

datos1<-cbind(datos1, marrmale, marrfem, singfem)
names(datos1)
##  [1] "wage"     "educ"     "exper"    "tenure"   "nonwhite" "female"  
##  [7] "married"  "numdep"   "smsa"     "northcen" "south"    "west"    
## [13] "construc" "ndurman"  "trcommpu" "trade"    "services" "profserv"
## [19] "profocc"  "clerocc"  "servocc"  "lwage"    "expersq"  "tenursq" 
## [25] "numdep_0" "numdep_1" "numdep_2" "numdep_3" "numdep_4" "numdep_5"
## [31] "numdep_6" "male"     "single"   "marrmale" "marrfem"  "singfem"

0.4.5 Estimamos modelos para el salario

El modelo mas básico estima el log del salario en función de la educación de las personas. \[\begin{equation} log(salario_i) = \beta _0+\beta _1 Educ_i + u_i \\ \end{equation}\]

modelo1 <- lm(lwage ~ educ, data = datos1)
summary(modelo1)
## 
## Call:
## lm(formula = lwage ~ educ, data = datos1)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -2.21158 -0.36393 -0.07263  0.29712  1.52339 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 0.583773   0.097336   5.998 3.74e-09 ***
## educ        0.082744   0.007567  10.935  < 2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4801 on 524 degrees of freedom
## Multiple R-squared:  0.1858, Adjusted R-squared:  0.1843 
## F-statistic: 119.6 on 1 and 524 DF,  p-value: < 2.2e-16

Nos preguntamos si el efecto de la educación es el mismo para hombres y para mujeres, para ello corremos los modelos considerando dos bases separadas por sexo (mujeres y hombres)

library(dplyr)
mujeres<- datos1%>% 
  select("lwage", "educ", "exper", "tenure") %>% 
  filter(female==1)

hombres<- datos1%>% 
  select("lwage", "educ", "exper", "tenure") %>% 
  filter(female==0)

dim(mujeres)
## [1] 252   4
dim(hombres)
## [1] 274   4

Corremos el mismo modelo por un lado con las mujeres y por otro con los hombres:

\[\begin{equation} log(salario_mi) = \beta _0+\beta _1 Educ_mi + u_mi \\ \end{equation}\]

\[\begin{equation} log(salario_hi) = \beta _0+\beta _1 Educ_hi + u_hi \\ \end{equation}\]

sal_mujeres<-lm(lwage~ educ, data = mujeres)
sal_hombres<-lm(lwage~ educ, data = hombres)

summary(sal_mujeres)
## 
## Call:
## lm(formula = lwage ~ educ, data = mujeres)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -2.02673 -0.24397 -0.06163  0.21415  1.21924 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.46589    0.12890   3.614 0.000364 ***
## educ         0.07716    0.01026   7.520 9.82e-13 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.402 on 250 degrees of freedom
## Multiple R-squared:  0.1845, Adjusted R-squared:  0.1812 
## F-statistic: 56.55 on 1 and 250 DF,  p-value: 9.824e-13
summary(sal_hombres)
## 
## Call:
## lm(formula = lwage ~ educ, data = hombres)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.11585 -0.34240 -0.01708  0.32659  1.34740 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 0.825955   0.127813   6.462 4.75e-10 ***
## educ        0.077228   0.009731   7.936 5.47e-14 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4828 on 272 degrees of freedom
## Multiple R-squared:  0.188,  Adjusted R-squared:  0.185 
## F-statistic: 62.99 on 1 and 272 DF,  p-value: 5.471e-14
library(jtools)
library(huxtable)

export_summs(modelo1, sal_mujeres, sal_hombres, 
              modelo.nombres = c("General","Mujeres","Hombres"))
Model 1Model 2Model 3
(Intercept)0.58 ***0.47 ***0.83 ***
(0.10)   (0.13)   (0.13)   
educ0.08 ***0.08 ***0.08 ***
(0.01)   (0.01)   (0.01)   
N526       252       274       
R20.19    0.18    0.19    
*** p < 0.001; ** p < 0.01; * p < 0.05.

El retorno de la educación parece ser el mismo para ambos géneros. Eso no significa que el salario sea el mismo.

Etiquetamos las categorías para los gráficos:

datos1$Sexo_biológico<- factor(datos1$female, 
                       labels = c("hombre", "mujer"))

Graficamos los diagramas de dispersión del salario con respecto al sexo biológico para ver si existen diferencias:

library(ggplot2)
ggplot(datos1, aes(x=educ, y=lwage, color=Sexo_biológico)) + 
  geom_point() +
  geom_smooth(method='lm', formula=y~x, se=FALSE, col='dodgerblue1') +
      theme_light() +
  facet_wrap(~Sexo_biológico)+
     ggtitle("Relación entre salario y educación según el sexo biológico")

Para contrastar si el retorno de la educación es el mismo según el sexo biológico estimamos el siguiente modelo: \[\begin{equation} log(salario_i) = \beta _0+\beta _1 Mujer_i + u_i \\ \end{equation}\]

model_sexo<-lm(lwage~ female, data = datos1)
summary(model_sexo)
## 
## Call:
## lm(formula = lwage ~ female, data = datos1)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -2.05123 -0.31774 -0.04889  0.35548  1.65773 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  1.81357    0.02981  60.830   <2e-16 ***
## female      -0.39722    0.04307  -9.222   <2e-16 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4935 on 524 degrees of freedom
## Multiple R-squared:  0.1396, Adjusted R-squared:  0.138 
## F-statistic: 85.04 on 1 and 524 DF,  p-value: < 2.2e-16
boxplot(wage ~ Sexo_biológico, data = datos1)

library(tidyverse)
tabla_salarios <- datos1 %>%
  group_by(female) %>%
  summarise("cantidad"=n(), "Promedio"= mean(lwage))

print(tabla_salarios)
## # A tibble: 2 x 3
##   female cantidad Promedio
##    <int>    <int>    <dbl>
## 1      0      274     1.81
## 2      1      252     1.42

La diferencia en los valores promedio es igual a 1.416353-1.813570= -0.397217, que coincide con la estimación del coeficiente de regresión simple del logaritmo del salario contra el sexo biológico (variable female).

0.4.6 Modelos para el salario que consideran interacciones entre regresores

Estimaremos un modelo con interacciones para dar la posibilidad de que el retorno de la educación sea diferente para hombres que para mujeres.

\[\begin{equation} log(salario_i) = \beta _0+\beta _1 Female_i + \beta _2 Educ_i+\beta _3 Female_i*Educ_i+u_i \\ \end{equation}\]

modelo4<-lm(lwage ~ female+educ+female*educ, data=datos1)
summary(modelo4)
## 
## Call:
## lm(formula = lwage ~ female + educ + female * educ, data = datos1)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -2.02673 -0.27468 -0.03721  0.26221  1.34740 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  8.260e-01  1.181e-01   6.997 8.08e-12 ***
## female      -3.601e-01  1.854e-01  -1.942   0.0527 .  
## educ         7.723e-02  8.988e-03   8.593  < 2e-16 ***
## female:educ -6.408e-05  1.450e-02  -0.004   0.9965    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4459 on 522 degrees of freedom
## Multiple R-squared:  0.3002, Adjusted R-squared:  0.2962 
## F-statistic: 74.65 on 3 and 522 DF,  p-value: < 2.2e-16

Se observa que la interacción no es significativa, esto es equivalente al resultado visto antes que mostraba que el rendimiento de la educación era el mismo para mujeres que para hombres.

Estimamos un modelo para el logaritmo del salario, en función del sexo biológico, la educación, la experiencia, la antiguedad y sus efectos cuadráticos. \[\begin{equation} log(salario_i) = \beta _0+\beta _1 Female_i + \beta _2 Educ_i+ \beta _3 exper_i+ \beta _4 expersq_i+\beta _5 tenure_i+\beta _6 tenuresq_i+ \beta _7 Female_i*Educ_i+u_i \\ \end{equation}\]

modelo5 <- lm(lwage ~ female+ educ+ exper+ expersq+ tenure  + tenursq + educ*female, data = datos1)
summary(modelo5)
## 
## Call:
## lm(formula = lwage ~ female + educ + exper + expersq + tenure + 
##     tenursq + educ * female, data = datos1)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.83265 -0.25261 -0.02374  0.25396  1.13584 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.3888060  0.1186871   3.276  0.00112 ** 
## female      -0.2267886  0.1675394  -1.354  0.17644    
## educ         0.0823692  0.0084699   9.725  < 2e-16 ***
## exper        0.0293366  0.0049842   5.886 7.11e-09 ***
## expersq     -0.0005804  0.0001075  -5.398 1.03e-07 ***
## tenure       0.0318967  0.0068640   4.647 4.28e-06 ***
## tenursq     -0.0005900  0.0002352  -2.509  0.01242 *  
## female:educ -0.0055645  0.0130618  -0.426  0.67028    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.4001 on 518 degrees of freedom
## Multiple R-squared:  0.441,  Adjusted R-squared:  0.4334 
## F-statistic: 58.37 on 7 and 518 DF,  p-value: < 2.2e-16

¿Cuál es el rendimiento estimado de la educación para los hombres? ¿Y para las mujeres?

El rendimiento estimado de la educación para los hombres es \(\widehat{\beta_2}=0.082\), es decir, \(8.2\%\) Para las mujeres es (\(\widehat{\beta_2}-\widehat{\beta_7}=0.082-0.0056=0.0764\)), aprox \(7.6\%\)

La diferencia en el rendimiento de la educación para ambos sexos se estima en poco mas de medio punto porcentual (\(-0.56\%\)) Como vimos, esta diferencia no es grande desde el punto de vista económico ni tampoco es estadísticamente significativa (observar el p-valor) En consecuencia, no se ha encontrado suficiente evidencia empírica para rechazar la hipótesis de que el rendimiento en la educación sea el mismo para hombres y para mujeres.

Estimamos un modelo para un último modelo para el salario incluyendo las interacciones entre el sexo biológico y el estado civil.

names(datos1)

\[\begin{equation} log(salario_i) = \beta _0+\beta _1 marrmale_i + \beta_2 marrfem_i+ \beta _3 singfem_i+ \beta _4 educ_i+ \beta _5 exper_i+\beta _6 expersq_i+\beta _7 tenure_i+\beta _8 tenuresq_i+ +u_i \\ \end{equation}\]

modelo6 <- lm(lwage ~ marrmale+ marrfem+ singfem+ educ+exper+expersq+ tenure  + tenursq, data = datos1)
summary(modelo6)
## 
## Call:
## lm(formula = lwage ~ marrmale + marrfem + singfem + educ + exper + 
##     expersq + tenure + tenursq, data = datos1)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -1.89697 -0.24060 -0.02689  0.23144  1.09197 
## 
## Coefficients:
##               Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  0.3213781  0.1000090   3.213 0.001393 ** 
## marrmale     0.2126757  0.0553572   3.842 0.000137 ***
## marrfem     -0.1982676  0.0578355  -3.428 0.000656 ***
## singfem     -0.1103502  0.0557421  -1.980 0.048272 *  
## educ         0.0789103  0.0066945  11.787  < 2e-16 ***
## exper        0.0268006  0.0052428   5.112 4.50e-07 ***
## expersq     -0.0005352  0.0001104  -4.847 1.66e-06 ***
## tenure       0.0290875  0.0067620   4.302 2.03e-05 ***
## tenursq     -0.0005331  0.0002312  -2.306 0.021531 *  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.3933 on 517 degrees of freedom
## Multiple R-squared:  0.4609, Adjusted R-squared:  0.4525 
## F-statistic: 55.25 on 8 and 517 DF,  p-value: < 2.2e-16

La categoría de referencia en este caso son los hombres solteros, por lo tanto la estimación de los coeficientes de las dummies debe interpretarse como diferencial respecto a esta categoría.\

Manteniendo constantes los niveles de educación, experiencia y antiguedad (tenure), se estima que los hombres casados, en promedio, ganan aproximadamente 21,3% más que los solteros Se estima que, a igualdad del resto de las características, una mujer casada, gana en promedio 19,8% menos que un hombre soltero Se estima que, a igualdad del resto de las características, una mujer soltera, gana en promedio 11% menos que un hombre soltero.