1. Chuẩn bị dữ liệu

suppressMessages(library(vars))
library(thongke.club)
dulieu <-price_tom 

matranA <- diag(4)
diag(matranA) <- NA
matranA[2, 1] <- NA
matranA[3, 1] <- NA
matranA[3, 2] <- NA
matranA[4, 1] <- NA
matranA[4, 2] <- NA
matranA[4, 3] <- NA
print(matranA)
##      [,1] [,2] [,3] [,4]
## [1,]   NA    0    0    0
## [2,]   NA   NA    0    0
## [3,]   NA   NA   NA    0
## [4,]   NA   NA   NA   NA
matranB <- diag(4)
diag(matranB) <- NA
print(matranB)
##      [,1] [,2] [,3] [,4]
## [1,]   NA    0    0    0
## [2,]    0   NA    0    0
## [3,]    0    0   NA    0
## [4,]    0    0    0   NA

2 Ươcs lượng VAR

chonvar <- VAR(dulieu, p = 1, type = "both", season = 4)
summary(chonvar)
## 
## VAR Estimation Results:
## ========================= 
## Endogenous variables: price_ao, price_bs, price_bl, price_xk 
## Deterministic variables: both 
## Sample size: 119 
## Log Likelihood: -1948.211 
## Roots of the characteristic polynomial:
## 0.2234 0.08241 0.08241 0.002315
## Call:
## VAR(y = dulieu, p = 1, type = "both", season = 4L)
## 
## 
## Estimation results for equation price_ao: 
## ========================================= 
## price_ao = price_ao.l1 + price_bs.l1 + price_bl.l1 + price_xk.l1 + const + trend + sd1 + sd2 + sd3 
## 
##             Estimate Std. Error t value Pr(>|t|)    
## price_ao.l1 -0.03265    0.12815  -0.255    0.799    
## price_bs.l1  0.14200    0.12826   1.107    0.271    
## price_bl.l1  0.14896    0.12422   1.199    0.233    
## price_xk.l1 -0.04288    0.10253  -0.418    0.677    
## const       95.81012   21.21310   4.517 1.59e-05 ***
## trend       -0.04377    0.04810  -0.910    0.365    
## sd1         -3.90670    4.71354  -0.829    0.409    
## sd2         -3.22234    4.71152  -0.684    0.495    
## sd3         -5.04026    4.65980  -1.082    0.282    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## 
## Residual standard error: 17.89 on 110 degrees of freedom
## Multiple R-Squared: 0.04609, Adjusted R-squared: -0.02328 
## F-statistic: 0.6644 on 8 and 110 DF,  p-value: 0.7217 
## 
## 
## Estimation results for equation price_bs: 
## ========================================= 
## price_bs = price_ao.l1 + price_bs.l1 + price_bl.l1 + price_xk.l1 + const + trend + sd1 + sd2 + sd3 
## 
##               Estimate Std. Error t value Pr(>|t|)    
## price_ao.l1   0.047709   0.100151   0.476    0.635    
## price_bs.l1  -0.217667   0.100237  -2.172    0.032 *  
## price_bl.l1   0.053554   0.097081   0.552    0.582    
## price_xk.l1  -0.070068   0.080132  -0.874    0.384    
## const       145.547896  16.578441   8.779 2.41e-14 ***
## trend        -0.004617   0.037594  -0.123    0.902    
## sd1           1.474939   3.683724   0.400    0.690    
## sd2           0.107636   3.682139   0.029    0.977    
## sd3          -3.031705   3.641718  -0.832    0.407    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## 
## Residual standard error: 13.98 on 110 degrees of freedom
## Multiple R-Squared: 0.06256, Adjusted R-squared: -0.005622 
## F-statistic: 0.9175 on 8 and 110 DF,  p-value: 0.5049 
## 
## 
## Estimation results for equation price_bl: 
## ========================================= 
## price_bl = price_ao.l1 + price_bs.l1 + price_bl.l1 + price_xk.l1 + const + trend + sd1 + sd2 + sd3 
## 
##               Estimate Std. Error t value Pr(>|t|)    
## price_ao.l1  -0.032390   0.123174  -0.263    0.793    
## price_bs.l1   0.095157   0.123279   0.772    0.442    
## price_bl.l1   0.009707   0.119398   0.081    0.935    
## price_xk.l1   0.017139   0.098552   0.174    0.862    
## const       131.379271  20.389417   6.444 3.18e-09 ***
## trend        -0.005355   0.046236  -0.116    0.908    
## sd1           4.900592   4.530521   1.082    0.282    
## sd2           5.907655   4.528572   1.305    0.195    
## sd3           0.946433   4.478860   0.211    0.833    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## 
## Residual standard error: 17.19 on 110 degrees of freedom
## Multiple R-Squared: 0.0303,  Adjusted R-squared: -0.04022 
## F-statistic: 0.4296 on 8 and 110 DF,  p-value: 0.9011 
## 
## 
## Estimation results for equation price_xk: 
## ========================================= 
## price_xk = price_ao.l1 + price_bs.l1 + price_bl.l1 + price_xk.l1 + const + trend + sd1 + sd2 + sd3 
## 
##              Estimate Std. Error t value Pr(>|t|)    
## price_ao.l1   0.03752    0.13258   0.283    0.778    
## price_bs.l1   0.16149    0.13269   1.217    0.226    
## price_bl.l1  -0.02474    0.12851  -0.193    0.848    
## price_xk.l1   0.02100    0.10608   0.198    0.843    
## const       146.43220   21.94604   6.672 1.06e-09 ***
## trend        -0.03802    0.04977  -0.764    0.447    
## sd1          -4.25048    4.87640  -0.872    0.385    
## sd2           0.07519    4.87430   0.015    0.988    
## sd3           0.82891    4.82080   0.172    0.864    
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## 
## Residual standard error: 18.51 on 110 degrees of freedom
## Multiple R-Squared: 0.03662, Adjusted R-squared: -0.03344 
## F-statistic: 0.5227 on 8 and 110 DF,  p-value: 0.8373 
## 
## 
## 
## Covariance matrix of residuals:
##          price_ao price_bs price_bl price_xk
## price_ao    319.9   101.76   184.49   146.34
## price_bs    101.8   195.42    73.95    63.51
## price_bl    184.5    73.95   295.58    99.41
## price_xk    146.3    63.51    99.41   342.44
## 
## Correlation matrix of residuals:
##          price_ao price_bs price_bl price_xk
## price_ao   1.0000   0.4070   0.5999   0.4421
## price_bs   0.4070   1.0000   0.3077   0.2455
## price_bl   0.5999   0.3077   1.0000   0.3125
## price_xk   0.4421   0.2455   0.3125   1.0000

