Econometria: exemplos 10.5 e 12.3 do Wooldridge (2016, p.390)

Abstract

This is an undergrad student level exercise for class use. We analyse data on Barium market, from 1978m02 to 1988m12.

Licença

This work is licensed under the Creative Commons Attribution-ShareAlike 4.0 International License. To view a copy of this license, visit http://creativecommons.org/licenses/by-sa/4.0/ or send a letter to Creative Commons, PO Box 1866, Mountain View, CA 94042, USA.

License: CC BY-SA 4.0

Citação

Sugestão de citação: FIGUEIREDO, Adriano Marcos Rodrigues. Econometria: exemplos 10.5 e 12.3 do Wooldridge (2016, p.390). Campo Grande-MS,Brasil: RStudio/Rpubs, 2022. Disponível em http://rpubs.com/amrofi/ex_wooldridge_12_3_autocorrelation.

12.3. Wooldridge Source: C.M. Krupp and P.S. Pollard (1999), Market Responses to Antidumpting Laws: Some Evidence from the U.S. Chemical Industry, Canadian Journal of Economics 29, 199-227. They are monthly data covering February 1978 through December 1988.

A data.frame with 131 observations on 31 variables:

• chnimp: Chinese imports, bar. chl.

• bchlimp: total imports bar. chl.

• befile6: =1 for all 6 mos before filing

• affile6: =1 for all 6 mos after filing

• afdec6: =1 for all 6 mos after decision

• befile12: =1 all 12 mos before filing

• affile12: =1 all 12 mos after filing

• afdec12: =1 all 12 mos after decision

• chempi: chemical production index

• gas: gasoline production

• rtwex: exchange rate index

• spr: =1 for spring months

• sum: =1 for summer months

• fall: =1 for fall months

• lchnimp: log(chnimp)

• lgas: log(gas)

• lrtwex: log(rtwex)

• lchempi: log(chempi)

• t: time trend

• feb: =1 if month is feb

• mar: =1 if month is march

• apr:

• may:

• jun:

• jul:

• aug:

• sep:

• oct:

• nov:

• dec:

• percchn: percent imports from china

Ver também https://justinmshea.github.io/wooldridge/articles/Introductory-Econometrics-Examples.html e http://www.urfie.net/downloads/12/Example-12-3.R.

library(dynlm)
library(stargazer)
library(wooldridge)
data(barium, package = "wooldridge")
dados = barium
attach(dados)

Vamos ver como está a tabela importada:

library(knitr)
kable(dados[, 1:5], caption = "Dados do exemplo")

Dados do exemplo

