Tarea 3

Respuestas

Respuestas

Pregunta 1

En este ejercicio usaremos los datos de Blackburn & Neumark (1992), que pueden accederse instalando el paquete wooldridge y luego llamando el objeto wooldridge::wage2. La pregunta de investigación de este trabajo es la misma que la del estudio de Card (1993) visto en clase, pero con otra estrategia de identificación. El propósito de este ejercicio es estimar los retornos a la educación.

  1. [3 puntos] Estime una regresión por MCO para explicar el log del salario (lwage) en función de la educación educ y los siguientes controles: experiencia (exper), tiempo en el empleo actual (tenure), edad (age), casado (married), número de hermanos (sibs) y urbano (urban). Reporte errores clásicos y errores robustos. ¿Qué problema encuentra en la estimación de esta relación? ¿El coeficiente sobre educ tiene una interpretación causal del efecto de la educación en el salario?

    Estimamos por MCO la relación entre salarios y educación, controlando por un conjunto de regresores:

    library(wooldridge)
    
    db <- wooldridge::wage2
    
    regmco <- lm(lwage ~ educ + exper + tenure + age + married + sibs + urban,
                data = db)
    
    modelsummary(models = list(regmco, regmco),
              output = 'gt',
              stars = c('***' = .01, '**' = .05, '*' = 0.1),
              vcov = c("classical", "HC3"),
              fmt = 5,
              coef_map = 'educ',
              gof_map = "nobs")
    (1) (2)
    educ 0.06730*** 0.06730***
    (0.00668) (0.00676)
    Num.Obs. 935 935
    * p < 0.1, ** p < 0.05, *** p < 0.01

    Hay una relación de 6.7% mayor ingreso por cada año de educación adicional. Sin embargo, esta no es una relación causal pues es muy probable que la educación no sea exógena en la ecuación de salarios. Esto puede deberse, por ejemplo, a una variable omitida de habilidad que afecta tanto al número de años de educación alcanzados como al desempeño en el mercado de trabajo.

  2. [3 puntos] Se propone usar como instrumento de los años de educación a la educación del padre. ¿Qué condiciones debe cumplir la variable propuesta para funcionar como instrumento válido?

    El instrumento debe cumplir dos condiciones:

    Exogeneidad: el instrumento no debe pertenecer a la ecuación de salarios. Es decir, la educación del padre no debe afectar el salario contemporáneo de forma directa.

    Relevancia: el instrumento debe estar correlacionado con la variable endógena. En este caso, la educación del padre debe estar correlacionada con el número de años de educación completados.

  3. [4 puntos] ¿Cómo juzga la propuesta de usar la variable antes descrita como instrumento?

    Este argumento fue usado por Blackburn y Neumark (1992) para estudiar los diferneciales en los salarios entre industrias en los Estados Unidos.

    No hay una respuesta correcta o incorrecta. Quiero leer sus argumentos.

  4. [3 puntos] Estime la relación entre el logaritmo del salario y la educación usando la educación del padre, feduc, como instrumento. Emplee las mismas variables de control que en el modelo de MCO. Reporte errores clásicos y errores robustos.

    regvi <- ivreg(lwage ~ educ + exper + tenure + age + married + sibs + urban  |
                     feduc + exper + tenure + age + married + sibs + urban,
                   data=db)
    
    modelsummary(models = list(regmco, regmco,regvi, regvi),
              output = 'gt',
              stars = c('***' = .01, '**' = .05, '*' = 0.1),
              vcov = c("classical", "HC3", "classical", "HC3"),
              fmt = 5,
              coef_map = 'educ',
              gof_map = "nobs")
    (1) (2) (3) (4)
    educ 0.06730*** 0.06730*** 0.13283*** 0.13283***
    (0.00668) (0.00676) (0.02413) (0.02411)
    Num.Obs. 935 935 741 741
    * p < 0.1, ** p < 0.05, *** p < 0.01
  5. [4 puntos] Interprete la primera etapa en términos del coeficiente sobre el instrumento. Obtenga el estadístico \(F\) del instrumento excluido e interprete su magnitud.

    En la primera etapa, cada año de educación del padre se asocia con 0.19 años de escolaridad acumulados. Este efecto estadísticamente significativo al 1%. El estadístico F es de una magnitud de 85, por encima de 10, la regla de dedo comúnmente empleada para juzgar la presencia de instrumentos débiles.

    regpe <- lm(educ ~ feduc + exper + tenure + age + married + sibs + urban,
                data=db)
    
    
    modelsummary(models = list(regpe),
              output = 'gt',
              stars = c('***' = .01, '**' = .05, '*' = 0.1),
              vcov = c("HC3"),
              fmt = 5,
              coef_map = 'feduc',
              gof_map = "nobs")
    (1)
    feduc 0.19398***
    (0.02093)
    Num.Obs. 741
    * p < 0.1, ** p < 0.05, *** p < 0.01
    linearHypothesis(regpe, "feduc=0")
    Linear hypothesis test
    
    Hypothesis:
    feduc = 0
    
    Model 1: restricted model
    Model 2: educ ~ feduc + exper + tenure + age + married + sibs + urban
    
      Res.Df    RSS Df Sum of Sq      F    Pr(>F)    
    1    734 2543.1                                  
    2    733 2278.2  1    264.91 85.235 < 2.2e-16 ***
    ---
    Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
  6. [3 puntos] Interprete el coeficiente sobre la variable de educación en el modelo estructural. Compare la magnitud del efecto estimado con el resultado de MCO.

    El coeficiente estimado sobre los años de educación indica que un año adicional de escolaridad incrementa en 13.3% el salario. Este efecto es casi el doble del estimado por MCO y estadísticamente significativo al 1%.

    modelsummary(models = list(regmco, regvi),
              output = 'gt',
              stars = c('***' = .01, '**' = .05, '*' = 0.1),
              vcov = c("HC3", "HC3"),
              fmt = 5,
              coef_map = 'educ',
              gof_map = "nobs")
    (1) (2)
    educ 0.06730*** 0.13283***
    (0.00676) (0.02411)
    Num.Obs. 935 741
    * p < 0.1, ** p < 0.05, *** p < 0.01
  7. [3 puntos] Realice ahora el siguiente procedimiento. Primero, estime la primera etapa usando una regresión por MCO. Obtenga los valores ajustados de educación y llámelos educ_hat. Luego, estime la segunda etapa empleando educ_hat como variable independiente, además del resto de variables de control. ¿Cómo cambian sus resultados en comparación con la parte d.?

    La magnitud de los coeficientes estimados es la misma. Esto es lo que esperábamos pues sabemos que el estimador de MC2E puede entenderse como un procedimiento donde primero se estiman los valores ajustados de la variable endógena usando el instrumento y las variables de control y luego se usan estos valores ajustados en la ecuación estructural. En cambio, los errores estándar son algo distintos.

    db <- db %>% 
      mutate(educ_hat = predict(regpe, .))
    
    reg2e <- lm(lwage ~ educ_hat + exper + tenure + age + married + sibs + urban,
                data=db)
    
    #Comparamos
    
    modelsummary(models = list(regvi, reg2e),
              output = 'gt',
              stars = c('***' = .01, '**' = .05, '*' = 0.1),
              vcov = c("HC3", "HC3"),
              fmt = 5,
              coef_map = c('educ', 'educ_hat'),
              gof_map = "nobs")
    (1) (2)
    educ 0.13283***
    (0.02411)
    educ_hat 0.13283***
    (0.02417)
    Num.Obs. 741 741
    * p < 0.1, ** p < 0.05, *** p < 0.01
  8. [3 puntos] ¿A qué se deben las discrepancias que encuentra? ¿Cuál de las dos estrategias prefiere para estimar el modelo de variables instrumentales?

    Los coeficientes estimados son exactamente iguales, pero los errores estándar no. El problema es que nuestro procedimiento de MC2E a mano no toma en cuenta que en la ecuación estructural estamos usando valores ajustados de la variable endógena. Las funciones en la mayoría de los paquetes utilizados en econometría calculan los errores estándar de manera correcta. Preferimos usar las funciones previamente ya programadas cuando sea posible, aunque este ejercicio nos ayuda a reforzar la intuición del estimador de MC2E.

  9. [3 puntos] Reestime el modelo de variables instrumentales añadiendo un segundo instrumento, la educación de la madre, meduc, y reporte errores robustos (no es necesario usar gmm, puede seguir con ivreg, por lo que no estaría obteniendo el estimador de MGM óptimo). ¿Cómo cambian sus resultados para la ecuación estructural con respecto al caso exactamente identificado?

    El efecto estimado es significativo al 1% y en magnitud se incrementa ligeramente hasta 0.137.

    regvi2 <- ivreg(lwage ~ educ + exper + tenure + age + married + sibs + urban  |
                     feduc + meduc +  exper + tenure + age + married + sibs + urban,
                   data=db)
    
    
    modelsummary(models = list(regmco, regvi, regvi2),
              output = 'gt',
              stars = c('***' = .01, '**' = .05, '*' = 0.1),
              vcov = c("HC3", "HC3", "HC3"),
              fmt = 5,
              coef_map = c('educ'),
              gof_map = "nobs")
    (1) (2) (3)
    educ 0.06730*** 0.13283*** 0.13719***
    (0.00676) (0.02411) (0.02198)
    Num.Obs. 935 741 722
    * p < 0.1, ** p < 0.05, *** p < 0.01
  10. [3 puntos] Con el objeto que resulta de la estimación del modelo sobreidentificado, realice summary(OBJETO, vcov = sandwich, diagnostics = TRUE) para obtener las tres pruebas diagnóstico más usadas en variables instrumentales: prueba de instrumentos débiles, prueba de Hausman y prueba de Sargan. Interprete cada una de las pruebas.

    summary(regvi2, vcov = sandwich, diagnostics = TRUE)
    
    Call:
    ivreg(formula = lwage ~ educ + exper + tenure + age + married + 
        sibs + urban | feduc + meduc + exper + tenure + age + married + 
        sibs + urban, data = db)
    
    Residuals:
         Min       1Q   Median       3Q      Max 
    -1.82265 -0.22364  0.01904  0.24623  1.29960 
    
    Coefficients:
                 Estimate Std. Error t value Pr(>|t|)    
    (Intercept)  4.418980   0.268857  16.436  < 2e-16 ***
    educ         0.137194   0.021691   6.325 4.47e-10 ***
    exper        0.035400   0.007952   4.452 9.88e-06 ***
    tenure       0.007832   0.003069   2.552   0.0109 *  
    age         -0.008260   0.007499  -1.102   0.2710    
    married      0.210486   0.048810   4.312 1.84e-05 ***
    sibs         0.003731   0.007297   0.511   0.6093    
    urban        0.170560   0.032662   5.222 2.32e-07 ***
    
    Diagnostic tests:
                     df1 df2 statistic  p-value    
    Weak instruments   2 713    51.717  < 2e-16 ***
    Wu-Hausman         1 713    14.588 0.000145 ***
    Sargan             1  NA     0.342 0.558820    
    ---
    Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
    
    Residual standard error: 0.3949 on 714 degrees of freedom
    Multiple R-Squared: 0.1221,  Adjusted R-squared: 0.1134 
    Wald test: 25.22 on 7 and 714 DF,  p-value: < 2.2e-16 

    La prueba de instrumentos débiles rechaza la \(H_0\) de que los instrumentos son débiles, por lo que tenemos primera etapa.

    La prueba de Hausman rechaza al 10% que los estimadores de VI y de MCO sean iguales, por lo que se prefiere el de VI.

    La prueba de Sargan no rechaza la \(H_0\) de que el modelo esté mal especificado.

  11. [4 puntos] Considere la primera etapa del modelo sobreidentificado. Compruebe que si realiza una prueba de significancia conjunta para los instrumentos obtiene la prueba de instrumentos débiles que se reporta en el resumen que obtuvo con summary.

    En el apartado anterior, obtenemos un valor \(p<0.0000000000000002\) para la prueba de instrumentos débiles. Noten que se especificó una matriz de varianzas robustas. Entonces, tenemos que usar la misma matriz en la prueba \(F\):

    regpe2 <- lm(educ ~ feduc + meduc + exper + tenure + age + married + sibs + urban,
                data=db)
    
    linearHypothesis(regpe2, c("feduc=0", "meduc=0"), white.adjust = "hc0")
    Linear hypothesis test
    
    Hypothesis:
    feduc = 0
    meduc = 0
    
    Model 1: restricted model
    Model 2: educ ~ feduc + meduc + exper + tenure + age + married + sibs + 
        urban
    
    Note: Coefficient covariance matrix supplied.
    
      Res.Df Df      F    Pr(>F)    
    1    715                        
    2    713  2 51.717 < 2.2e-16 ***
    ---
    Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

    Obtenemos el mismo valor \(p\) y maginitud del estadístico \(F\).

  12. [4 puntos] Compruebe que si realiza el procedimiento de regresión auxiliar para la prueba de Hausman obtiene el mismo valor \(p\) que se reporta en el resumen que obtuvo con summary.

    De la primera etapa, obtenemos los residuales:

    db <- db %>% 
      filter(!is.na(lwage)) %>% 
      mutate(vhat = educ - predict(regpe2, newdata=db))

    Corremos la regresión auxiliar con los residuales:

    regaux <- lm(lwage ~ educ + exper + tenure + age + married + sibs + urban + vhat,
                   data=db)
    
    
    modelsummary(models = list(regaux),
              output = 'gt',
              stars = c('***' = .01, '**' = .05, '*' = 0.1),
              vcov = c("HC0"),
              statistic = c("p-val.={p.value}, e.e.={std.error}"),
              fmt = 6,
              coef_map = c('vhat'),
              gof_map = "nobs")
    (1)
    vhat -0.080976***
    p-val.=0.000145, e.e.=0.021201
    Num.Obs. 722
    * p < 0.1, ** p < 0.05, *** p < 0.01

    El valor \(p\) asociado al coeficiente de los residuales es idéntico al obtenido para la prueba de Hausman con summary.

