移動平均自己回帰モデル/ARMA(Auto Regressive Moving Average Model)

Equation
ARMA(p,q)
y_{t}=c+\sum^{p}_{i=1}a_{i}y_{t-i}+\varepsilon_{t}+\sum^{q}_{j=1}b_{j}\varepsilon_{t-j}

R
サンプル:2009年2月から2014年2月までのゴールド価格(月次平均、米ドル/トロイオンス、ロンドン市場、AM10:30)
※Jarque–Bera test(ジャック-ベラ検定)の帰無仮説:「分布は正規分布である」
※Ljung-Box test(リュング・ボックス検定)の帰無仮説:「自己相関無し」

goldlondongoldlondon_1stdiffgoldarimaqqplotgoldarima

> library(tseries)
> ts.gold <- ts(dataset$gold,start=c(2009,2),frequency=12)
> ts.gold
          Jan      Feb      Mar      Apr      May      Jun      Jul      Aug      Sep      Oct      Nov      Dec
2009           943.163  924.273  890.200  928.645  945.670  934.228  949.375  996.591 1043.159 1127.036 1134.724
2010 1117.963 1095.413 1113.337 1148.688 1205.434 1232.920 1192.966 1215.810 1270.977 1342.024 1369.886 1390.553
2011 1356.400 1372.725 1424.011 1473.806 1510.438 1528.659 1572.810 1755.807 1771.852 1665.214 1738.977 1652.306
2012 1656.119 1742.619 1673.773 1650.066 1585.505 1596.697 1593.909 1626.034 1744.450 1747.011 1721.136 1688.529
2013 1670.955 1627.588 1592.863 1485.083 1413.500 1342.363 1286.724 1347.095 1348.798 1316.185 1275.821 1225.403
2014 1244.795 1300.975                                                                                          
> adf.test(ts.gold)

        Augmented Dickey-Fuller Test

data:  ts.gold
Dickey-Fuller = -0.7502, Lag order = 3, p-value = 0.9609
alternative hypothesis: stationary

> ts.diff.gold <- diff(ts.gold)
> ts.diff.gold
          Jan      Feb      Mar      Apr      May      Jun      Jul      Aug      Sep      Oct      Nov      Dec
2009                    -18.890  -34.073   38.445   17.025  -11.442   15.147   47.216   46.568   83.877    7.688
2010  -16.761  -22.550   17.924   35.351   56.746   27.486  -39.954   22.844   55.167   71.047   27.862   20.667
2011  -34.153   16.325   51.286   49.795   36.632   18.221   44.151  182.997   16.045 -106.638   73.763  -86.671
2012    3.813   86.500  -68.846  -23.707  -64.561   11.192   -2.788   32.125  118.416    2.561  -25.875  -32.607
2013  -17.574  -43.367  -34.725 -107.780  -71.583  -71.137  -55.639   60.371    1.703  -32.613  -40.364  -50.418
2014   19.392   56.180                                                                                          
> adf.test(ts.diff.gold)

        Augmented Dickey-Fuller Test

data:  ts.diff.gold
Dickey-Fuller = -4.5173, Lag order = 3, p-value = 0.01
alternative hypothesis: stationary

 警告メッセージ: 
In adf.test(ts.diff.gold) : p-value smaller than printed p-value

> for(i in 0:3)
+ {
+ for(j in 0:3)
+ {
+ tmp <- arima(ts.diff.gold,order=c(i,0,j))
+ cat ("order=(",i,",0,",j,"),AIC=",tmp$aic,"\n")
+ }
+ }

order=( 0 ,0, 0 ),AIC= 652.4788 
order=( 0 ,0, 1 ),AIC= 650.8186 
order=( 0 ,0, 2 ),AIC= 651.7375 
order=( 0 ,0, 3 ),AIC= 648.4721 
order=( 1 ,0, 0 ),AIC= 651.1776 
order=( 1 ,0, 1 ),AIC= 645.007 
order=( 1 ,0, 2 ),AIC= 646.9539 
order=( 1 ,0, 3 ),AIC= 647.6591 
order=( 2 ,0, 0 ),AIC= 653.1097 
order=( 2 ,0, 1 ),AIC= 646.9268 
order=( 2 ,0, 2 ),AIC= 648.0994 
order=( 2 ,0, 3 ),AIC= 649.9996 
order=( 3 ,0, 0 ),AIC= 654.2283 
order=( 3 ,0, 1 ),AIC= 647.6391 
order=( 3 ,0, 2 ),AIC= 649.5061 
order=( 3 ,0, 3 ),AIC= 651.6075 

> ts.diff.gold.arima <- arima(ts.diff.gold,order=c(1,0,1))
> plot(ts.gold,type="l",main="Gold.London.Monthly",ylab="$/ozt")
> plot(ts.diff.gold,type="l",main="1stDiff.Gold.London.Monthly",ylab="$/ozt",ylim=c(-200,200))
> par(new=t)
> plot(ts.diff.gold-ts.diff.gold.arima$resid,main="",ylab="",xlab="",col=2,ylim=c(-200,200))
> ts.diff.gold.arima

Call:
arima(x = ts.diff.gold, order = c(1, 0, 1))

Coefficients:
          ar1     ma1  intercept
      -0.5898  1.0000     6.1301
s.e.   0.1077  0.0511     7.7406

sigma^2 estimated as 2281:  log likelihood = -318.5,  aic = 645.01
> ts.diff.gold.arima$resid
              Jan          Feb          Mar          Apr          May
2009                            -22.3073501  -31.9481555   34.8560834
2010   -0.9830779  -39.9124186   32.5050271    4.5179629   61.7571541
2011  -45.7138365   30.7245513   20.6338465   49.1806340    7.7763153
2012   42.3195451   36.7620535  -63.0856295  -11.6084397  -75.9194346
2013  -21.0142401  -42.2582099  -27.9322123 -109.3116071  -36.2710964
2014   24.3552981   33.4463051                                       
              Jun          Jul          Aug          Sep          Oct
2009   -2.1184722   -8.6350393    6.4340433   38.4070614   26.8641898
2010   -8.6628475  -24.4390527   13.0703113   45.1075341   48.6763790
2011   22.0771921   23.0652197  173.9551078  -56.3107538  -50.7013578
2012   37.9565032  -42.9671746   62.5279287   65.0702941   -1.6964590
2013  -86.3806565  -21.5459918   38.8250985  -10.8304540  -30.3595850
2014                                                                 
              Nov          Dec
2009   72.8829916  -21.7529612
2010   12.1231172   15.1756339
2011   50.4135673 -101.2723873
2012  -32.0981872  -25.5856001
2013  -38.9166717  -45.0070796
2014                          
> qqnorm(ts.diff.gold.arima$resid)
> qqline(ts.diff.gold.arima$resid,col=2)
> jarque.bera.test(ts.diff.gold.arima$resid)

        Jarque Bera Test

data:  ts.diff.gold.arima$resid
X-squared = 7.9733, df = 2, p-value = 0.01856

> Box.test(ts.diff.gold.arima$resid,lag=30,type="Ljung")

        Box-Ljung test

data:  ts.diff.gold.arima$resid
X-squared = 22.9319, df = 30, p-value = 0.818

> 

参考文献
大崎秀一,吉川大介(2013).『ファイナンスのためのRプログラミング』.共立出版.179pp.

アプリケーション
R Core Team (2013). R: A language and environment for statistical computing.
R Foundation for Statistical Computing, Vienna, Austria.
URL http://www.R-project.org/.