¿Cómo dejar de preocuparse y empezar a utilizar eficazmente los paquetes de R para la econometría?

Actualización: GrantMcDermott ha proporcionado una solución más fácil para los errores estándar agrupados usando R base.

Dedicado a todos los que escriben comentarios estúpidos sobre las personas con discapacidades, LGBT+, no blancos y mujeres, y especialmente a los que promueven la falsa inclusión. Como siempre, este artículo es una opinión personal y no involucra a las personas que trabajan conmigo a menos que estén explícitamente de acuerdo con estas palabras.

Contexto

Cuando se hace econometría, usar R no es inmediatamente sencillo y tiene una curva de aprendizaje empinada que podemos aplanar con algunas adiciones a la configuración de R.

Si usted recién está probando R y viene de Stata es posible que ya se haya sorprendido con al menos uno de estos (me pasó a mí en 2015):

  • No **es fácil calcular errores estándar clusterizados para modelos lineales
  • La obtención de un Pseudo-R2 para modelos lineales generalizados requiere trabajo adicional
  • Hay diferentes paquetes para usos similares
  • No todos los paquetes son compatibles entre sí

Probablemente piense que R es difícil de usar, y que su transparencia viene con esa dificultad adicional. Pero si utiliza paquetes creados por la comunidad, como el magnífico tidyverse, entonces tendrá acceso a funciones que no vienen con R cuando lo instala, pero que facilitan mucho su trabajo en R (¡piense en esttout en Stata!).

Pero, ¿qué pasa con otros paquetes? a veces es difícil encontrar tutoriales paso a paso del tipo ‘En Stata se hace X, en R se hace Y’. Aquí voy a explicar un poco de econometría usando R y cómo crear una colección de paquetes.

Piense en los paquetes de R como en las extensiones del navegador, estos sirven para un propósito particular para un grupo particular de usuarios. Por ejemplo, no todos los usuarios de Firefox/Chrome necesitan la extensión Enhancer for YouTube, ya que no todos los usuarios ven YouTube. La instalación por defecto de R está pensada para ser ligera y de uso general.

Puede que haya encontrado Stackoverflow y el hashtag #rstats en Twitter donde hay gente que ayuda mucho y hay buenas comunidades como rOpenSci. Pero, siendo R un grupo muy grande (y creciente) de personas, no olvide que es una muestra de la raza humana, y por lo tanto encontrará tutoriales escritos en un lenguaje que no es diferente de lo que la intersección de muchos cultos y pseudociencia (no tan mala como la de los criptobros por supuesto) escribiría, y habrá gente que responda “Búscalo en Google” cuando tenga un problema que no pueda resolver. No te decepciones, por favor. Pregunta a las personas adecuadas.

Solo manda un poco de sombra a los que no contribuyen a una comunidad sana. sólo un GIF de RuPaul en referencia al concepto de 'sombra'

Afortunadamente, si perteneces a un grupo sub-representado o estás empezando con R, hay una fracción importante de la comunidad R que es realmente inclusiva y que no tiene una baja autoestima hasta el punto de sentirse asustada por personas que son diferentes a ellos. Por favor, cuenten con los valiosos y si ves algo, di algo.

Instalación de paquetes R desde CRAN (e.g. la fuente oficial)

CRAN, el servidor oficial de paquetes de R, hace un muy buen trabajo al proporcionar algunos paquetes como

  • dplyr: Una herramienta rápida y consistente para trabajar con objetos tipo data frame, tanto en memoria y fuera de memoria.
  • ggplot2: Un sistema para crear gráficos “declarativamente”, basado en “La gramática de gráficos”.
  • msm: Modelos de Markov multiestado y de Markov latente en tiempo continuo

Pero hay más paquetes fuera de CRAN, o que fueron eliminados de él debido a cambios de política, etc. Uno de ellos es RVAideMemoire, que proporciona ‘Procedimientos de prueba y trazado para la bioestadística’ y que es muy útil para trabajar con modelos lineales generalizados.

Dplyr y ggplot2 están ampliamente cubiertos en muchos otros blogs, ambos son parte del Tidyverse, que es una especie de mochila que contiene diferentes paquetes.

Puedes instalar los paquetes de R con install.packages("tidyverse") o similar. Cuando se necesita instalar desde otras fuentes, la cosa se complica.

Recomiendo instalar paquetes desde CRAN, a menos que necesites una nueva funcionalidad, porque a veces las fuentes no oficiales introducen cambios experimentales.

Instalación de paquetes desde fuentes no oficiales (e.g. GitHub)

