tải các gói dữ liệu

library(cointReg)
## Warning: package 'cointReg' was built under R version 4.0.5
## cointReg (v0.2.0): Parameter Estimation and Inference in a Cointegrating Regression.
library(readxl)

Kết nối dữ liệu

setwd("d:/DATA2021/ThuySan.VECM/data5.5.21")
solieu <-read_excel("data.ca.5.5.21.xlsx")

congthuc1 <-LnGDP ~ LnFEX + LnLAB + REER + LnOPEN + LnFDI
congthuc2 <-LnFGDP ~ LnFEX + LnLAB + REER + LnOPEN + LnFDI
lm(data=solieu,congthuc1) -> hoiquy1
lm(data=solieu,congthuc2) -> hoiquy2
summary(hoiquy1)
## 
## Call:
## lm(formula = congthuc1, data = solieu)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.210940 -0.058072  0.003733  0.039472  0.269841 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 23.35330    4.89950   4.766 9.14e-06 ***
## LnFEX        0.40039    0.06716   5.962 7.86e-08 ***
## LnLAB       -2.96920    0.63490  -4.677 1.28e-05 ***
## REER        -1.00061    0.53801  -1.860  0.06688 .  
## LnOPEN       0.53708    0.08649   6.210 2.81e-08 ***
## LnFDI        0.13661    0.04067   3.359  0.00124 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.09556 on 74 degrees of freedom
## Multiple R-squared:  0.9455, Adjusted R-squared:  0.9418 
## F-statistic: 256.6 on 5 and 74 DF,  p-value: < 2.2e-16
summary(hoiquy2)
## 
## Call:
## lm(formula = congthuc2, data = solieu)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -0.8529 -0.1832  0.0617  0.2316  1.0591 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  -0.3754    21.8932  -0.017 0.986364    
## LnFEX        -0.5312     0.3001  -1.770 0.080796 .  
## LnLAB        -1.4912     2.8370  -0.526 0.600712    
## REER         -2.4703     2.4041  -1.028 0.307505    
## LnOPEN        1.3723     0.3865   3.551 0.000672 ***
## LnFDI         0.4846     0.1817   2.667 0.009406 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.427 on 74 degrees of freedom
## Multiple R-squared:   0.87,  Adjusted R-squared:  0.8612 
## F-statistic: 99.03 on 5 and 74 DF,  p-value: < 2.2e-16

Tạo dữ liệu để chạy Cointreg

attach(solieu)
doclap <-cbind(LnFEX,  LnLAB,  REER,  LnOPEN,  LnFDI)
phuthuoc1 <-cbind(LnGDP)
phuthuoc2 <-cbind(LnFGDP)

Chạy FMOLS

cointRegFM(x=doclap, y=phuthuoc1)
## 
## ### FM-OLS model ###
## 
## Model:       phuthuoc1 ~ doclap
## 
## Parameters:  Kernel = "ba"  //  Bandwidth = 5.278472 ("Andrews")
## 
## Coefficients:
##          Estimate    Std.Err t value Pr(|t|>0)    
## LnFEX   0.4291942  0.0904645  4.7443 9.769e-06 ***
## LnLAB   0.1567290  0.1474494  1.0629  0.291222    
## REER   -2.0591037  0.7135396 -2.8858  0.005096 ** 
## LnOPEN  0.1919875  0.0410366  4.6784 1.253e-05 ***
## LnFDI   0.0054448  0.0530039  0.1027  0.918455    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
cointRegFM(x=doclap, y=phuthuoc2)
## 
## ### FM-OLS model ###
## 
## Model:       phuthuoc2 ~ doclap
## 
## Parameters:  Kernel = "ba"  //  Bandwidth = 10.90026 ("Andrews")
## 
## Coefficients:
##        Estimate  Std.Err t value Pr(|t|>0)    
## LnFEX  -1.11228  0.59090 -1.8824   0.06367 .  
## LnLAB  -0.48315  0.96311 -0.5017   0.61738    
## REER   -5.02120  4.66071 -1.0773   0.28478    
## LnOPEN  1.46460  0.26804  5.4640 5.823e-07 ***
## LnFDI   0.58440  0.34621  1.6880   0.09557 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Chạy kiểm định