chnimp	bchlimp	befile6	affile6	afdec6
220.4620	9578.376	0	0	0
94.7980	11219.480	0	0	0
219.3575	9719.900	0	0	0
317.4215	12920.950	0	0	0
114.6390	9790.446	0	0	0
129.5240	11020.470	0	0	0
110.0920	8302.168	0	0	0
39.6830	9583.325	0	0	0
262.8230	9880.729	0	0	0
173.9980	9003.923	0	0	0
548.9455	10153.390	0	0	0
433.9755	10527.060	0	0	0
213.7855	9238.292	0	0	0
788.6440	11740.390	0	0	0
287.0395	13699.070	0	0	0
278.9545	17716.330	0	0	0
429.8970	15868.760	0	0	0
344.8820	15512.320	0	0	0
538.9150	13831.560	0	0	0
535.7180	13669.760	0	0	0
445.3290	16691.211	0	0	0
649.2550	16617.330	0	0	0
639.3350	15912.760	0	0	0
209.4370	9401.297	0	0	0
512.5915	9195.309	0	0	0
605.3425	11727.520	0	0	0
144.2360	20584.260	0	0	0
405.5360	8929.146	0	0	0
621.1075	10912.110	0	0	0
893.7185	13581.910	0	0	0
140.2300	9218.409	0	0	0
243.4765	20751.770	0	0	0
330.5750	12935.290	0	0	0
46.5860	26984.199	0	0	0
215.1960	18886.381	0	0	0
214.8385	25679.510	0	0	0
238.2085	15502.940	0	0	0
265.6545	18247.141	0	0	0
315.2515	12060.010	0	0	0
426.5905	17777.609	0	0	0
124.7860	14406.150	0	0	0
387.8990	20570.359	0	0	0
400.5400	14210.780	0	0	0
876.8340	15190.780	0	0	0
535.9930	14593.670	0	0	0
815.9555	14044.220	0	0	0
771.8130	10431.430	0	0	0
352.8730	13228.840	0	0	0
431.9340	12801.710	0	0	0
1211.4520	16403.340	0	0	0
575.9175	13591.700	0	0	0
1491.6219	16239.280	0	0	0
380.2920	20010.461	0	0	0
391.5965	16873.131	0	0	0
505.5120	15379.570	0	0	0
341.0220	10551.650	0	0	0
372.5765	20755.381	0	0	0
483.3590	18066.199	0	0	0
224.0470	17886.760	0	0	0
532.3940	19134.020	0	0	0
239.1405	17455.141	0	0	0
689.4535	18215.721	0	0	0
318.0120	20764.830	1	0	0
1048.0360	28414.789	1	0	0
945.7295	18819.391	1	0	0
295.3720	11737.220	1	0	0
319.4580	12520.280	1	0	0
855.1370	16988.230	1	0	0
671.8070	15835.780	0	0	0
657.9315	17282.949	0	1	0
847.9595	26960.510	0	1	0
849.8735	25137.689	0	1	0
852.9305	19807.490	0	1	0
435.8370	32134.600	0	1	0
204.7075	22954.551	0	1	0
241.5675	24438.789	0	0	0
244.1380	19057.439	0	0	0
656.4015	26784.340	0	0	0
588.0975	37244.379	0	0	0
1128.4440	28129.359	0	0	0
1328.9670	37253.469	0	0	0
227.0635	18289.359	0	0	1
583.2610	14795.910	0	0	1
499.8935	25700.260	0	0	1
165.3460	25240.311	0	0	1
1057.9720	24854.939	0	0	1
465.8225	24119.289	0	0	1
349.4040	27244.980	0	0	0
466.8240	30912.420	0	0	0
717.0470	23173.529	0	0	0
826.7580	29031.439	0	0	0
1571.7360	19400.510	0	0	0
542.9925	29887.711	0	0	0
670.7935	20670.631	0	0	0
430.2755	34375.250	0	0	0
524.8625	19555.820	0	0	0
317.7500	27362.230	0	0	0
442.8195	31952.500	0	0	0
322.6325	19578.160	0	0	0
742.0795	26605.789	0	0	0
665.3660	20034.420	0	0	0
737.2670	33408.488	0	0	0
1065.2111	20967.311	0	0	0
547.3285	26526.840	0	0	0
1258.4070	28304.020	0	0	0
976.5165	17168.150	0	0	0
1192.6840	12874.300	0	0	0
1047.2390	16153.360	0	0	0
944.7875	15638.680	0	0	0
840.0295	18391.830	0	0	0
578.9175	18786.789	0	0	0
828.1195	17661.711	0	0	0
1010.0570	17623.770	0	0	0
1546.4830	18482.400	0	0	0
1160.1580	21132.420	0	0	0
862.6800	16421.520	0	0	0
692.2970	16128.660	0	0	0
451.3840	14769.440	0	0	0
766.4015	34192.180	0	0	0
1794.5160	17537.230	0	0	0
972.8010	25518.910	0	0	0
738.5085	36649.191	0	0	0
786.1805	18668.000	0	0	0
849.0650	19710.660	0	0	0
895.4205	38521.531	0	0	0
814.3615	19679.381	0	0	0
589.4220	31686.859	0	0	0
529.7335	16370.340	0	0	0
968.0375	15653.960	0	0	0
794.0600	45716.609	0	0	0
1087.8060	29580.500	0	0	0

Estimando o modelo linear de regressao multipla fazendo conforme a expressão do enunciado.

1 Resultados

1.1 Estimação

Fazendo as regressoes com uso de logaritmos.

tsdata <- ts(dados, start = c(1978, 2), frequency = 12)
set.seed(12345)
mod1 <- dynlm(log(chnimp) ~ log(chempi) + log(gas) + log(rtwex) + befile6 + affile6 +
    afdec6, data = tsdata)

Vamos utilizar o pacote stargazer para organizar as saídas de resultados. Se a saída fosse apenas pelo comando summary, sairia da forma:

summary(mod1)

Time series regression with "ts" data:
Start = 1978(2), End = 1988(12)

Call:
dynlm(formula = log(chnimp) ~ log(chempi) + log(gas) + log(rtwex) + 
    befile6 + affile6 + afdec6, data = tsdata)

Residuals:
     Min       1Q   Median       3Q      Max 
-2.03356 -0.39080  0.03048  0.40248  1.51720 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept) -17.80277   21.04537  -0.846   0.3992    
log(chempi)   3.11719    0.47920   6.505 1.72e-09 ***
log(gas)      0.19634    0.90662   0.217   0.8289    
log(rtwex)    0.98302    0.40015   2.457   0.0154 *  
befile6       0.05957    0.26097   0.228   0.8198    
affile6      -0.03241    0.26430  -0.123   0.9026    
afdec6       -0.56524    0.28584  -1.978   0.0502 .  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.5974 on 124 degrees of freedom
Multiple R-squared:  0.3049,    Adjusted R-squared:  0.2712 
F-statistic: 9.064 on 6 and 124 DF,  p-value: 3.255e-08