Como forma de facilitar la instalación de paquetes desde otras fuentes que no sean CRAN, tenemos el paquete ‘remotes’. Por ejemplo, el paquete gravity no ha sido actualizado en CRAN debido a algunos cambios de política que me impiden añadir algunas mejoras sin tomar desvíos. Estos desvíos son irrelevantes para el usuario final, que puede instalar la versión mejorada de gravity desde el repositorio del paquete.

La forma más fácil sería:

install.packages("remotes")
remotes::install_github("pachadotdev/gravity")

Para los paquetes retirados de CRAN es lo mismo:

install.packages("remotes")
remotes::install_github("pachadotdev/RVAideMemoire")

Por supuesto, puede instalar desde diferentes repositorios como lrberge/fixest en lugar de install.packages("fixest"). Por cierto, fixest es más que excelente y es muy útil para ajustar modelos con grandes efectos.

Ajuste de modelos con el paquete fixest versus R base

Al principio del post mencioné cosas que no son fáciles de hacer con R estándar.

Como ejemplo, digamos que necesitamos replicar el primer modelo OLS de la página 42 de An Advanced Guide to Trade Policy Analysis (véa el manual de soluciones).

La especificación del modelo es: \[\begin{align*} \log X_{ij,t} =& \:\beta_0 + \beta_1 DIST_{i,j} + \beta_2 CNTG_{i,j} + \beta_3 LANG_{i,j} + \beta_4 CLNY_{i,j} + \beta_5 \log Y_{i,t} +\\ \text{ }& \:\beta_6 \log E_{j,t} + \varepsilon_{ij,t} \end{align*}\]

donde todas las variables, excepto \(Y,E\) están en el dataset incluido en el paquete

> tradepolicy::agtpa_applications
# A tibble: 99,981 x 17
   exporter importer pair_id  year    trade   dist  cntg  lang  clny   rta rta_lag3 rta_lag4 rta_lag6 rta_lag8
   <chr>    <chr>      <dbl> <dbl>    <dbl>  <dbl> <dbl> <dbl> <dbl> <dbl>    <dbl>    <dbl>    <dbl>    <dbl>
 1 ARG      ARG        12339  1986 61289.     534.     0     0     0     0        0        0        0        0
 2 ARG      AUS            1  1986    27.8  12045.     0     0     0     0        0        0        0        0
 3 ARG      AUT            2  1986     3.56 11751.     0     0     0     0        0        0        0        0
 4 ARG      BEL            4  1986    96.1  11305.     0     0     0     0        0        0        0        0
 5 ARG      BGR            3  1986     3.13 12116.     0     0     0     0        0        0        0        0
 6 ARG      BOL            6  1986    52.7   1866.     1     1     0     0        0        0        0        0
 7 ARG      BRA            8  1986   405.    2392.     1     0     0     0        0        0        0        1
 8 ARG      CAN           10  1986    48.3   9391.     0     0     0     0        0        0        0        0
 9 ARG      CHE           12  1986    23.6  11233.     0     0     0     0        0        0        0        0
10 ARG      CHL           14  1986   109.    1157.     1     1     0     0        0        0        0        1
# … with 99,971 more rows, and 3 more variables: rta_lag9 <dbl>, rta_lag12 <dbl>, rta_lead4 <dbl>

Lo que se haría para replicar los resultados es seguir los pasos del libro pero en R:

  1. Cargar el paquete tradepolicy para tener los datos del libro
  2. Filtrar las observaciones para un rango de años (1986, 1990, 1994, 1998, 2002 y 2006)
  3. Transformar algunas variables a escala logarítmica (comercio y dist) y crear nuevas variables a partir de las del conjunto de datos original
  4. Eliminar los casos en los que el exportador y el importador son el mismo
  5. Eliminar las observaciones en las que el flujo comercial es cero

El primer punto consistiría en

install.packages("tradepolicy")
library(tradepolicy)

Tradepolicy también llamará a dplyr y a otros cuando se cargue, y dplyr proporciona el operador %>% para encadenar operaciones. Para completar el punto 2 con dplyr (mucho más fácil que usando R base) puede hacer esto:

ch1_application1 <- agtpa_applications %>% 
  filter(year %in% seq(1986, 2006, 4))

Luego el paso 3 es similar (consulte el libro para la definición de \(Y,E\)):

ch1_application1 <- ch1_application1 %>% 
  mutate(
    log_trade = log(trade),
    log_dist = log(dist)
  )
  
ch1_application1 <- ch1_application1 %>% 
  # Create Yit
  group_by(exporter, year) %>% 
  mutate(
    y = sum(trade),
    log_y = log(y)
  ) %>% 
  
  # Create Eit
  group_by(importer, year) %>% 
  mutate(
    e = sum(trade),
    log_e = log(e)
  )

Los pasos 4 y 5 son más de lo mismo:

ch1_application1 <- ch1_application1 %>% 
  filter(exporter != importer, trade > 0)