library(tseries) # Có hàm lag L()
## Registered S3 method overwritten by 'quantmod':
##   method            from
##   as.zoo.data.frame zoo
mohinh1 <-cbind(LnGDP, LnFEX,  LnLAB,  REER,  LnOPEN,  LnFDI)
mohinh2 <-cbind(LnFGDP, LnFEX,  LnLAB,  REER,  LnOPEN,  LnFDI)

#Null hypothesis: Series are not cointegrated       
po.test(mohinh1)
## Warning in po.test(mohinh1): p-value smaller than printed p-value
## 
##  Phillips-Ouliaris Cointegration Test
## 
## data:  mohinh1
## Phillips-Ouliaris demeaned = -61.168, Truncation lag parameter = 0,
## p-value = 0.01
po.test(mohinh2)
## Warning in po.test(mohinh2): p-value greater than printed p-value
## 
##  Phillips-Ouliaris Cointegration Test
## 
## data:  mohinh2
## Phillips-Ouliaris demeaned = -20.277, Truncation lag parameter = 0,
## p-value = 0.15
library(aTSA)
## Warning: package 'aTSA' was built under R version 4.0.3
## 
## Attaching package: 'aTSA'
## The following objects are masked from 'package:tseries':
## 
##     adf.test, kpss.test, pp.test
## The following object is masked from 'package:graphics':
## 
##     identify
#Null hypothesis: Series are not cointegrated       
coint.test(phuthuoc1,doclap, d=2, nlag=3)
## Response: diff(phuthuoc1,2) 
## Input: diff(doclap,2) 
## Number of inputs: 5 
## Model: y ~ X - 1 
## ------------------------------- 
## Engle-Granger Cointegration Test 
## alternative: cointegrated 
## 
## Type 1: no trend 
##     lag      EG p.value 
##    3.00  -10.14    0.01 
## ----- 
##  Type 2: linear trend 
##     lag      EG p.value 
##    3.00    6.52    0.10 
## ----- 
##  Type 3: quadratic trend 
##     lag      EG p.value 
##    3.00    4.97    0.10 
## ----------- 
## Note: p.value = 0.01 means p.value <= 0.01 
##     : p.value = 0.10 means p.value >= 0.10
coint.test(phuthuoc2, doclap,d=1, nlag=3)
## Response: diff(phuthuoc2,1) 
## Input: diff(doclap,1) 
## Number of inputs: 5 
## Model: y ~ X - 1 
## ------------------------------- 
## Engle-Granger Cointegration Test 
## alternative: cointegrated 
## 
## Type 1: no trend 
##     lag      EG p.value 
##    3.00   -6.88    0.01 
## ----- 
##  Type 2: linear trend 
##     lag      EG p.value 
##    3.00    1.09    0.10 
## ----- 
##  Type 3: quadratic trend 
##     lag      EG p.value 
##    3.00    0.39    0.10 
## ----------- 
## Note: p.value = 0.01 means p.value <= 0.01 
##     : p.value = 0.10 means p.value >= 0.10

Chạy ECM

