Content

\[ \Delta y_{t}= \alpha + \beta t + \gamma y_{t-1} + \delta_{1} \Delta y_{t-1}+ \delta_{2} \Delta y_{t-2} + \epsilon_{t} \]

head(USDJPY)
cat("------------------------\n")
tail(USDJPY)
cat("------------------------\n")
ts_data <- USDJPY$`USD/JPY`
ts_data
cat("------------------------\n")
diff_ts_data <- diff(ts_data)
diff_ts_data
cat("------------------------\n")
length(diff_ts_data)
        Date  USD/JPY
1 2017-05-01 112.2436
2 2017-06-01 110.9141
3 2017-07-01 112.4170
4 2017-08-01 109.8270
5 2017-09-01 110.7760
6 2017-10-01 112.9148
------------------------
         Date  USD/JPY
19 2018-11-01 113.3380
20 2018-12-01 112.1994
21 2019-01-01 108.9605
22 2019-02-01 110.4400
23 2019-03-01 111.1443
24 2019-04-01 111.6073
------------------------
 [1] 112.2436 110.9141 112.4170 109.8270 110.7760 112.9148 112.8190 112.9405 110.8710 107.9700 106.0468 107.6562 109.6882 110.0638 111.5210 110.9965 112.0974 112.7218 113.3380 112.1994 108.9605 110.4400 111.1443 111.6073
------------------------
 [1] -1.3295  1.5029 -2.5900  0.9490  2.1388 -0.0958  0.1215 -2.0695 -2.9010 -1.9232  1.6094  2.0320  0.3756  1.4572 -0.5245  1.1009  0.6244  0.6162 -1.1386 -3.2389  1.4795  0.7043  0.4630
------------------------
[1] 23
# ここではラグを2次まで取ります。
lag <- 2
# 自己回帰でお世話になるembed{stats}
M_diff <- embed(x = diff_ts_data, dimension = lag + 1)
M_diff
         [,1]    [,2]    [,3]
 [1,] -2.5900  1.5029 -1.3295
 [2,]  0.9490 -2.5900  1.5029
 [3,]  2.1388  0.9490 -2.5900
 [4,] -0.0958  2.1388  0.9490
 [5,]  0.1215 -0.0958  2.1388
 [6,] -2.0695  0.1215 -0.0958
 [7,] -2.9010 -2.0695  0.1215
 [8,] -1.9232 -2.9010 -2.0695
 [9,]  1.6094 -1.9232 -2.9010
[10,]  2.0320  1.6094 -1.9232
[11,]  0.3756  2.0320  1.6094
[12,]  1.4572  0.3756  2.0320
[13,] -0.5245  1.4572  0.3756
[14,]  1.1009 -0.5245  1.4572
[15,]  0.6244  1.1009 -0.5245
[16,]  0.6162  0.6244  1.1009
[17,] -1.1386  0.6162  0.6244
[18,] -3.2389 -1.1386  0.6162
[19,]  1.4795 -3.2389 -1.1386
[20,]  0.7043  1.4795 -3.2389
[21,]  0.4630  0.7043  1.4795
# 原系列からラグ分を除くための添字の作成(トレンド項(linear trend)として流用)
trend_term <- c((lag + 1):(length(ts_data) - 1))
trend_term
 [1]  3  4  5  6  7  8  9 10 11 12 13 14 15 16 17 18 19 20 21 22 23
# 始めにトレンド項、ドリフト項(constant、定数項)ともに含む場合
result_lm <- summary(lm(M_diff[, 1] ~ ts_data[trend_term] + M_diff[, 2] + M_diff[, 3] + trend_term))
result_lm
cat("------------------------\n")
# ADF検定統計量 = 自己回帰係数(ここではts_dataの係数γ)のt値(但しt分布には従わない)
result_lm$coefficient[2, 1]/result_lm$coefficient[2, 2]
cat("------------------------\n")
# Rには拡張ディッキー=フラー検定関数としてCADFtest{CADFtest}があります
result_adftest <- CADFtest(model = ts_data, max.lag.y = lag, type = "trend")
result_adftest
cat("------------------------\n")
result_adftest %>>% (statistic)

Call:
lm(formula = M_diff[, 1] ~ ts_data[trend_term] + M_diff[, 2] + 
    M_diff[, 3] + trend_term)

Residuals:
    Min      1Q  Median      3Q     Max 
-2.2350 -0.8228  0.1731  1.0604  2.5038 