Agora, criando uma tabela com as várias saídas de modelos, com o pacote stargazer tem-se, com a geração de AIC e BIC:

library(stargazer)
mod1$AIC <- AIC(mod1)

mod1$BIC <- BIC(mod1)

library(stargazer)
star.1 <- stargazer(mod1, title = "Título: Resultado da Regressão", align = TRUE,
    type = "text", style = "all", keep.stat = c("aic", "bic", "rsq", "adj.rsq", "n"))

Título: Resultado da Regressão
===============================================
                        Dependent variable:    
                    ---------------------------
                            log(chnimp)        
-----------------------------------------------
log(chempi)                  3.117***          
                              (0.479)          
                             t = 6.505         
                             p = 0.000         
log(gas)                       0.196           
                              (0.907)          
                             t = 0.217         
                             p = 0.829         
log(rtwex)                    0.983**          
                              (0.400)          
                             t = 2.457         
                             p = 0.016         
befile6                        0.060           
                              (0.261)          
                             t = 0.228         
                             p = 0.820         
affile6                       -0.032           
                              (0.264)          
                            t = -0.123         
                             p = 0.903         
afdec6                        -0.565*          
                              (0.286)          
                            t = -1.978         
                             p = 0.051         
Constant                      -17.803          
                             (21.045)          
                            t = -0.846         
                             p = 0.400         
-----------------------------------------------
Observations                    131            
R2                             0.305           
Adjusted R2                    0.271           
Akaike Inf. Crit.             245.573          
Bayesian Inf. Crit.           268.575          
===============================================
Note:               *p<0.1; **p<0.05; ***p<0.01

2 Autocorrelação dos resíduos

library(car)
library(lmtest)
library(sandwich)
n <- nrow(dados)  #armazenando o valor de n=numero de observacoes
k <- length(coefficients(mod1)) - 1  #obtendo o numero de regressores menos o intercepto
c("número de observações: n=", n)

[1] "número de observações: n=" "131"                      

c("número de variáveis explicativas k=", k)

[1] "número de variáveis explicativas k=" "6"                                  

2.1 Breusch-Godfrey LM Serial Correlation test:

# Exemplo Wooldridge 12.3:
set.seed(12345)
library(car)
library(dynlm)
residual <- resid(mod1)

resreg <- dynlm::dynlm(residual ~ L(residual) + L(residual, 2) + L(residual, 3) +
    log(chempi) + log(gas) + log(rtwex) + befile6 + affile6 + afdec6, data = tsdata)
summary(resreg)

Time series regression with "ts" data:
Start = 1978(5), End = 1988(12)

Call:
dynlm::dynlm(formula = residual ~ L(residual) + L(residual, 2) + 
    L(residual, 3) + log(chempi) + log(gas) + log(rtwex) + befile6 + 
    affile6 + afdec6, data = tsdata)

Residuals:
     Min       1Q   Median       3Q      Max 
-1.89072 -0.32250  0.05873  0.36376  1.19650 

Coefficients:
                Estimate Std. Error t value Pr(>|t|)  
(Intercept)    -14.36906   20.65566  -0.696   0.4880  
L(residual)      0.22149    0.09166   2.417   0.0172 *
L(residual, 2)   0.13404    0.09216   1.454   0.1485  
L(residual, 3)   0.12554    0.09112   1.378   0.1709  
log(chempi)     -0.14316    0.47203  -0.303   0.7622  
log(gas)         0.62330    0.88597   0.704   0.4831  
log(rtwex)       0.17867    0.39103   0.457   0.6486  
befile6         -0.08592    0.25101  -0.342   0.7327  
affile6         -0.12212    0.25470  -0.479   0.6325  
afdec6          -0.06683    0.27437  -0.244   0.8080  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.5704 on 118 degrees of freedom
Multiple R-squared:  0.1159,    Adjusted R-squared:  0.04848 
F-statistic: 1.719 on 9 and 118 DF,  p-value: 0.09202

2.1.0.1 Teste F e Qui-quadrado pela função linearHypothesis

Atentar que a saída em html é formatada pelo tidyr e mostra-se um pouco diferente do formato original obtido no console R.

linearHypothesis(resreg, test = c("Chisq"), c("L(residual)", "L(residual, 2)", "L(residual, 3)"))

linearHypothesis(resreg, test = c("F"), c("L(residual)", "L(residual, 2)", "L(residual, 3)"))

r2aux <- summary(resreg)$r.squared
r2aux  # R2 da equação de teste

[1] 0.1159108

nr2 <- n * r2aux  # Estatística de teste Qui-Quadrado estimado
nr2

[1] 15.18431