library(dynlm)
## Loading required package: zoo
## 
## Attaching package: 'zoo'
## The following objects are masked from 'package:base':
## 
##     as.Date, as.Date.numeric
library(collapse) #lấy hàm sai phân D()
## Warning: package 'collapse' was built under R version 4.0.5
## collapse 1.5.3, see ?`collapse-package` or ?`collapse-documentation`
## Note: stats::D  ->  D.expression, D.call, D.name
## 
## Attaching package: 'collapse'
## The following object is masked from 'package:zoo':
## 
##     is.regular
## The following object is masked from 'package:stats':
## 
##     D
solieu <- ts(solieu,start = c(2000,1),frequency = 4)
uocluong1 <-dynlm(congthuc1)
summary(uocluong1)
## 
## Time series regression with "numeric" data:
## Start = 1, End = 80
## 
## Call:
## dynlm(formula = congthuc1)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.210940 -0.058072  0.003733  0.039472  0.269841 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept) 23.35330    4.89950   4.766 9.14e-06 ***
## LnFEX        0.40039    0.06716   5.962 7.86e-08 ***
## LnLAB       -2.96920    0.63490  -4.677 1.28e-05 ***
## REER        -1.00061    0.53801  -1.860  0.06688 .  
## LnOPEN       0.53708    0.08649   6.210 2.81e-08 ***
## LnFDI        0.13661    0.04067   3.359  0.00124 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.09556 on 74 degrees of freedom
## Multiple R-squared:  0.9455, Adjusted R-squared:  0.9418 
## F-statistic: 256.6 on 5 and 74 DF,  p-value: < 2.2e-16
phandu1 <- residuals(uocluong1)
uocluong2 <-dynlm(congthuc2)
summary(uocluong2)
## 
## Time series regression with "numeric" data:
## Start = 1, End = 80
## 
## Call:
## dynlm(formula = congthuc2)
## 
## Residuals:
##     Min      1Q  Median      3Q     Max 
## -0.8529 -0.1832  0.0617  0.2316  1.0591 
## 
## Coefficients:
##             Estimate Std. Error t value Pr(>|t|)    
## (Intercept)  -0.3754    21.8932  -0.017 0.986364    
## LnFEX        -0.5312     0.3001  -1.770 0.080796 .  
## LnLAB        -1.4912     2.8370  -0.526 0.600712    
## REER         -2.4703     2.4041  -1.028 0.307505    
## LnOPEN        1.3723     0.3865   3.551 0.000672 ***
## LnFDI         0.4846     0.1817   2.667 0.009406 ** 
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.427 on 74 degrees of freedom
## Multiple R-squared:   0.87,  Adjusted R-squared:  0.8612 
## F-statistic: 99.03 on 5 and 74 DF,  p-value: < 2.2e-16
phandu2 <- residuals(uocluong2)
dulieu <-cbind(solieu,phandu1,phandu2)
head(dulieu)
##      solieu.STT solieu.QUY solieu.LnGDP solieu.LnFGDP solieu.LnFEX solieu.LnLAB
## [1,]          1          1     12.32272      7.204922     19.51234     10.55482
## [2,]          2          2     12.57529      7.457496     19.76492     10.56007
## [3,]          3          3     12.50762      7.389821     19.69724     10.56939
## [4,]          4          4     12.70953      7.591730     19.89915     10.55959
## [5,]          5          5     12.39648      7.320969     19.72526     10.58206
## [6,]          6          6     12.63004      7.554526     19.95881     10.58733
##      solieu.REER solieu.LnOPEN solieu.LnFDI      phandu1    phandu2
## [1,]   0.9962832      22.62871     10.34393 -0.073536826 0.08117174
## [2,]   0.9997200      22.73124     10.59650  0.007360259 0.22113928
## [3,]   1.0071360      22.71304     10.52882  0.020915954 0.20751552
## [4,]   1.0072360      22.87906     10.73074 -0.003762970 0.17665428
## [5,]   0.9931406      22.57954     10.45481  0.003977180 0.35692495
## [6,]   0.9989467      22.81070     10.68836  0.009430951 0.30637390
EC1 <- D(LnGDP)~L(D(LnFEX),1)+L(D(LnLAB),1)+ L(D(REER),1) + L(D(LnOPEN),1)+ L(D(LnFDI),1) + L(phandu1,1)