Pregunta 2

Considere los datos comportamiento_wide.csv, que contienen información individual de niñas y niños, incluyendo su género, edad, raza e información de sus madres. Además, se incluye una medida auto reportada de autoestima (self) y una evaluación de comportamiento antisocial (anti). Se quiere conocer cómo influye la autoestima en el comportamiento antisocial. Para cada niño o niña hay tres observaciones en el tiempo. Se busca explicar el comportamiento antisocial en función de la autoestima y la condición de pobreza (pov):

\[anti_{it}=\alpha_i+\beta_1 self_{it}+\beta_2 pov_{it}+\varepsilon_{it}\]

  1. [2 puntos] La base se encuentra en formato wide. Ponga la base en formato long, donde haya una columna para cada variable y donde las filas representen a un individuo en un periodo.

    Hay muchas formas de hacer esto. Podemos usar las funciones pivot_longer y pivot_wider, por ejemplo.

    data.comp <-read_csv("../files/comportamiento_wide.csv",
                          locale = locale(encoding = "latin1")) %>%
      pivot_longer(c(anti90:anti94,self90:self94,pov90:pov94),
                   names_to = c("measure", "year"),
                   names_pattern = "(.*)(..)")  %>%
      pivot_wider(names_from = measure,
                  values_from = value)
    
    colnames(data.comp)
     [1] "id"       "momage"   "gender"   "childage" "hispanic" "black"   
     [7] "momwork"  "married"  "year"     "anti"     "self"     "pov"     
  2. [2 puntos] Estime la ecuación de comportamiento antisocial empleando MCO pooled. ¿Cuáles son los supuestos que se deben cumplir para que \(\hat{\beta}_1^{MCO}\) sea consistente?

    m.mco <- plm( anti ~ self + pov,
                          data=data.comp,
                          model="pooling",
                          index = c("id", "year"))
    
    
    
    modelsummary(models = list("MCO"=m.mco),
              output = 'gt',
              stars = c('***' = .01, '**' = .05, '*' = 0.1),
              vcov = c("iid"),
              fmt = 5,
              coef_map = c('self', 'pov'),
              gof_map = "nobs")
    MCO
    self -0.06510***
    (0.01108)
    pov 0.51581***
    (0.07873)
    Num.Obs. 1743
    * p < 0.1, ** p < 0.05, *** p < 0.01

    La variable self tiene un efecto negativo y estadísticamente significativo sobre anti. La variable pov tiene un efecto positivo y estadísticamente significativo. El estimador de MCO será consistente solo si las variables self y pov no están correlacionadas con el error. Además, para estimar este modelo, asumimos que la heterogeneidad no observada \(\alpha_i\) puede escribirse simplemente como \(\alpha\). Otra forma de pensar sobre este modelo es si el mismo modelo es válido para todos los periodos como para asumir una ordenada al origen y una pendiente común.

    El modelo pooled ignora la naturaleza en panel de los datos. Sin embargo, como tenemos a los mismos individuos en varios puntos del tiempo, los errores están agrupados, así que se deben de estimar errores con esta estructura. En este caso, al tomar en cuenta esta correlación entre individuos, los errores estándar son más grandes, pero los resultados siguen siendo significativos. En muchos casos, no tomar en cuenta la estructura agrupada de los errores puede llevar a rechazar hipótesis nulas que son ciertas.

    modelsummary(models = list("MCO"=m.mco,
                               "MCO, errores agrupados"=m.mco),
              output = 'gt',
              stars = c('***' = .01, '**' = .05, '*' = 0.1),
              vcov=list(NULL, clubSandwich::vcovCR(m.mco, type='CR1', cluster=data.comp$id)),
              fmt = 5,
              coef_map = c('self', 'pov'),
              gof_map = "nobs")
    MCO MCO, errores agrupados
    self -0.06510*** -0.06510***
    (0.01108) (0.01369)
    pov 0.51581*** 0.51581***
    (0.07873) (0.10496)
    Num.Obs. 1743 1743
    * p < 0.1, ** p < 0.05, *** p < 0.01
  3. [3 puntos] Estime la ecuación de comportamiento antisocial empleando el estimador within. ¿Cuáles son los supuestos que se deben cumplir para que \(\hat{\beta}_1^{FE}\) sea consistente?

    Si asumimos que la heterogeneidad no observada y el error están potencialmente correlacionados, entonces podemos usar un estimador de efectos fijos para deshacernos de la heterogeneidad no observada y estimar consistentemente los parámetros sobre self y pov.

    m.fe <- plm( anti ~ self + pov,
                 data=data.comp,
                 model="within",
                 index = c("id", "year"))
    
    modelsummary(models = list("MCO"=m.mco,
                               "MCO, errores agrupados"=m.mco,
                               "Efectos fijos"=m.fe),
              output = 'gt',
              stars = c('***' = .01, '**' = .05, '*' = 0.1),
              vcov=list(NULL, clubSandwich::vcovCR(m.mco, type='CR1', cluster=data.comp$id), clubSandwich::vcovCR(m.fe, type='CR1', cluster=data.comp$id)),
              fmt = 5,
              coef_map = c('self', 'pov'),
              gof_map = "nobs")
    MCO MCO, errores agrupados Efectos fijos
    self -0.06510*** -0.06510*** -0.05150***
    (0.01108) (0.01369) (0.01131)
    pov 0.51581*** 0.51581*** 0.10490
    (0.07873) (0.10496) (0.09922)
    Num.Obs. 1743 1743 1743
    * p < 0.1, ** p < 0.05, *** p < 0.01
  4. [3 puntos] Estime la ecuación de comportamiento antisocial empleando efectos aleatorios. ¿Cuáles son los supuestos que se deben cumplir para que \(\hat{\beta}_1^{RE}\) sea consistente?

    Si estamos dispuestos a asumir que la heterogeneidad no observada y el error son independientes, podemos emplear el estimador de efectos aleatorios. MCO pooled también es consistente pero no es eficiente.

    m.re <- plm( anti ~ self + pov,
                 data=data.comp,
                 model="random",
                 index = c("id", "year"))
    
    
    modelsummary(models = list("MCO"=m.mco,
                               "MCO, errores agrupados"=m.mco,
                               "Efectos fijos"=m.fe,
                               "Efectos aleatorios"=m.re),
              output = 'gt',
              stars = c('***' = .01, '**' = .05, '*' = 0.1),
              vcov=list(NULL,
                        clubSandwich::vcovCR(m.mco, type='CR1',cluster=data.comp$id),
                        clubSandwich::vcovCR(m.fe, type='CR1', cluster=data.comp$id),
                        clubSandwich::vcovCR(m.re, type='CR1', cluster=data.comp$id)),
              fmt = 5,
              coef_map = c('self', 'pov'),
              gof_map = "nobs")
    MCO MCO, errores agrupados Efectos fijos Efectos aleatorios
    self -0.06510*** -0.06510*** -0.05150*** -0.05673***
    (0.01108) (0.01369) (0.01131) (0.01022)
    pov 0.51581*** 0.51581*** 0.10490 0.29241***
    (0.07873) (0.10496) (0.09922) (0.08196)
    Num.Obs. 1743 1743 1743 1743
    * p < 0.1, ** p < 0.05, *** p < 0.01
  5. [3 puntos] Se desea incorporar en el análisis el género (gender) y una variable dicotómica para los hispanos (hispanic). Indique qué modelo usaría y estime dicho modelo.

    No es posible estimar los coeficientes sobre variables que no varían en el tiempo usando efectos fijos, por lo que este modelo queda descartado. Podríamos usar MCO pooled, que impone supuestos muy fuertes. La otra alternativa es un modelo de efectos aleatorios, que asume que la heterogeneidad no observada y el error no están correlacionados.

    m.sex <- plm( anti ~ self + pov + gender,
                  data=data.comp,
                  model="random",
                  index = c("id", "year"))
    
    
    modelsummary(models = list("MCO"=m.mco,
                               "MCO, errores agrupados"=m.mco,
                               "Efectos fijos"=m.fe,
                               "Efectos aleatorios"=m.re,
                               "Efectos aleatorios con género"=m.sex),
              output = 'gt',
              stars = c('***' = .01, '**' = .05, '*' = 0.1),
              vcov=list(NULL,
                        clubSandwich::vcovCR(m.mco, type='CR1',cluster=data.comp$id),
                        clubSandwich::vcovCR(m.fe, type='CR1', cluster=data.comp$id),
                        clubSandwich::vcovCR(m.re, type='CR1', cluster=data.comp$id),
                        clubSandwich::vcovCR(m.sex, type='CR1', cluster=data.comp$id)),
              fmt = 5,
              coef_map = c('self', 'pov', 'gender'),
              gof_map = "nobs")
    MCO MCO, errores agrupados Efectos fijos Efectos aleatorios Efectos aleatorios con género
    self -0.06510*** -0.06510*** -0.05150*** -0.05673*** -0.05856***
    (0.01108) (0.01369) (0.01131) (0.01022) (0.01022)
    pov 0.51581*** 0.51581*** 0.10490 0.29241*** 0.30500***
    (0.07873) (0.10496) (0.09922) (0.08196) (0.08146)
    gender -0.48047***
    (0.10710)
    Num.Obs. 1743 1743 1743 1743 1743
    * p < 0.1, ** p < 0.05, *** p < 0.01
  6. [2 puntos] Regrese al modelo que incluye solo la autoestima y el estado de pobreza como covariables. Realice una prueba de Hausman para determinar si se prefiere un modelo de efectos fijos o uno de efectos aleatorios.

    La implementación de la prueba de Hausman indica que se rechaza la H0 de que los coeficientes estimados son iguales (y que el modelo de efectos aleatorios es el adecuado). Hay evidencia de que se prefiere un modelo de efectos fijos, aunque tendremos que vivir con el hecho de no poder estimar el coeficiente asociado a las variables que no varían en el tiempo en este caso.

    phtest(m.fe, m.re)
    
     Hausman Test
    
    data:  anti ~ self + pov
    chisq = 13.578, df = 2, p-value = 0.001126
    alternative hypothesis: one model is inconsistent