2.1.1 Opção do Breusch-Godfrey LM Test pela função lmtest::bgtest

bgtest3.chi <- bgtest(mod1, order = 3, type = c("Chisq"), data = tsdata)
bgtest3.chi  # padrao do teste de BG, com distribuição qui-quadrado

    Breusch-Godfrey test for serial correlation of order up to 3

data:  mod1
LM test = 14.768, df = 3, p-value = 0.002026

bgtest3.F <- bgtest(mod1, order = 3, type = c("F"), data = tsdata)
bgtest3.F  #alternativa do teste de BG com versão para amostra finita

    Breusch-Godfrey test for serial correlation of order up to 3

data:  mod1
LM test = 5.1247, df1 = 3, df2 = 121, p-value = 0.002264

Enxergando a equação de teste do bgtest:

lmtest::coeftest(bgtest3.chi)

z test of coefficients:

               Estimate Std. Error z value Pr(>|z|)  
(Intercept)  -13.140690  20.350979 -0.6457   0.5185  
log(chempi)   -0.043471   0.457134 -0.0951   0.9242  
log(gas)       0.543595   0.875736  0.6207   0.5348  
log(rtwex)     0.202386   0.385555  0.5249   0.5996  
befile6       -0.081533   0.250628 -0.3253   0.7449  
affile6       -0.124293   0.254311 -0.4887   0.6250  
afdec6        -0.077584   0.273697 -0.2835   0.7768  
lag(resid)_1   0.218596   0.090487  2.4158   0.0157 *
lag(resid)_2   0.134116   0.091995  1.4579   0.1449  
lag(resid)_3   0.125677   0.090972  1.3815   0.1671  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

2.2 Durbin-Watson Test (AR(1))

dw.mod1 <- dwtest(mod1)
dw.mod1

    Durbin-Watson test

data:  mod1
DW = 1.4584, p-value = 0.0001688
alternative hypothesis: true autocorrelation is greater than 0

Fiz uma rotina para rodar vários BGtest até ordem 12.

# padrao do teste de BG, com distribuição qui-quadrado
library(lmtest)
bgorder = 1:12  # definindo até a máxima ordem do bgtest
d = NULL
for (p in bgorder) {
    bgtest.chi <- bgtest(mod1, order = p, type = c("Chisq"), data = dados)
    print(bgtest.chi)
    d = rbind(d, data.frame(bgtest.chi$statistic, bgtest.chi$p.value))
    print(coeftest(bgtest.chi))
}

    Breusch-Godfrey test for serial correlation of order up to 1

data:  mod1
LM test = 9.8291, df = 1, p-value = 0.001718


z test of coefficients:

              Estimate Std. Error z value Pr(>|z|)   
(Intercept)  -7.918439  20.476585 -0.3867 0.698973   
log(chempi)  -0.016122   0.462772 -0.0348 0.972209   
log(gas)      0.329089   0.881656  0.3733 0.708953   
log(rtwex)    0.103582   0.387799  0.2671 0.789391   
befile6      -0.026860   0.252150 -0.1065 0.915167   
affile6      -0.065503   0.256061 -0.2558 0.798096   
afdec6       -0.039794   0.276305 -0.1440 0.885484   
lag(resid)_1  0.277166   0.087746  3.1587 0.001585 **
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1


    Breusch-Godfrey test for serial correlation of order up to 2

data:  mod1
LM test = 12.935, df = 2, p-value = 0.001553


z test of coefficients:

               Estimate Std. Error z value Pr(>|z|)   
(Intercept)  -11.090299  20.372219 -0.5444 0.586178   
log(chempi)   -0.031910   0.458756 -0.0696 0.944545   
log(gas)       0.461069   0.876944  0.5258 0.599049   
log(rtwex)     0.154146   0.385398  0.4000 0.689182   
befile6       -0.044374   0.250107 -0.1774 0.859178   
affile6       -0.111486   0.255087 -0.4371 0.662074   
afdec6        -0.048485   0.273900 -0.1770 0.859495   
lag(resid)_1   0.236378   0.089900  2.6294 0.008555 **
lag(resid)_2   0.161523   0.090164  1.7914 0.073224 . 
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1


    Breusch-Godfrey test for serial correlation of order up to 3

data:  mod1
LM test = 14.768, df = 3, p-value = 0.002026


z test of coefficients:

               Estimate Std. Error z value Pr(>|z|)  
(Intercept)  -13.140690  20.350979 -0.6457   0.5185  
log(chempi)   -0.043471   0.457134 -0.0951   0.9242  
log(gas)       0.543595   0.875736  0.6207   0.5348  
log(rtwex)     0.202386   0.385555  0.5249   0.5996  
befile6       -0.081533   0.250628 -0.3253   0.7449  
affile6       -0.124293   0.254311 -0.4887   0.6250  
afdec6        -0.077584   0.273697 -0.2835   0.7768  
lag(resid)_1   0.218596   0.090487  2.4158   0.0157 *
lag(resid)_2   0.134116   0.091995  1.4579   0.1449  
lag(resid)_3   0.125677   0.090972  1.3815   0.1671  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1


    Breusch-Godfrey test for serial correlation of order up to 4