Coefficients:
                    Estimate Std. Error t value Pr(>|t|)  
(Intercept)         67.26657   23.97459   2.806   0.0127 *
ts_data[trend_term] -0.61067    0.21598  -2.828   0.0121 *
M_diff[, 2]          0.41852    0.21915   1.910   0.0743 .
M_diff[, 3]          0.20474    0.22639   0.904   0.3792  
trend_term           0.02984    0.05305   0.562   0.5816  
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 1.47 on 16 degrees of freedom
Multiple R-squared:  0.3658,    Adjusted R-squared:  0.2073 
F-statistic: 2.307 on 4 and 16 DF,  p-value: 0.1028

------------------------
[1] -2.827502
------------------------

    ADF test

data:  ts_data
ADF(2) = -2.8275, p-value = 0.2036
alternative hypothesis: true delta is less than 0
sample estimates:
     delta 
-0.6106707 

------------------------
   ADF(2) 
-2.827502 
# 続いてドリフト項を含みトレンド項は含まない場合(ドリフト付きランダムウォークモデル仮定)
result_lm <- summary(lm(M_diff[, 1] ~ ts_data[trend_term] + M_diff[, 2] + M_diff[, 3]))
result_lm
cat("------------------------\n")
result_lm$coefficient[2, 1]/result_lm$coefficient[2, 2]
cat("------------------------\n")
result_adftest <- CADFtest(model = ts_data, max.lag.y = lag, type = "drift")
result_adftest
cat("------------------------\n")
result_adftest %>>% (statistic)

Call:
lm(formula = M_diff[, 1] ~ ts_data[trend_term] + M_diff[, 2] + 
    M_diff[, 3])

Residuals:
    Min      1Q  Median      3Q     Max 
-2.0230 -0.7796  0.1067  1.0745  2.2787 

Coefficients:
                    Estimate Std. Error t value Pr(>|t|)   
(Intercept)          67.9926    23.4535   2.899  0.00998 **
ts_data[trend_term]  -0.6137     0.2115  -2.901  0.00993 **
M_diff[, 2]           0.4217     0.2146   1.965  0.06600 . 
M_diff[, 3]           0.2114     0.2215   0.955  0.35318   
---
Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1

Residual standard error: 1.44 on 17 degrees of freedom
Multiple R-squared:  0.3533,    Adjusted R-squared:  0.2391 
F-statistic: 3.095 on 3 and 17 DF,  p-value: 0.05471

------------------------
[1] -2.90144
------------------------

    ADF test

data:  ts_data
ADF(2) = -2.9014, p-value = 0.06203
alternative hypothesis: true delta is less than 0
sample estimates:
     delta 
-0.6137171 

------------------------
  ADF(2) 
-2.90144 
# 最後にトレンド項、ドリフト項ともに含まない場合(ランダムウォークモデル仮定)
result_lm <- summary(lm(M_diff[, 1] ~ ts_data[trend_term] + M_diff[, 2] + M_diff[, 3] - 1))
result_lm
cat("------------------------\n")
result_lm$coefficient[1, 1]/result_lm$coefficient[1, 2]
cat("------------------------\n")
result_adftest <- CADFtest(model = ts_data, max.lag.y = lag, type = "none")
result_adftest
cat("------------------------\n")
result_adftest %>>% (statistic)

Call:
lm(formula = M_diff[, 1] ~ ts_data[trend_term] + M_diff[, 2] + 
    M_diff[, 3] - 1)

Residuals:
    Min      1Q  Median      3Q     Max 
-2.9382 -1.0912  0.4493  1.5630  1.8780 

Coefficients:
                      Estimate Std. Error t value Pr(>|t|)
ts_data[trend_term] -0.0005631  0.0033713  -0.167    0.869
M_diff[, 2]          0.1452452  0.2284271   0.636    0.533
M_diff[, 3]         -0.1170646  0.2260907  -0.518    0.611

Residual standard error: 1.711 on 18 degrees of freedom
Multiple R-squared:  0.03409,   Adjusted R-squared:  -0.1269 
F-statistic: 0.2117 on 3 and 18 DF,  p-value: 0.8869

------------------------
[1] -0.1670245
------------------------

    ADF test

data:  ts_data
ADF(2) = -0.16702, p-value = 0.614
alternative hypothesis: true delta is less than 0
sample estimates:
        delta 
-0.0005630961 

------------------------
    ADF(2) 
-0.1670245