Ahora estoy listo para ajustar el modelo, pero en lugar de los pasos del manual de soluciones usaré fixest para hacerlo, sin usar las funciones de tradepolicy:

model1 <- feols(
  log_trade ~ log_dist + cntg + lang + clny + log_y + log_e,
  data = ch1_application1,
  cluster = ~ pair_id
)

Como puede observar en el código, he especificado una variable de clusterización para los errores estándar. R tiene una función por defecto para los modelos lineales, lm, pero requiere el uso de paquetes adicionales para la clsuterización.

Vamos a explorar la salida del modelo:

> model1
OLS estimation, Dep. Var.: log_trade
Observations: 25,689 
Standard-errors: Clustered (pair_id) 
              Estimate Std. Error  t value Pr(>|t|))    
(Intercept) -11.283000   0.295827 -38.1410 < 2.2e-16 ***
log_dist     -1.001600   0.027340 -36.6350 < 2.2e-16 ***
cntg          0.573805   0.184710   3.1065  0.001916 ** 
lang          0.801548   0.082102   9.7629 < 2.2e-16 ***
clny          0.734853   0.144193   5.0963  3.74e-07 ***
log_y         1.190200   0.009456 125.8700 < 2.2e-16 ***
log_e         0.907588   0.009910  91.5850 < 2.2e-16 ***
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
RMSE: 1.7427   Adj. R2: 0.758469

¡Este es exactamente el mismo resultado del libro!

Ahora hagámoslo con lm y algunos refuerzos, sándwich para los errores estándar y lmtest para la prueba de coeficientes:

install.packages(c("sandwich","lmtest"))
library(sandwich)
library(lmtest)

model2 <- lm(
  log_trade ~ log_dist + cntg + lang + clny + log_y + log_e,
  data = ch1_application1
)

coef_test <- coeftest(
  model2,
  vcov = vcovCL, 
  cluster = ~ pair_id
)

Ahora veamos coef_test:

> coef_test

t test of coefficients:

               Estimate  Std. Error  t value  Pr(>|t|)    
(Intercept) -11.2830801   0.2958274 -38.1408 < 2.2e-16 ***
log_dist     -1.0016075   0.0273400 -36.6353 < 2.2e-16 ***
cntg          0.5738051   0.1847098   3.1065  0.001895 ** 
lang          0.8015482   0.0821017   9.7629 < 2.2e-16 ***
clny          0.7348535   0.1441929   5.0963 3.488e-07 ***
log_y         1.1902357   0.0094560 125.8716 < 2.2e-16 ***
log_e         0.9075884   0.0099098  91.5846 < 2.2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Por supuesto, se puede continuar el ejemplo y construir sobre el segundo modelo hasta obtener un resumen completo como este pero informando de los errores estándar clusterizados:

> summary(model2)

Call:
lm(formula = log_trade ~ log_dist + cntg + lang + clny + log_y + 
    log_e, data = ch1_application1 %>% filter(trade > 0))

Residuals:
     Min       1Q   Median       3Q      Max 
-14.5421  -0.8281   0.1578   1.0476   7.6585 

Coefficients:
              Estimate Std. Error t value Pr(>|t|)    
(Intercept) -11.283080   0.151732  -74.36  < 2e-16 ***
log_dist     -1.001607   0.014159  -70.74  < 2e-16 ***
cntg          0.573805   0.074427    7.71 1.31e-14 ***
lang          0.801548   0.033748   23.75  < 2e-16 ***
clny          0.734853   0.070387   10.44  < 2e-16 ***
log_y         1.190236   0.005402  220.32  < 2e-16 ***
log_e         0.907588   0.005577  162.73  < 2e-16 ***
---
Signif. codes:  0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1

Residual standard error: 1.743 on 25682 degrees of freedom
Multiple R-squared:  0.7585,    Adjusted R-squared:  0.7585 
F-statistic: 1.345e+04 on 6 and 25682 DF,  p-value: < 2.2e-16

Como puede ver, fixest ayuda mucho, que es lo que hacen los paquetes.

Usando el R-Universe

El R-Universe, creado por Jeroen Ooms, proporciona una forma muy sencilla de crear repositorios personales tipo CRAN, lo que significa una forma de mostrar su colección de herramientas en uso a la comunidad.

Además, puede usarlo para publicar artículos usando rmarkdown, un paquete de R que le permite escribir texto y código, y generar un documento PDF incluso con las plantillas requeridas para la revista en la que se desea publicar.

Para instalar leontief o cualquier otro paquete de mi sección de R-Universe puede hacer esto:

# Install new leontief version
install.packages("leontief", repos = "https://pachadotdev.r-universe.dev")