data:  mod1
LM test = 17.463, df = 4, p-value = 0.001571


z test of coefficients:

               Estimate Std. Error z value Pr(>|z|)   
(Intercept)  -8.4647011 20.3864595 -0.4152 0.677987   
log(chempi)   0.0035169  0.4545352  0.0077 0.993826   
log(gas)      0.3380107  0.8776186  0.3851 0.700130   
log(rtwex)    0.1590660  0.3835037  0.4148 0.678310   
befile6      -0.0780067  0.2487446 -0.3136 0.753824   
affile6      -0.1479439  0.2527796 -0.5853 0.558367   
afdec6       -0.0716150  0.2716537 -0.2636 0.792068   
lag(resid)_1  0.2355629  0.0903646  2.6068 0.009139 **
lag(resid)_2  0.1524613  0.0919449  1.6582 0.097281 . 
lag(resid)_3  0.1568731  0.0921576  1.7022 0.088713 . 
lag(resid)_4 -0.1544626  0.0915235 -1.6877 0.091472 . 
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1


    Breusch-Godfrey test for serial correlation of order up to 5

data:  mod1
LM test = 17.465, df = 5, p-value = 0.003697


z test of coefficients:

               Estimate Std. Error z value Pr(>|z|)  
(Intercept)  -8.4989457 20.4840385 -0.4149  0.67821  
log(chempi)   0.0022123  0.4572368  0.0048  0.99614  
log(gas)      0.3398494  0.8821133  0.3853  0.70004  
log(rtwex)    0.1587066  0.3851800  0.4120  0.68032  
befile6      -0.0775195  0.2499892 -0.3101  0.75649  
affile6      -0.1474423  0.2540498 -0.5804  0.56167  
afdec6       -0.0704850  0.2737936 -0.2574  0.79684  
lag(resid)_1  0.2362415  0.0918264  2.5727  0.01009 *
lag(resid)_2  0.1517863  0.0933835  1.6254  0.10408  
lag(resid)_3  0.1561947  0.0936052  1.6687  0.09519 .
lag(resid)_4 -0.1554378  0.0941029 -1.6518  0.09858 .
lag(resid)_5  0.0044694  0.0926364  0.0482  0.96152  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1


    Breusch-Godfrey test for serial correlation of order up to 6

data:  mod1
LM test = 17.555, df = 6, p-value = 0.007445


z test of coefficients:

               Estimate Std. Error z value Pr(>|z|)  
(Intercept)  -8.7383883 20.5773853 -0.4247  0.67108  
log(chempi)  -0.0054375  0.4596686 -0.0118  0.99056  
log(gas)      0.3515545  0.8863181  0.3966  0.69163  
log(rtwex)    0.1604875  0.3866991  0.4150  0.67813  
befile6      -0.0780190  0.2509520 -0.3109  0.75588  
affile6      -0.1503270  0.2551971 -0.5891  0.55582  
afdec6       -0.0639994  0.2756586 -0.2322  0.81641  
lag(resid)_1  0.2361406  0.0921787  2.5618  0.01041 *
lag(resid)_2  0.1563668  0.0949291  1.6472  0.09952 .
lag(resid)_3  0.1518751  0.0950185  1.5984  0.10996  
lag(resid)_4 -0.1594817  0.0953836 -1.6720  0.09452 .
lag(resid)_5 -0.0017308  0.0951737 -0.0182  0.98549  
lag(resid)_6  0.0284533  0.0929982  0.3060  0.75964  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1


    Breusch-Godfrey test for serial correlation of order up to 7

data:  mod1
LM test = 18.881, df = 7, p-value = 0.008569


z test of coefficients:

               Estimate Std. Error z value Pr(>|z|)  
(Intercept)  -7.6947819 20.5632048 -0.3742  0.70825  
log(chempi)  -0.0185979  0.4590604 -0.0405  0.96768  
log(gas)      0.3157464  0.8854060  0.3566  0.72138  
log(rtwex)    0.1245003  0.3872833  0.3215  0.74785  
befile6      -0.0648778  0.2507945 -0.2587  0.79588  
affile6      -0.1413088  0.2548990 -0.5544  0.57932  
afdec6       -0.0036094  0.2799610 -0.0129  0.98971  
lag(resid)_1  0.2327739  0.0920739  2.5281  0.01147 *
lag(resid)_2  0.1560081  0.0947759  1.6461  0.09975 .
lag(resid)_3  0.1673345  0.0957709  1.7472  0.08060 .
lag(resid)_4 -0.1774415  0.0964456 -1.8398  0.06580 .
lag(resid)_5 -0.0168956  0.0958904 -0.1762  0.86014  
lag(resid)_6  0.0049764  0.0949691  0.0524  0.95821  
lag(resid)_7  0.1107054  0.0941282  1.1761  0.23955  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1


    Breusch-Godfrey test for serial correlation of order up to 8

