Ya hice otro tutorial de datos en panel en un post anterior. En este post llegaremos a modelos dinámicos haciendo uso de los paquetes plm y gmm.
Usaremos la base de datos de Grunfeld (1958) sobre las determinantes de la inversión corporativa. Esta base de datos consta de 20 observaciones anuales de tres variables: Inversión Bruta en términos reales (invest), valor de la empresa en términos reales (value) y valor del stock de capital en términos reales (capital). Esta base se extiende a 11 empresas grandes en Estados Unidos en el periodo de 1935 a 1954. Esta base de datos es clásica para el estudio de bases de datos en panel, proviene de la tesis doctoral en la Universidad de Chicago de Yehuda Grunfeld y la podemos encontrar en el paquete AER, del libro Applied Econometrics with R. Para cargar la base de datos, simplemente tecleamos las siguientes líneas de código. Además, para ejemplificar, haremos de la base de datos un subconjunto que solamente utilice cuatro empresas.
data("Grunfeld", package = "AER") #Para cargar los datos. Asegurate de que has instalado ya la biblioteca AER antes de correr esta linea.
library('plm')## Loading required package: Formulagr <- subset(Grunfeld, firm %in% c("General Electric", "General Motors", "US Steel", "Westinghouse"))
p.gr <- pdata.frame(gr, index = c("firm", "year"))
head(p.gr)## invest value capital firm year
## General Motors-1935 317.6 3078.5 2.8 General Motors 1935
## General Motors-1936 391.8 4661.7 52.6 General Motors 1936
## General Motors-1937 410.6 5387.1 156.9 General Motors 1937
## General Motors-1938 257.7 2792.2 209.2 General Motors 1938
## General Motors-1939 330.8 4313.2 203.4 General Motors 1939
## General Motors-1940 461.2 4643.9 207.2 General Motors 1940El modelo de regresi??n en panel se puede modelar como \[ Y_{it}=\beta_1V_{it}+\beta_2K_{it}+\alpha_i+\varepsilon_{it}, \]
donde \(Y_it\) representa la inversi??n de la empresa \(i\in[1,n]\) en el periodo \(t\in[1,T]\) y \(\alpha_i\) denota el efecto particular del modelo. Esta elecci??n de los ??ndices es una mera convenci??n que decidimos acatar en este texto.
Para realizar la estimaci??n del modelo descrito con anterioridad, dependemos de algunos supuestos respecto al intercepto, los coeficientes y los t??rminos de error. Para el tratamiento de datos en panel, se consideran las siguientes opciones:
Todos estos casos pueden darse y aumentan en complejidad y en realismo de los supuestos. Hay que recordar que al momento de realizar un modelo, siempre existe un costo entre la complejidad del modelo y su realismo: un modelo mas realista suele requerir de supuestos mas complejos, lo cual complica su an??lisis.
El supuesto inicial implica que no existe heterogeneidad entre individuos ni en el tiempo (para ver c??mo probar heterogeneidad, ver el otro tutorial de Panel) y por lo tanto, se puede realizar una regresi??n de M??nimos Cuadrados Ordinarios (MCO). Para hacer esta regresi??n podemos ejecutar la siguiente l??nea.
gr_ols <- lm(invest ~ value + capital, data = p.gr)
summary(gr_ols)##
## Call:
## lm(formula = invest ~ value + capital, data = p.gr)
##
## Residuals:
## Min 1Q Median 3Q Max
## -319.68 -99.95 1.96 65.99 336.21
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -62.83184 29.72538 -2.114 0.0378 *
## value 0.11052 0.01378 8.023 9.19e-12 ***
## capital 0.30046 0.04940 6.082 4.27e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 142.9 on 77 degrees of freedom
## Multiple R-squared: 0.7547, Adjusted R-squared: 0.7483
## F-statistic: 118.4 on 2 and 77 DF, p-value: < 2.2e-16Podemos observar que el resultado no es tan malo. De hecho, los coeficientes tienen sentido y son estad??sticamente significativos y el valor de \(R^2\) es razonable. El valor de los coeficientes es positivos, que es lo que a priori esper??bamos. El problema en estos resultados est?? en el valor del estad??stico Durbin-Watson muy baja, lo cual podr??a ser un indicador de la presencia de autocorrelaci??n o de problemas en la especificaci??n del modelo:
library(car)
durbinWatsonTest(gr_ols)## lag Autocorrelation D-W Statistic p-value
## 1 0.8315513 0.333299 0
## Alternative hypothesis: rho != 0Un posible problema de especificaci??n es que estamos suponiendo que el valor del intercepto es el mismo para las diferentes empresas. Este supuesto puede resultar irreal si suponemos algo de heterogeneidad entre ellas.
En este modelo, se considera que, a pesar de haber heterogeneidad entre individuos, no la hay as?? a trav??s del tiempo. Esto se puede lograr creando una variable dicot??mica para el efecto de cada uno de los individuos. Esto es, podemos expresar el modelo como \[ Y_{it} = \alpha_0 + \alpha_1\delta_{1i} + \alpha_{2}\delta_{2i}+\alpha_3\delta_{3i}+\beta_1V_{it}+\beta_2K_{it}+\varepsilon_{it}, \]
donde \(\delta_{1i}\) es una variable dicot??mica o dummy que denota que los datos son de General Motors cuando esta es igual a 1 y cuando se trata de cualquier otra empresa esta tiene valor de 0. De igual manera, cuando de US Steel y Westinghouse, las variables \(\delta_{2i}\) y \(\delta_{3i}\) se convierten en 1, respectivamente. Note que no hay una variable dicot??mica que sea indicador de General Electric, esto para evitar la existencia de colinealidad perfecta (v??ase la trampa de las variables dicot??micas). En otras palabras, \(\alpha_0\) representa el intercepto que tiene la regresi??n ??nicamente para General Electric, mientras que el resto de los coeficientes de \(\alpha_j\) representan el valor diferencial en las pendientes. Por supuesto, es posible elegir a cualquier empresa o individuo como el punto de comparaci??n; esta elecci??n puede ser arbitraria o puede servir con fines ilustrativos a lo que se quiere mostrar en la investigaci??n. Aqu?? un ejemplo de una regresi??n con efectos fijos en R:
gr_fe <- plm(invest ~ value + capital, data = p.gr, model = "within")
summary(gr_fe)## Oneway (individual) effect Within Model
##
## Call:
## plm(formula = invest ~ value + capital, data = p.gr, model = "within")
##
## Balanced Panel: n = 4, T = 20, N = 80
##
## Residuals:
## Min. 1st Qu. Median 3rd Qu. Max.
## -185.00 -48.30 9.33 40.50 198.00
##
## Coefficients:
## Estimate Std. Error t-value Pr(>|t|)
## value 0.108400 0.017566 6.1711 3.3e-08 ***
## capital 0.345058 0.026708 12.9195 < 2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Total Sum of Squares: 2171500
## Residual Sum of Squares: 422220
## R-Squared: 0.80556
## Adj. R-Squared: 0.79242
## F-statistic: 153.291 on 2 and 74 DF, p-value: < 2.22e-16summary(fixef(gr_fe))## Estimate Std. Error t-value Pr(>|t|)
## General Motors -85.515 73.490 -1.1636 0.248310
## US Steel 94.988 36.664 2.5907 0.011530 *
## General Electric -246.228 35.938 -6.8515 1.857e-09 ***
## Westinghouse -59.386 20.233 -2.9351 0.004439 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1Note c??mo R ya tiene correctamente identificados los valores de \(n = 4\) y \(T = 20\), es decir, el n??mero de empresas, en este caso y de periodos de estudio que abarca el panel. El estimador de efectos fijos (o estimador de efectos within o dentro de) sirve para estimar los coeficientes del modelo incluyendo los efectos fijos.
El resumen de la regresi??n nos indica que la informaci??n est?? en un panel balanceado (esto significa que cada unidad transversal –cada empresa –tiene el mismo n??mero de observaciones de series de tiempo. En este caso, 20), y nos proporciona informaci??n respecto a los coeficientes (sin el valor del efecto fijo) y medidas de bondad de ajuste. A priori, nosotros esperamos que exista un efecto positivo entre las variables de inter??s (el nivel de capital y el valor de la empresa con el nivel de inversi??n).
El m??todo fixef() nos permite examinar los efectos fijos individuales para cada individuo. Esta funci??n tiene un summary() asociado que nos permite anexar los efectos individuales.
Alternativamente a realizar un modelo de efectos fijos, se puede realizar una regresi??n agrupada, que es lo mismo que nuestro modelo de m??nimos cuadrados que hemos hecho previamente. Esta regresi??n es nuestro punto de referencia para comparar nuestro modelo con efectos fijos. Podemos hacer una prueba F entre estos dos modelos.
\[ F\frac{(R^2_{FE}-R^2_{OLS})/3}{(1-R^2_{FE})/74}= F\frac{(0.9345-0.7547)/3}{(1-0.9345)/74} = 67.71 \]
Adem??s, el paquete plm contiene su propio m??todo para comparar dos modelos en panel. A continuaci??n se presenta el m??todo para hacer una regresi??n agrupada con plm, que debe de generar el mismo resultado que hemos hecho anteriormente con una regresi??n por MCO.
gr_pool <- plm(invest ~ value + capital, data = p.gr, model = "pooling")
pFtest(gr_fe, gr_pool)##
## F test for individual effects
##
## data: invest ~ value + capital
## F = 67.215, df1 = 3, df2 = 74, p-value < 2.2e-16
## alternative hypothesis: significant effectsEsta prueba nos indica que hay suficiente variaci??n entre empresas. Esto implica que, efectivamente los interceptos de cada una de las empresas es distinto. Esto se puede deber a caracter??sticas inherentes a cada empresa que no est??n expresados en el modelo y pueden ir desde efectos del mercado en el que se encuentran hasta las distintas habilidades de direcci??n que tienen los gerentes de cada empresa.
Una de las desventajas del modelo de efectos fijos es que suele resultar costoso en t??rminos de grados de libertad, debido a la inclusi??n de variables dicot??micas. La inclusi??n de estas variables nos ayuda a compensar nuestra falta de conocimiento del modelo real que gobierna a nuestra relaci??n entre variables. En tal caso, ser??a m??s honesto expresar directamente nuestra ignorancia respecto a las fuerzas que hay detr??s del fen??meno que estudiamos y expresarla en t??rmino de variables aleatorias. A este modelo, le llamamos modelo de Efectos Fijos (EF).
Comenzamos con nuestro modelo est??ndar, con \[ Y_{it} = \beta_{0i} +\beta_{1i}X_{1it}+\beta_{2i}X_{2it}+\varepsilon_{it} \]
En este modelo, consideramos a \(\beta_{0i}\) como una variable aleatoria en el que \(\frac{1}{n}\sum_{i=1}^n\beta_{0i} = \beta_0\), en otras palabras \[ \beta_{0i} = \beta_0 + \eta_i \] donde \(\eta_i\sim\mathcal{N}(0,\sigma^2_\eta)\). Si incluimos este supuesto en el modelo original, tenemos entonces \[ Y_{it} = \beta_{0i} +\beta_{1i}X_{1it}+\beta_{2i}X_{2it}+\varepsilon_{it}+ \eta_i \] Suponemos que los t??rminos de error no est??n correlacionados entre si y no est??n correlacionados a trav??s del tiempo ni transversalmente entre empresas.
Existen diferentes t??cnicas para el an??lisis de datos en panel. R usa por default la t??cnica de Swamy-Arora, pero en este tutorial haremos uso ilustrativo de la t??cnica de Wallace-Hussain (cf. Baltagi, 2008).
gr_re <- plm(invest ~ value + capital, data = p.gr, model = "random", random.method = "walhus")
summary(gr_re)## Oneway (individual) effect Random Effect Model
## (Wallace-Hussain's transformation)
##
## Call:
## plm(formula = invest ~ value + capital, data = p.gr, model = "random",
## random.method = "walhus")
##
## Balanced Panel: n = 4, T = 20, N = 80
##
## Effects:
## var std.dev share
## idiosyncratic 5786.52 76.07 0.294
## individual 13872.63 117.78 0.706
## theta: 0.8571
##
## Residuals:
## Min. 1st Qu. Median 3rd Qu. Max.
## -194.00 -46.20 1.35 42.00 198.00
##
## Coefficients:
## Estimate Std. Error t-value Pr(>|t|)
## (Intercept) -72.632165 68.908322 -1.0540 0.2952
## value 0.107899 0.016658 6.4773 7.973e-09 ***
## capital 0.344258 0.026868 12.8128 < 2.2e-16 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Total Sum of Squares: 2258100
## Residual Sum of Squares: 446130
## R-Squared: 0.80243
## Adj. R-Squared: 0.7973
## F-statistic: 156.365 on 2 and 77 DF, p-value: < 2.22e-16Note que los resultados de ambos m??todos son similares, al menos para estos datos. Para comprobar si es realmente necesario este tipo de modelo, se realiza una prueba de multiplicadores de Lagrange. Esta es sencilla de hacer en R,
plmtest(gr_fe)##
## Lagrange Multiplier Test - (Honda) for balanced panels
##
## data: invest ~ value + capital
## normal = 19.454, p-value < 2.2e-16
## alternative hypothesis: significant effectsPor default, esta prueba usa la prueba de Honda (1985).
Los efectos aleatorios son m??s eficientes que los modelos de efectos fijos bajo condiciones de exogeneidad de los efectos individuales. La prueba de Hausman nos sirve para realizar pruebas de endogeneidad:
phtest(gr_re, gr_fe)##
## Hausman Test
##
## data: invest ~ value + capital
## chisq = 0.074639, df = 2, p-value = 0.9634
## alternative hypothesis: one model is inconsistentPresentamos el modelo dinámico de datos en panel: \[ Y_{it} = \alpha_i + \beta X_{it} + \sum_{k=1}^p\gamma_k Y_{i,t-k}+\varepsilon_{it} \]
Este modelo presenta además una serie de rezagos, que permiten expresar la idea de que el valor de \(Y_{it}\) depende de los valores previos que se hayan observado. La elecci??n de \(p\) implica expresar de cuántos valores de periodos anteriores depende nuestra variable de interéss. Considere los datos de Arellano-Bond (1991) sobre las determinantes del empleo en un panel de 140 empresas del Reino Unido para el periodo de 1976 a 1984, contenidos en el paquete plm con fines ilustrativos. Estos datos componen un panel no-balanceado, con entre siete y nueve observaciones por empresa. Comenzaremos con una f??rmula sencilla en el que se consideran variables relevantes: salario promedio anual por empleado wage, el valor en libros del capital fijo capital y un índice de valor agregado considerando un costo de factores constante output, todo en logaritmos:
data("EmplUK", package = "plm")
form <- log(emp) ~ log(wage) + log(capital) + log(output)
head(EmplUK)## firm year sector emp wage capital output
## 1 1 1977 7 5.041 13.1516 0.5894 95.7072
## 2 1 1978 7 5.600 12.3018 0.6318 97.3569
## 3 1 1979 7 5.015 12.8395 0.6771 99.6083
## 4 1 1980 7 4.715 13.8039 0.6171 100.5501
## 5 1 1981 7 4.093 14.2897 0.5076 99.5581
## 6 1 1982 7 3.166 14.8681 0.4229 98.6151Usaremos la función pgmm() (GMM proviene de Generalized Method of Moments) para realizar la estimación tomando como primer argumento una fórmula dinámica, expresada por medio de la función dynformula. Esta función toma como argumentos la fórmula que expresamos arriba y una lista con enteros que se interpreta como el número de rezagos de cada variable. En otras palabras,dynformula(form, list(2,1,0,1)) significa que se utilizará una fórmula dinámica como la expresada con la variable form que hicimos arriba, donde el empleo emp toma dos rezagos, los sueldos wage un rezago y así respectivamente. El estimador de Arellano y Bond (1991) y su aplicación dinámica a los datos de empleo se estiman de la siguiente manera:
empl_ab <- pgmm(dynformula(form, list(2,1,0,1)),
data = EmplUK, index = c("firm", "year"),
effect = "twoways", model = "twosteps",
gmm.inst = ~ log(emp), lag.gmm = list(c(2,99)))El argumento lag.gmm = list(c(2,99)) significa que todos los rezagos más allá de 1 se han de usar como instrumentos. Mostramos el resultado de esta regresión.
summary(empl_ab)## Twoways effects Two steps model
##
## Call:
## pgmm(formula = dynformula(form, list(2, 1, 0, 1)), data = EmplUK,
## effect = "twoways", model = "twosteps", index = c("firm",
## "year"), gmm.inst = ~log(emp), lag.gmm = list(c(2, 99)))
##
## Unbalanced Panel: n = 140, T = 7-9, N = 1031
##
## Number of Observations Used: 611
##
## Residuals:
## Min. 1st Qu. Median Mean 3rd Qu. Max.
## -0.6191000 -0.0255700 0.0000000 -0.0001339 0.0332000 0.6410000
##
## Coefficients:
## Estimate Std. Error z-value Pr(>|z|)
## lag(log(emp), c(1, 2))1 0.474151 0.185398 2.5575 0.0105437 *
## lag(log(emp), c(1, 2))2 -0.052967 0.051749 -1.0235 0.3060506
## log(wage) -0.513205 0.145565 -3.5256 0.0004225 ***
## lag(log(wage), 1) 0.224640 0.141950 1.5825 0.1135279
## log(capital) 0.292723 0.062627 4.6741 2.953e-06 ***
## log(output) 0.609775 0.156263 3.9022 9.530e-05 ***
## lag(log(output), 1) -0.446373 0.217302 -2.0542 0.0399605 *
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Sargan test: chisq(25) = 30.11247 (p-value = 0.22011)
## Autocorrelation test (1): normal = -1.53845 (p-value = 0.12394)
## Autocorrelation test (2): normal = -0.2796829 (p-value = 0.77972)
## Wald test for coefficients: chisq(7) = 142.0353 (p-value = < 2.22e-16)
## Wald test for time dummies: chisq(6) = 16.97046 (p-value = 0.0093924)Las pruebas muestran dos cosas: la primera es que los rezagos son útiles y mejoran el modelo sustancialmente. Las razones de esto pueden variar y dependen mucho de la dinámica existente entre los datos. Lo segundo que se puede observar es que este modelo es insuficiente, pero por motivos de exposición se deja como ejercicio al lector la especificación de Arellano y Bond, que toma a los salarios y al capital como variables endógenas.