Aquí tiene algunos universos que puede añadir permanentemente a su configuración:

# Enable some universes
options(repos = c(
    pachadotdev = 'https://pachadotdev.r-universe.dev',
    tidyverse = 'https://github.com/r-universe/tidyverse',
    rlib = 'https://github.com/r-universe/r-lib',
    tidymodels = 'https://github.com/r-universe/tidymodels',
    rspatial = 'https://github.com/r-universe/r-spatial'
    CRAN = 'https://cloud.r-project.org'))

Cree su propio R-Universe

Esto es particularmente útil si enseña cursos y proporciona, por ejemplo paquetes de datos para sus estudiantes, o si tiene paquetes que no le preocupa enviar a CRAN.

Para unirse al R-Universe, se necesita un perfil GitHub y una muy buena referencia para empezar con Git y GitHub es Happy Git and GitHub for the useR.

Por ejemplo, tradepolicy es un repositorio donde tengo todos los códigos para reproducir los resultados de An Advanced Guide to Trade Policy Analysis]. Es de mi interés listar los paquetes de R utilizados allí, y en otros repositorios, para que en el R-Universe otros usuarios puedan descubrir fácilmente las herramientas que utilizo.

Cuando visite el R-Universe verá esta página de inicio donde tiene que hacer click en ‘Configurar’.

página de inicio de r-universe, en el menú de la izquierda tiene que
hacer click en 'configuración'

A continuación, tiene que seleccionar su perfil personal (o una cuenta de organización si tiene autorización).

página de github en la que hay que pulsar la cuenta a utilizar con el
r-universe

Ahora puedes elegir todos tus repositorios o sólo algunos. Yo elegiré sólo gravity y agtpa, y algunos otros y luego hice click en ‘Instalar’.

si quiere añadir todos sus repositorios seleccione 'todos los repositorios', luego
haga click en instalar

para seleccionar unos pocos repositorios haga click en 'only select repositories' y añada lo que
que necesite de la lista, luego haga click en instalar

Se le pedirá que confirme el acceso de R-Universe a sus repositorios.

escribe tu contraseña y haz click en el botón 'confirmar contraseña'

Una vez que esté listo, verá esto

el sitio que muestra r-universe cuando confirmas el acceso desde github

Y la URL de su repositorio será de la forma

https://githubusername.r-universe.dev/
(https://pachadotdev.r-universe.dev/ en mi caso)

Unos minutos después de configurarlo, el sitio tendrá este aspecto:

el sitio que muestra r-universe cuando está listo o casi listo

Tenga en cuenta que los repositorios que he añadido son paquetes de R, puede simplemente añadir una lista de paquetes (por ejemplo, ‘uso dplyr, haven, etc’) a R-Universe si crea un repositorio con un archivo packages.json.

Así es como se añade una lista de tus paquetes favoritos al R-Universe:

  1. Vaya a github.com y cree un nuevo repositorio

haz click en el signo más de la parte superior derecha y luego en 'nuevo repositorio'

  1. Nombre su repositorio como “universe”, ponganlo como público y con una licencia adecuada (a mí me gusta la licencia Apache) y luego proceda con los pasos de la imagen. Los pasos específicos son

escribir un nombre de repositorio, seleccionar la opción 'pública', elegir la licencia 
(por ejemplo, apache), luego haga click en 'crear repositorio'

  1. Crea un nuevo archivo packages.json que contenga tus paquetes

haga click en 'añadir archivo' y luego en 'crear nuevo archivo'

escriba el nombre del archivo 'packages.json' en el primer cuadro, luego liste sus 
paquetes en el cuadro de texto más grande

haga un desplazamiento hacia abajo y haga click en 'commit changes'

Este es el texto tipo que he utilizado:

[
  {
    "paquete": "tradepolicy",
    "url": "https://github.com/pachadotdev/tradepolicy"
  },
  {
    "paquete": "leontief",
    "url": "https://github.com/pachadotdev/leontief"
  },
  {
    "paquete": "tradestatistics",
    "url": "https://github.com/ropensci/tradestatistics"
  },
  {
    "paquete": "gravedad",
    "url": "https://github.com/pachadotdev/gravity"
  },
  {
    "paquete": "flecha",
    "url": "https://github.com/apache/arrow"
  },
  {
    "paquete": "RVAideMemoire",
    "url": "https://github.com/pachadotdev/RVAideMemoire"
  }
]
  1. Si ha añadido algunos repositorios a R-Universe, añada el repositorio recién creado de universe, de lo contrario ya está listo

haga click en el icono del perfil en la parte superior derecha y luego en 'settings'

haga click en 'application' y luego en 'configure' a la derecha del icono 'r-universe'

añada el nuevo repositorio, y luego haga click en 'save'