Pregunta 3

Retome los datos de la pregunta 2 y el modelo del comportamiento antisocial en función de la autoestima y la pobreza. En esta pregunta mostraremos la equivalencia del estimador within con otros estimadores.

  1. [3 puntos] Compruebe que el estimador de efectos fijos es equivalente a MCO con dummies de individuos.

    Comprobamos:

    m.fe <- plm( anti ~ self + pov,
                 data=data.comp,
                 model="within",
                 index = c("id", "year"))
    
    m.dummy <- lm(anti ~ self + pov + factor(id),
                  data=data.comp)
    
    
    modelsummary(models = list("MCO"=m.fe,
                               "MCO con dummies"=m.dummy),
              output = 'gt',
              stars = c('***' = .01, '**' = .05, '*' = 0.1),
              vcov=list(clubSandwich::vcovCR(m.fe, type='CR1',cluster=data.comp$id),
                        clubSandwich::vcovCR(m.dummy, type='CR1',cluster=data.comp$id)),
              fmt = 5,
              coef_map = c('self', 'pov'),
              gof_map = "nobs")
    MCO MCO con dummies
    self -0.05150*** -0.05150***
    (0.01131) (0.01131)
    pov 0.10490 0.10490
    (0.09922) (0.09922)
    Num.Obs. 1743 1743
    * p < 0.1, ** p < 0.05, *** p < 0.01
  2. [2 puntos] Compruebe que en un modelo de efectos fijos las características que no varían en el tiempo no pueden ser identificadas. Añada la variable black para comprobarlo.

    Comprobamos que la variable simplemente es omitida del análisis:

    summary(plm( anti ~ self + pov + black,
                 data=data.comp,
                 model="within",
                 index = c("id", "year")))
    Oneway (individual) effect Within Model
    
    Call:
    plm(formula = anti ~ self + pov + black, data = data.comp, model = "within", 
        index = c("id", "year"))
    
    Balanced Panel: n = 581, T = 3, N = 1743
    
    Residuals:
          Min.    1st Qu.     Median    3rd Qu.       Max. 
    -3.7868224 -0.4706542 -0.0012721  0.4534891  3.2646729 
    
    Coefficients:
          Estimate Std. Error t-value  Pr(>|t|)    
    self -0.051495   0.010530 -4.8902 1.149e-06 ***
    pov   0.104899   0.093880  1.1174    0.2641    
    ---
    Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
    
    Total Sum of Squares:    1190.7
    Residual Sum of Squares: 1165.4
    R-Squared:      0.021182
    Adj. R-Squared: -0.46991
    F-statistic: 12.5514 on 2 and 1160 DF, p-value: 4.0471e-06
  3. [5 puntos] Compruebe que el estimador de efectos fijos es equivalente a MCO sobre el modelo en diferencias con respecto a la media. Para esto, conserve dos periodos consecutivos de datos y solo observaciones que tengan datos para las variables dependientes e independientes en los dos años que elija. Luego estime por MCO el modelo con variables transformadas.

    Nos quedamos con un subconjunto de datos:

    data.comp.sub <- data.comp %>% 
      dplyr::select(id, year, anti, self, pov) %>% 
      filter(year==90 | year==92)
    
    #Nos quedamos con los que no son NA
    data.comp.sub <- data.comp.sub[complete.cases(data.comp.sub), ]

    Creamos las variables como diferencias respecto a la media y estimamos el modelo within y el modelo de MCO en las variables transformadas:

    data.comp.sub <- data.comp.sub %>%
      group_by(id) %>% 
      mutate(m.anti = mean(anti),
             m.self = mean(self),
             m.pov = mean(pov)) %>% 
      mutate(dm.anti = anti - m.anti,
             dm.self = self - m.self,
             dm.pov = pov - m.pov)
    
    m.fe.sub <- plm( anti ~ self + pov,
                     data=data.comp.sub,
                     model="within",
                     index = c("id", "year"))
    
    
    m.demean <- lm(dm.anti ~ dm.self + dm.pov,
                   data.comp.sub)
    
    
    modelsummary(models = list("Efectos fijos"=m.fe.sub,
                               "MCO demean"=m.demean),
              output = 'gt',
              stars = c('***' = .01, '**' = .05, '*' = 0.1),
              vcov=list(clubSandwich::vcovCR(m.fe.sub, type='CR1',cluster=data.comp.sub$id),
                        clubSandwich::vcovCR(m.demean, type='CR1',cluster=data.comp.sub$id)),
              fmt = 5,
              coef_map = c('self', 'pov', 'dm.self', 'dm.pov'),
              gof_map = "nobs")
    Efectos fijos MCO demean
    self -0.03835**
    (0.01487)
    pov 0.19469
    (0.14506)
    dm.self -0.03835**
    (0.01487)
    dm.pov 0.19469
    (0.14506)
    Num.Obs. 1162 1162
    * p < 0.1, ** p < 0.05, *** p < 0.01
  4. [5 puntos] Compruebe que el estimador de efectos fijos es equivalente a MCO sobre el modelo en primeras diferencias. Parta de la muestra con dos años de la parte d. para estimar por MCO el modelo con variables transformadas.

    Usando el mismo subconjunto, calculamos ahora las primeras diferencias y estimamos:

    data.comp.sub <- data.comp.sub %>%
      group_by(id) %>% 
      mutate(d.anti = anti-dplyr::lag(anti, order_by = year),
             d.self = self-dplyr::lag(self, order_by = year),
             d.pov = pov-dplyr::lag(pov, order_by = year)) %>% 
      ungroup()
    
    
    m.difs <- lm(d.anti ~ -1 + d.self + d.pov,
                 data=data.comp.sub)
    
    
    
    modelsummary(models = list("Efectos fijos"=m.fe.sub,
                               "MCO demean"=m.demean,
                               "MCO primeras diferencias"=m.difs),
              output = 'gt',
              stars = c('***' = .01, '**' = .05, '*' = 0.1),
              vcov=list(clubSandwich::vcovCR(m.fe.sub, type='CR1',cluster=data.comp.sub$id),
                        clubSandwich::vcovCR(m.demean, type='CR1',cluster=data.comp.sub$id),
                        clubSandwich::vcovCR(m.difs, type='CR1',cluster=data.comp.sub$id)),
              fmt = 5,
              coef_map = c('self', 'pov', 'dm.self', 'dm.pov', 'd.self', 'd.pov'),
              gof_map = "nobs")
    Efectos fijos MCO demean MCO primeras diferencias
    self -0.03835**
    (0.01487)
    pov 0.19469
    (0.14506)
    dm.self -0.03835**
    (0.01487)
    dm.pov 0.19469
    (0.14506)
    d.self -0.03835**
    (0.01487)
    d.pov 0.19469
    (0.14506)
    Num.Obs. 1162 1162 581
    * p < 0.1, ** p < 0.05, *** p < 0.01