data:  mod1
LM test = 18.9, df = 8, p-value = 0.0154


z test of coefficients:

               Estimate Std. Error z value Pr(>|z|)  
(Intercept)  -7.7001231 20.6498876 -0.3729  0.70923  
log(chempi)  -0.0206594  0.4612232 -0.0448  0.96427  
log(gas)      0.3169018  0.8891741  0.3564  0.72154  
log(rtwex)    0.1219837  0.3893186  0.3133  0.75403  
befile6      -0.0653532  0.2518735 -0.2595  0.79527  
affile6      -0.1402735  0.2560769 -0.5478  0.58384  
afdec6        0.0038994  0.2860676  0.0136  0.98912  
lag(resid)_1  0.2312244  0.0931031  2.4835  0.01301 *
lag(resid)_2  0.1560350  0.0951754  1.6394  0.10112  
lag(resid)_3  0.1674233  0.0961765  1.7408  0.08172 .
lag(resid)_4 -0.1752485  0.0980749 -1.7869  0.07396 .
lag(resid)_5 -0.0190656  0.0974989 -0.1955  0.84496  
lag(resid)_6  0.0030772  0.0963021  0.0320  0.97451  
lag(resid)_7  0.1081650  0.0962020  1.1244  0.26086  
lag(resid)_8  0.0135157  0.0951539  0.1420  0.88705  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1


    Breusch-Godfrey test for serial correlation of order up to 9

data:  mod1
LM test = 22.435, df = 9, p-value = 0.007598


z test of coefficients:

              Estimate Std. Error z value Pr(>|z|)  
(Intercept)  -4.940925  20.459635 -0.2415  0.80917  
log(chempi)   0.049523   0.457303  0.1083  0.91376  
log(gas)      0.183517   0.881539  0.2082  0.83509  
log(rtwex)    0.112900   0.384822  0.2934  0.76923  
befile6      -0.067963   0.248950 -0.2730  0.78485  
affile6      -0.149620   0.253147 -0.5910  0.55449  
afdec6       -0.066031   0.285043 -0.2317  0.81681  
lag(resid)_1  0.233588   0.092029  2.5382  0.01114 *
lag(resid)_2  0.174727   0.094564  1.8477  0.06464 .
lag(resid)_3  0.168145   0.095059  1.7688  0.07692 .
lag(resid)_4 -0.179427   0.096959 -1.8505  0.06424 .
lag(resid)_5 -0.051431   0.097806 -0.5258  0.59900  
lag(resid)_6  0.032043   0.096353  0.3326  0.73947  
lag(resid)_7  0.132182   0.095891  1.3785  0.16806  
lag(resid)_8  0.048850   0.095804  0.5099  0.61013  
lag(resid)_9 -0.180887   0.093482 -1.9350  0.05299 .
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1


    Breusch-Godfrey test for serial correlation of order up to 10

data:  mod1
LM test = 22.438, df = 10, p-value = 0.01302


z test of coefficients:

                Estimate Std. Error z value Pr(>|z|)  
(Intercept)   -4.9894613 20.5679780 -0.2426  0.80833  
log(chempi)    0.0512483  0.4603762  0.1113  0.91136  
log(gas)       0.1853512  0.8860171  0.2092  0.83430  
log(rtwex)     0.1125301  0.3865603  0.2911  0.77097  
befile6       -0.0673127  0.2503171 -0.2689  0.78800  
affile6       -0.1482142  0.2555418 -0.5800  0.56191  
afdec6        -0.0667997  0.2866303 -0.2331  0.81572  
lag(resid)_1   0.2326780  0.0939092  2.4777  0.01322 *
lag(resid)_2   0.1749524  0.0950652  1.8403  0.06572 .
lag(resid)_3   0.1688346  0.0962987  1.7532  0.07956 .
lag(resid)_4  -0.1792313  0.0974476 -1.8393  0.06588 .
lag(resid)_5  -0.0516759  0.0983351 -0.5255  0.59923  
lag(resid)_6   0.0310887  0.0983269  0.3162  0.75187  
lag(resid)_7   0.1330110  0.0974882  1.3644  0.17245  
lag(resid)_8   0.0496699  0.0973797  0.5101  0.61001  
lag(resid)_9  -0.1797681  0.0960855 -1.8709  0.06136 .
lag(resid)_10 -0.0052265  0.0953616 -0.0548  0.95629  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1


    Breusch-Godfrey test for serial correlation of order up to 11