EC2 <-D(LnFGDP)~L(D(LnFEX),1)+L(D(LnLAB),1)+ L(D(REER),1) + L(D(LnOPEN),1)+ L(D(LnFDI),1) + L(phandu2,1)
ECM1 <-dynlm(EC1, data=dulieu[2:80])
ECM2 <-dynlm(EC2, data=dulieu[2:80])
summary(ECM1)
## 
## Time series regression with "numeric" data:
## Start = 1, End = 79
## 
## Call:
## dynlm(formula = EC1, data = dulieu[2:80])
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.17580 -0.03513 -0.01002  0.03479  0.19111 
## 
## Coefficients:
##                  Estimate Std. Error t value Pr(>|t|)    
## (Intercept)      0.004337   0.009717   0.446  0.65674    
## L(D(LnFEX), 1)   0.426100   0.070154   6.074 5.35e-08 ***
## L(D(LnLAB), 1)  -1.992476   0.694725  -2.868  0.00542 ** 
## L(D(REER), 1)   -0.166709   1.094446  -0.152  0.87936    
## L(D(LnOPEN), 1)  0.287304   0.110339   2.604  0.01119 *  
## L(D(LnFDI), 1)   0.378410   0.057574   6.573 6.72e-09 ***
## L(phandu1, 1)    0.541213   0.100791   5.370 9.21e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.07045 on 72 degrees of freedom
##   (1 observation deleted due to missingness)
## Multiple R-squared:  0.9341, Adjusted R-squared:  0.9286 
## F-statistic:   170 on 6 and 72 DF,  p-value: < 2.2e-16
summary(ECM2)
## 
## Time series regression with "numeric" data:
## Start = 1, End = 79
## 
## Call:
## dynlm(formula = EC2, data = dulieu[2:80])
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.33482 -0.08772 -0.03993  0.04108  1.25331 
## 
## Coefficients:
##                 Estimate Std. Error t value Pr(>|t|)    
## (Intercept)     0.009182   0.027852   0.330 0.742611    
## L(D(LnFEX), 1)  0.709869   0.200586   3.539 0.000708 ***
## L(D(LnLAB), 1)  1.869965   1.998832   0.936 0.352643    
## L(D(REER), 1)   4.451317   3.030098   1.469 0.146179    
## L(D(LnOPEN), 1) 0.021020   0.315028   0.067 0.946987    
## L(D(LnFDI), 1)  0.411386   0.147801   2.783 0.006866 ** 
## L(phandu2, 1)   0.104218   0.058341   1.786 0.078248 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2031 on 72 degrees of freedom
##   (1 observation deleted due to missingness)
## Multiple R-squared:  0.7041, Adjusted R-squared:  0.6795 
## F-statistic: 28.56 on 6 and 72 DF,  p-value: < 2.2e-16

Chạy ECM bằng gói aTSA

library(aTSA)
ECM1 <-ecm(phuthuoc1,doclap)
## 
## Call:
## lm(formula = dy ~ dX + ECM - 1)
## 
## Residuals:
##       Min        1Q    Median        3Q       Max 
## -0.172016 -0.024001 -0.004248  0.029151  0.219716 
## 
## Coefficients:
##          Estimate Std. Error t value Pr(>|t|)    
## dXLnFEX   0.43767    0.07065   6.195 3.11e-08 ***
## dXLnLAB  -1.05283    0.63447  -1.659   0.1013    
## dXREER    0.14292    1.06249   0.135   0.8934    
## dXLnOPEN  0.19695    0.09711   2.028   0.0462 *  
## dXLnFDI   0.38412    0.05672   6.772 2.74e-09 ***
## ECM      -0.48433    0.09056   5.348 9.75e-07 ***
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.07055 on 73 degrees of freedom
## Multiple R-squared:  0.9334, Adjusted R-squared:  0.9279 
## F-statistic: 170.5 on 6 and 73 DF,  p-value: < 2.2e-16
ECM2 <-ecm(phuthuoc2,doclap)
## 
## Call:
## lm(formula = dy ~ dX + ECM - 1)
## 
## Residuals:
##      Min       1Q   Median       3Q      Max 
## -0.32746 -0.08091 -0.03058  0.04549  1.25992 
## 
## Coefficients:
##          Estimate Std. Error t value Pr(>|t|)    
## dXLnFEX   0.70178    0.19791   3.546 0.000687 ***
## dXLnLAB   2.17686    1.75347   1.241 0.218409    
## dXREER    4.30587    2.97986   1.445 0.152740    
## dXLnOPEN  0.06761    0.28043   0.241 0.810165    
## dXLnFDI   0.40753    0.14643   2.783 0.006849 ** 
## ECM      -0.10441    0.05798   1.801 0.075886 .  
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.2019 on 73 degrees of freedom
## Multiple R-squared:  0.7076, Adjusted R-squared:  0.6836 
## F-statistic: 29.44 on 6 and 73 DF,  p-value: < 2.2e-16