Pregunta 4

Considere los datos mlbook1.csv con información sobre 2287 estudiantes en 131 escuelas. Nos interesa la relación entre una medida de aptitud verbal, (iq_vert) y el resultado de un examen de inglés (langpost). Las variables schoolnr y pupilnr identifican a las escuelas y los estudiantes, respectivamente. El modelo a estimar es el siguiente:

\[langpost_{i}=\alpha+\beta iqvert_{i}+BX_{i}+\varepsilon_{i}\]

donde \(i\) indexa y \(X_i\) son tres características usadas como control: el sexo, sex, si el estudiante es de una población minoritaria, minority y el número de años repetidos, repeatgr.

  1. [3 puntos] ¿Por qué es posible que estemos frente a una situación de errores agrupados?

    Los datos están agrupados a nivel escuela. Los estudiantes en una misma escuela comparten características observadas y no observadas que hacen altamente probable que los factores no observables estén correlacionados entre los individuos, rompiendo el supuesto de independencia.

  2. [2 puntos] Estime la ecuación de calificación usando MCO ignorando la agrupación de datos. ¿Qué concluye respecto a la relación entre la aptitud verbal y la prueba de inglés?

    Se concluye que una hora más en la prueba de aptitud incrementa en 2.49 puntos la calificación del examen. El error estándar es 0.072.

    data.examen<-read_csv("../files/mlbook1.csv",
                          locale = locale(encoding = "latin1")) 
    
    m.mco <- lm(langpost ~ iq_verb + sex + minority + repeatgr,
                data=data.examen)
    
    
    
    modelsummary(models = list("MCO"=m.mco),
              output = 'gt',
              stars = c('***' = .01, '**' = .05, '*' = 0.1),
              vcov=list(NULL),
              fmt = 5,
              coef_map = c('iq_verb', 'sex', 'minority', 'repeatgr'),
              gof_map = "nobs")
    MCO
    iq_verb 2.48635***
    (0.07233)
    sex 2.42228***
    (0.28871)
    minority -0.03701
    (0.62762)
    repeatgr -4.40860***
    (0.43222)
    Num.Obs. 2287
    * p < 0.1, ** p < 0.05, *** p < 0.01
  3. [3 puntos] Estime ahora los errores robustos a heteroscedasticidad del tipo HC1. ¿Qué cambia y por qué en la interpretación de la relación entre la prueba de aptitud y el examen?

    El coeficiente estimado es el mismo. La fórmula empleada para calcular la varianza es una en forma de sándwich, que toma en cuenta la posible heterocedasticidad. El error estándar es apromximadamente 5% más grande, 0.076.

    modelsummary(models = list("MCO"=m.mco,
                               "MCO, errores robustos"=m.mco),
              output = 'gt',
              stars = c('***' = .01, '**' = .05, '*' = 0.1),
              vcov=list(NULL,
                        "HC1"),
              fmt = 5,
              coef_map = c('iq_verb', 'sex', 'minority', 'repeatgr'),
              gof_map = "nobs")
    MCO MCO, errores robustos
    iq_verb 2.48635*** 2.48635***
    (0.07233) (0.07587)
    sex 2.42228*** 2.42228***
    (0.28871) (0.28852)
    minority -0.03701 -0.03701
    (0.62762) (0.61245)
    repeatgr -4.40860*** -4.40860***
    (0.43222) (0.44862)
    Num.Obs. 2287 2287
    * p < 0.1, ** p < 0.05, *** p < 0.01
  4. [2 puntos] Estime la ecuación de calificación usando MCO y efectos fijos de escuela. ¿Qué resuelve este procedimiento?

    Al incluir efectos fijos a nivel escuela controlamos por características no observadas a nivel escuela. Estas diferencias se incorporan en el modelo como desplazamientos de la ordenada al origen. Este procedimiento no tiene nada que ver con la agrupación de errores.

    m.mco.ef <- lm(langpost ~ iq_verb + sex + minority + repeatgr + factor(schoolnr),
                   data=data.examen)
    
    modelsummary(models = list("MCO"=m.mco,
                               "MCO, errores robustos"=m.mco,
                               "Efectos fijos"=m.mco.ef),
              output = 'gt',
              stars = c('***' = .01, '**' = .05, '*' = 0.1),
              vcov=list(NULL,
                        "HC1",
                        "HC1"),
              fmt = 5,
              coef_map = c('iq_verb', 'sex', 'minority', 'repeatgr'),
              gof_map = "nobs")
    MCO MCO, errores robustos Efectos fijos
    iq_verb 2.48635*** 2.48635*** 2.25997***
    (0.07233) (0.07587) (0.07310)
    sex 2.42228*** 2.42228*** 2.41900***
    (0.28871) (0.28852) (0.26990)
    minority -0.03701 -0.03701 0.23262
    (0.62762) (0.61245) (0.70803)
    repeatgr -4.40860*** -4.40860*** -4.43503***
    (0.43222) (0.44862) (0.42731)
    Num.Obs. 2287 2287 2287
    * p < 0.1, ** p < 0.05, *** p < 0.01
  5. [5 puntos] Estime la ecuación de calificación usando MCO y con errores agrupados a nivel escuela (sin efectos fijos de escuela). ¿Qué resuelve este procedimiento?

    Al estimar los errores agrupados y robustos a heterocedasticidad se toma en cuenta la correlación que existe en los errores dentro de cada escuela. Los errores agrupados estimados con la opción cluster asumen correlación de errores dentro del grupo, pero no entre grupos. Con respecto a las partes b. y c., el error estándar asociado al tiempo dedicado a la tarea es aproximadamente 20% mayor. Este es un ejemplo típico en el que los errores agrupados se inflan con respecto a los errores de MCO clásicos y los errores robustos.

    Nota: es posible que los errores agrupados sean menores que los errores de MCO. Para ver eso, considere un modelo simple con datos agrupados de la forma siguiente: \[y_{ig=\alpha+\beta x_{ig}+u_{ig}}\] donde \(x_{ig}\) es un regresor escalar.

    Se asume que el tamaño promedio de los grupos es \(\bar{N}_g\). Moulton (1990) muestra que el error estándar de MCO esta sesgado hacia abajo por una cantidad igual a la raíz de \(\tau \approx 1 +\rho_x \rho_u (\bar{N}_g-1)\), donde \(\rho_x\) es la correlación dentro de los grupos de \(x\) y \(\rho_u\) es la correlación dentro de los grupos de los errores. Esto implica que para obtener el error correcto que toma en cuenta la agrupación hay que multiplicar el error de MCO por la raíz de \(\tau\). Sin embargo, note que dependiendo del signo y la magnitud de \(\rho_x\) y \(\rho_u\), la raíz de \(\tau\) puede llegar a ser menor que 1 y, por tanto, el error agrupado puede llegar a ser menor que el de MCO. \(\tau\) se conoce como el factor de Moulton y puede ser extendido para un modelo más complicado. La intuición funciona de manera similar para un modelo más complicado: todo depende de las correlaciones entre grupos de los regresores y la correlación de los errores.

    modelsummary(models = list("MCO"=m.mco,
                               "MCO, errores robustos"=m.mco,
                               "Efectos fijos"=m.mco.ef,
                               "MCO, errores agrupados"=m.mco),
              output = 'gt',
              stars = c('***' = .01, '**' = .05, '*' = 0.1),
              vcov=list(NULL,
                        "HC1",
                        "HC1",
                        clubSandwich::vcovCR(m.mco, type='CR1',cluster=data.examen$schoolnr)),
              fmt = 5,
              coef_map = c('iq_verb', 'sex', 'minority', 'repeatgr'),
              gof_map = "nobs")
    MCO MCO, errores robustos Efectos fijos MCO, errores agrupados
    iq_verb 2.48635*** 2.48635*** 2.25997*** 2.48635***
    (0.07233) (0.07587) (0.07310) (0.08986)
    sex 2.42228*** 2.42228*** 2.41900*** 2.42228***
    (0.28871) (0.28852) (0.26990) (0.28436)
    minority -0.03701 -0.03701 0.23262 -0.03701
    (0.62762) (0.61245) (0.70803) (0.85396)
    repeatgr -4.40860*** -4.40860*** -4.43503*** -4.40860***
    (0.43222) (0.44862) (0.42731) (0.39547)
    Num.Obs. 2287 2287 2287 2287
    * p < 0.1, ** p < 0.05, *** p < 0.01
  6. [5 puntos] Estime la ecuación de calificación usando MCO, efectos fijos de escuela y con errores agrupados a nivel escuela. ¿Qué resuelve este procedimiento?

    Al controlar por características no observadas de las escuelas empleando efectos fijos por escuela y además estimando los errores que toman en cuenta la estructura agrupada de los errores obtenemos un coeficiente estimado de 2.26, pero con un error estándar mayor, 0.08879.

    modelsummary(models = list("MCO"=m.mco,
                               "MCO, errores robustos"=m.mco,
                               "Efectos fijos"=m.mco.ef,
                               "MCO, errores agrupados"=m.mco,
                               "Efectos fijos, errores agrupados"=m.mco.ef),
              output = 'gt',
              stars = c('***' = .01, '**' = .05, '*' = 0.1),
              vcov=list(NULL,
                        "HC1",
                        "HC1",
                        clubSandwich::vcovCR(m.mco, type='CR1',cluster=data.examen$schoolnr),
                        clubSandwich::vcovCR(m.mco.ef, type='CR1',cluster=data.examen$schoolnr)),
              fmt = 5,
              coef_map = c('iq_verb', 'sex', 'minority', 'repeatgr'),
              gof_map = "nobs")
    MCO MCO, errores robustos Efectos fijos MCO, errores agrupados Efectos fijos, errores agrupados
    iq_verb 2.48635*** 2.48635*** 2.25997*** 2.48635*** 2.25997***
    (0.07233) (0.07587) (0.07310) (0.08986) (0.08879)
    sex 2.42228*** 2.42228*** 2.41900*** 2.42228*** 2.41900***
    (0.28871) (0.28852) (0.26990) (0.28436) (0.27519)
    minority -0.03701 -0.03701 0.23262 -0.03701 0.23262
    (0.62762) (0.61245) (0.70803) (0.85396) (0.69585)
    repeatgr -4.40860*** -4.40860*** -4.43503*** -4.40860*** -4.43503***
    (0.43222) (0.44862) (0.42731) (0.39547) (0.39617)
    Num.Obs. 2287 2287 2287 2287 2287
    * p < 0.1, ** p < 0.05, *** p < 0.01