data:  mod1
LM test = 23.335, df = 11, p-value = 0.01585


z test of coefficients:

               Estimate Std. Error z value Pr(>|z|)  
(Intercept)   -1.609972  20.866053 -0.0772  0.93850  
log(chempi)    0.108207   0.464222  0.2331  0.81569  
log(gas)       0.033154   0.900021  0.0368  0.97061  
log(rtwex)     0.074641   0.388627  0.1921  0.84769  
befile6       -0.034812   0.252612 -0.1378  0.89039  
affile6       -0.141119   0.255712 -0.5519  0.58104  
afdec6        -0.067337   0.286705 -0.2349  0.81431  
lag(resid)_1   0.230696   0.093956  2.4554  0.01407 *
lag(resid)_2   0.157690   0.096740  1.6300  0.10309  
lag(resid)_3   0.172454   0.096396  1.7890  0.07361 .
lag(resid)_4  -0.168216   0.098132 -1.7142  0.08649 .
lag(resid)_5  -0.049656   0.098382 -0.5047  0.61375  
lag(resid)_6   0.025333   0.098531  0.2571  0.79710  
lag(resid)_7   0.117270   0.098854  1.1863  0.23550  
lag(resid)_8   0.065279   0.098724  0.6612  0.50847  
lag(resid)_9  -0.165277   0.097264 -1.6993  0.08927 .
lag(resid)_10  0.014987   0.097635  0.1535  0.87800  
lag(resid)_11 -0.093949   0.096832 -0.9702  0.33193  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1


    Breusch-Godfrey test for serial correlation of order up to 12

data:  mod1
LM test = 23.935, df = 12, p-value = 0.02076


z test of coefficients:

                Estimate Std. Error z value Pr(>|z|)  
(Intercept)   -1.7226416 20.9009656 -0.0824   0.9343  
log(chempi)    0.1286140  0.4657002  0.2762   0.7824  
log(gas)       0.0355733  0.9015113  0.0395   0.9685  
log(rtwex)     0.0652576  0.3894481  0.1676   0.8669  
befile6       -0.0017098  0.2564538 -0.0067   0.9947  
affile6       -0.1387581  0.2561514 -0.5417   0.5880  
afdec6        -0.0525709  0.2877815 -0.1827   0.8551  
lag(resid)_1   0.2231724  0.0945882  2.3594   0.0183 *
lag(resid)_2   0.1582820  0.0969025  1.6334   0.1024  
lag(resid)_3   0.1587020  0.0981016  1.6177   0.1057  
lag(resid)_4  -0.1620733  0.0985987 -1.6438   0.1002  
lag(resid)_5  -0.0388890  0.0994769 -0.3909   0.6958  
lag(resid)_6   0.0268346  0.0987117  0.2718   0.7857  
lag(resid)_7   0.1136053  0.0991247  1.1461   0.2518  
lag(resid)_8   0.0530359  0.1000868  0.5299   0.5962  
lag(resid)_9  -0.1508785  0.0991041 -1.5224   0.1279  
lag(resid)_10  0.0269707  0.0989581  0.2725   0.7852  
lag(resid)_11 -0.0791453  0.0987739 -0.8013   0.4230  
lag(resid)_12 -0.0774090  0.0976794 -0.7925   0.4281  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

d

A solução passa por alternativas na modelagem de AR e MA ou procedimentos de Cochrane-Orcutt.

2.3 Cochrane-Orcutt MQG Factível

library(orcutt)

# OLS estimation

olsres <- dynlm(log(chnimp) ~ log(chempi) + log(gas) + log(rtwex) + befile6 + affile6 +
    afdec6, data = tsdata)

# Cochrane-Orcutt estimation
cochrane.orcutt(olsres)

Cochrane-orcutt estimation for first order autocorrelation 
 
Call:
dynlm(formula = log(chnimp) ~ log(chempi) + log(gas) + log(rtwex) + 
    befile6 + affile6 + afdec6, data = tsdata)

 number of interaction: 8
 rho 0.293362

Durbin-Watson statistic 
(original):    1.45841 , p-value: 1.688e-04
(transformed): 2.06330 , p-value: 4.91e-01
 
 coefficients: 
(Intercept) log(chempi)    log(gas)  log(rtwex)     befile6     affile6 
 -37.322241    2.947434    1.054858    1.136918   -0.016372   -0.033082 
     afdec6 
  -0.577158 

2.4 Prais-Winsten MQG Factível

Carregar o pacote prais e usar a função prais_winsten para estimar.