Chạy SVAR

Không dùng hàm summary

Phân tích irf,fevd, & test (shortsvar + longsvar)

# Chạy short-run
shortsvar <- SVAR(chonvar, Amat = matranA, Bmat = matranB, max.iter = 10000, hessian = TRUE)
## Warning in SVAR(chonvar, Amat = matranA, Bmat = matranB, max.iter = 10000, : The
## AB-model is just identified. No test possible.
shortsvar
## 
## SVAR Estimation Results:
## ======================== 
## 
## 
## Estimated A matrix:
##          price_ao price_bs price_bl price_xk
## price_ao   1.0000  0.00000  0.00000        0
## price_bs  -0.3181  1.00000  0.00000        0
## price_bl  -0.5469 -0.09363  1.00000        0
## price_xk  -0.3851 -0.09737 -0.07155        1
## 
## Estimated B matrix:
##          price_ao price_bs price_bl price_xk
## price_ao    17.89     0.00      0.0     0.00
## price_bs     0.00    12.77      0.0     0.00
## price_bl     0.00     0.00     13.7     0.00
## price_xk     0.00     0.00      0.0    16.52
#Chay long-run 
longsvar <- BQ(chonvar)
summary(longsvar)
## 
## SVAR Estimation Results:
## ======================== 
## 
## Call:
## BQ(x = chonvar)
## 
## Type: Blanchard-Quah 
## Sample size: 119 
## Log Likelihood: -1966.928 
## 
## Estimated contemporaneous impact matrix:
##          price_ao price_bs price_bl price_xk
## price_ao   17.707  -1.5618  -1.8498   0.7361
## price_bs    6.708  12.1900  -0.6172   1.2028
## price_bl   11.763   0.2568  12.5323  -0.2942
## price_xk    7.636  -0.5904   1.1713  16.8051
## 
## Estimated identified long run impact matrix:
##          price_ao price_bs price_bl price_xk
## price_ao   19.356    0.000   0.0000     0.00
## price_bs    6.262   10.007   0.0000     0.00
## price_bl   12.008    1.239  12.6703     0.00
## price_xk    9.271    1.016   0.8762    17.17
## 
## Covariance matrix of reduced form residuals (*100):
##          price_ao price_bs price_bl price_xk
## price_ao    31995    10176    18449    14634
## price_bs    10176    19542     7395     6351
## price_bl    18449     7395    29558     9941
## price_xk    14634     6351     9941    34244

Chạy SVEC

suppressMessages(library(urca)) # lấy hàm ca.jo
vecm <- ca.jo(dulieu, type = "trace", ecdet = "trend", K = 3, spec = "transitory")

shortrun <- matrix(NA, nrow = 4, ncol = 4) # 4 số biến trong mô hình
shortrun[4, 2] <- 0
shortrun
##      [,1] [,2] [,3] [,4]
## [1,]   NA   NA   NA   NA
## [2,]   NA   NA   NA   NA
## [3,]   NA   NA   NA   NA
## [4,]   NA    0   NA   NA
longrun <- matrix(NA, nrow = 4, ncol = 4)
longrun[1, 2:4] <- 0
longrun[2:4, 4] <- 0
longrun
##      [,1] [,2] [,3] [,4]
## [1,]   NA    0    0    0
## [2,]   NA   NA   NA    0
## [3,]   NA   NA   NA    0
## [4,]   NA   NA   NA    0
SVEC(vecm, LR = longrun, SR = shortrun, r = 1, lrtest = TRUE, boot = FALSE)
## Warning in SVEC(vecm, LR = longrun, SR = shortrun, r = 1, lrtest = TRUE, : The
## SVEC is just identified. No test possible.
## 
## SVEC Estimation Results:
## ======================== 
## 
## 
## Estimated contemporaneous impact matrix:
##          price_ao price_bs price_bl price_xk
## price_ao   19.306   0.1526   0.4445  -4.2140
## price_bs    8.862   6.9261  -1.5892   8.1572
## price_bl   14.029 -11.9952   0.6142   5.4248
## price_xk    9.171   0.0000  18.0323   0.2427
## 
## Estimated long run impact matrix:
##          price_ao price_bs price_bl price_xk
## price_ao   10.544   0.0000   0.0000        0
## price_bs    3.629   3.4383  -0.7269        0
## price_bl    7.233  -5.1947   0.9056        0
## price_xk    4.841  -0.3375   9.1126        0