library(prais)
barium_prais_winsten <- prais_winsten(log(chnimp) ~ log(chempi) + log(gas) + log(rtwex) +
    befile6 + affile6 + afdec6, data = tsdata, index = "t")

Iteration 0: rho = 0
Iteration 1: rho = 0.2708
Iteration 2: rho = 0.291
Iteration 3: rho = 0.293
Iteration 4: rho = 0.2932
Iteration 5: rho = 0.2932
Iteration 6: rho = 0.2932
Iteration 7: rho = 0.2932

summary(barium_prais_winsten)

Call:
prais_winsten(formula = log(chnimp) ~ log(chempi) + log(gas) + 
    log(rtwex) + befile6 + affile6 + afdec6, data = tsdata, index = "t")

Residuals:
     Min       1Q   Median       3Q      Max 
-1.99386 -0.32219  0.03747  0.40226  1.50281 

AR(1) coefficient rho after 7 iterations: 0.2932

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept) -37.07742   22.77831  -1.628   0.1061    
log(chempi)   2.94095    0.63284   4.647 8.46e-06 ***
log(gas)      1.04637    0.97734   1.071   0.2864    
log(rtwex)    1.13279    0.50666   2.236   0.0272 *  
befile6      -0.01648    0.31938  -0.052   0.9589    
affile6      -0.03316    0.32181  -0.103   0.9181    
afdec6       -0.57681    0.34199  -1.687   0.0942 .  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.5733 on 124 degrees of freedom
Multiple R-squared:  0.2021,    Adjusted R-squared:  0.1635 
F-statistic: 5.235 on 6 and 124 DF,  p-value: 7.764e-05

Durbin-Watson statistic (original): 1.458 
Durbin-Watson statistic (transformed): 2.087

2.5 Inferência Robusta - Newey-West

Conforme a ordem da autocorelação, pode ser mais indicado considerar o tratamento para a estrutura geral de Newey e West (1987) e Wooldridge (1989), citados por Wooldridge (2016, seção 12.5).

# library(sandwich) OLS regression foi feita a mod1
mod1 <- dynlm(log(chnimp) ~ log(chempi) + log(gas) + log(rtwex) + befile6 + affile6 +
    afdec6, data = tsdata)
# resultados com erros-padrões usuais
coeftest(mod1)

t test of coefficients:

              Estimate Std. Error t value  Pr(>|t|)    
(Intercept) -17.802768  21.045367 -0.8459   0.39922    
log(chempi)   3.117194   0.479202  6.5050 1.724e-09 ***
log(gas)      0.196341   0.906617  0.2166   0.82890    
log(rtwex)    0.983016   0.400154  2.4566   0.01541 *  
befile6       0.059574   0.260970  0.2283   0.81980    
affile6      -0.032406   0.264297 -0.1226   0.90261    
afdec6       -0.565245   0.285835 -1.9775   0.05020 .  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

# resultados com Erros-padrões HAC - robusto a heterocedasticidade e
# autocorrelação
coeftest(mod1, vcovHAC)

t test of coefficients:

              Estimate Std. Error t value  Pr(>|t|)    
(Intercept) -17.802768  26.497047 -0.6719   0.50291    
log(chempi)   3.117194   0.654191  4.7650 5.188e-06 ***
log(gas)      0.196341   1.196616  0.1641   0.86994    
log(rtwex)    0.983016   0.450206  2.1835   0.03088 *  
befile6       0.059574   0.153259  0.3887   0.69815    
affile6      -0.032406   0.233788 -0.1386   0.88998    
afdec6       -0.565245   0.249559 -2.2650   0.02525 *  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

summary(mod1, vcov. = vcovHAC)

Time series regression with "ts" data:
Start = 1978(2), End = 1988(12)

Call:
dynlm(formula = log(chnimp) ~ log(chempi) + log(gas) + log(rtwex) + 
    befile6 + affile6 + afdec6, data = tsdata)

Residuals:
     Min       1Q   Median       3Q      Max 
-2.03356 -0.39080  0.03048  0.40248  1.51720 

Coefficients:
             Estimate Std. Error t value Pr(>|t|)    
(Intercept) -17.80277   26.49705  -0.672   0.5029    
log(chempi)   3.11719    0.65419   4.765 5.19e-06 ***
log(gas)      0.19634    1.19662   0.164   0.8699    
log(rtwex)    0.98302    0.45021   2.183   0.0309 *  
befile6       0.05957    0.15326   0.389   0.6982    
affile6      -0.03241    0.23379  -0.139   0.8900    
afdec6       -0.56524    0.24956  -2.265   0.0253 *  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 0.5974 on 124 degrees of freedom
Multiple R-squared:  0.3049,    Adjusted R-squared:  0.2712 
F-statistic:   5.8 on 6 and 124 DF,  p-value: 2.366e-05