GRGと金融工学・統計解析

趣味で金融工学や統計解析を勉強しており，このブログでは自分が面白いと思ったことの紹介などをしていきます．色々な方との議論などできたらいいなぁと思います．

S&P500と日経225の時系列分析（コピュラ ＆ DCDモデル）

時系列 金融工学 統計

前回 –> コピュラとは（時系列分析・DCDモデルの前準備として） - GRGと金融工学・統計解析

はじめに

前回はコピュラについて簡単に触れましたが，今回はそのコピュラを導入した多変量時系列モデルであるDCDモデルについて見ていきたいと思います．DCDモデルはCopula-Based DCCモデルやCopula GARCHモデルなど様々な呼び方があるそうです．先ずはコピュラを実際のデータに適用して，どのような結果を得ることができるのかについてみていきます．そして，その後にコピュラを組み込んだ多変量の時系列モデルについて説明していきます．

コピュラ

今回は，コピュラの中で基本となる正規コピュラとtコピュラ，そしてアルキメディアンコピュラと呼ばれるクラスに属するクレイトンコピュラ，グンベルコピュラ（ガンベルと呼ぶときもあります）の全部で4つのコピュラを対象にします．前回にも説明をしましたが，自身が持つデータに最適なコピュラを選ぶことはとても難しい問題です．今回は多く用いられている方法の一つであるAICによる比較により，どのコピュラが相応しいかを決定したいと思います．コピュラの推定を行う前に，簡単にそれぞれのコピュラについて簡単な説明をします．

  1. 正規コピュラ：相関パラメーター\( \rho \)によって相関を表現します．裾部分の相関についてはあまりうまく表現できません．

  2. tコピュラ：自由度\( \nu \)と相関パラメーター\( \rho \)によって相関を表現します．上側と下側両方の裾部分の相関について表現することができます．しかし，その相関の仕方は上側と下側で同じ強さであることに注意する必要があります．自由度が無限大に近づくことで正規コピュラも表現することができます．

  3. クレイトンコピュラ：下側の裾部分に着目したようなコピュラです．下側の裾部分の相関を捉えることはできますが，上側は捉えることができません．

  4. グンベルコピュラ：上側の裾部分に着目したようなコピュラです．上側の裾部分の相関を捉えることはできますが，下側は捉えることができません．

コピュラのパラメーター推定

RのVineCopulaパッケージを用いて，4つのコピュラのパラメーターの推定と同時にAICの比較を行います．ここで用いるのは最尤推定法ですが，一つ問題があります．それは「それぞれの変数の周辺分布をどのように決めるのか」という問題です．解決の方法としては色々あります．一つがパラメトリックに周辺分布を決定する方法．他には，ノンパラメトリックに周辺分布を推定する方法もあります（間をとったセミパラメトリックな方法というものもあります）．立場によって使うべき方法が異なりますが，近年の研究では「コピュラを決めるために自分でコピュラをはじめに選択しているのに（今回の場合，4つのコピュラのこと），周辺分布までパラメトリックに決めてしまってはモデルリスクが高すぎやしないか？」という意見が多いような気がします．そのため，ノンパラメトリックな方法で周辺分布を推定する研究が増えているように感じます．ここでもこの意見に賛成し，ノンパラメトリックな方法で進んで行きましょう．（ちなみに，最尤推定法ではなくケンドールのタウなどを用いてパラメーターを逆算する方法もあります．）．周辺分布はRのksパッケージを用いて推定を行います．

d_SP500 <- kcde(Dat[,1])  #S&P500の累積分布関数を推定
d_NIKKEI <- kcde(Dat[,2])  #NIKKEI225の累積分布関数を推定
uni_SP500_data <- predict(d_SP500,x=Dat[,1])  #データを周辺分布を用いて変換する
uni_NIKKEI_data <- predict(d_NIKKEI,x=Dat[,2])  #データを周辺分布を用いて変換する

copula_select <- BiCopSelect(uni_SP500_data,uni_NIKKEI_data,indeptest = TRUE,familyset = 0:4) #コピュラの推定
summary(copula_select)

Family
------ 
No:    2
Name:  t

Parameter(s)
------------
par:  0.12  (SE = 0.03)
par2: 6.38  (SE = 1.39)

Dependence measures
-------------------
Kendall's tau:    0.08 (empirical = 0.08, p value < 0.01)
Upper TD:         0.05 
Lower TD:         0.05 

Fit statistics
--------------
logLik:  23.9 
AIC:    -43.8 
BIC:    -33.64

VineCopulaパッケージのBiCopSelect関数はコピュラのパラメーター推定と，AICによる比較を自動で行ってくれます．その結果より今回選択した4つのコピュラの中で一番AICが低かったコピュラは，tコピュラであることがわかります．パラメータのpar=0.12が相関パラメーターであり，par2=6.38が自由度を表しています．相関パラメーターは正であり，自由度も大きくないことから，S&P500のlog-returnとNIKKEI225のlog-returnは正の関係性があり，片方が大きく変動する時には，もう片方も大きく変動する傾向がある，ということがわかりました．同時に大きく変動する傾向の強さを表しているのが，Upper TDとLower TDです．TDとはTail Dependenceのことで，日本語では裾従属性と呼ばれます．この数値が高いほど，同時に大きく変動する傾向が強いことを意味しています．上述しましたが，今回選ばれたtコピュラは上側と下側の裾の相関の強さは同じとなるコピュラであり，結果からも0.05と同じになっていることが確認できます．

DCDモデル

では次にDCDモデルについてみていきます．DCDモデルはDCCモデルと似ていますが，コピュラを含んでいるので残差項が従う分布が異なります．DCCモデルは相関行列と分散が時間によって変化することで，全体の関係性の時間変化について表現しています．DCDモデルは分散は時間変化しますが，相関自体は変化しません．しかし，その相関の表現の幅をコピュラを用いることでとても広くしていることが特徴です．DCCモデルについては，前の記事を参考にしてください．S&P500と日経225の時系列分析（DCCモデル・多変量GARCH） - GRGと金融工学・統計解析

DCDモデルは次のように表せます．

{ \displaystyle {\bf y(t) = \mu (t) + u (t) =\mu (t) + D(t) \epsilon (t)} \dots (1) }

{\bf D(t)} = diag(h_{1,1} (t), h_{2,2} (t), \dots, h_{n,n} (t) )^{\frac{1}{2}} \dots (3)

{ \displaystyle \epsilon (t) \sim iid C(F_1 , F_2) }

ここで，\( \epsilon (t) = {\bf D(t)^{-1}u(t)} \)です．では，rmgarchパッケージを用いてパラメーター推定を行なってみましょう．今回は代表として，先ほど選ばれたtコピュラを用いて推定を行なっていきます．

spec_garch = ugarchspec(mean.model = list(armaOrder = c(1,1)), variance.model = list(garchOrder = c(1,1), model = "eGARCH"), distribution.model = "norm") #前回フィットが良かったEGARCHモデルを用いる
spec_copula = cgarchspec(uspec = multispec( replicate(2, spec_garch) ),  distribution.model = list(copula = "mvt", method = "ML", time.varying = FALSE, transformation = "empirical")) #最尤推定で，時間変化なし，周辺分布はノンパラメトリック（経験分布）
fit_dcd_garch = cgarchfit(spec_copula, data = Dat,solver.control=list(trace=1)) #DCDモデルのパラメーター推定
fit_dcd_garch

*-------------------------------------------------*
*                  Copula GARCH Fit               *
*-------------------------------------------------*

Distribution        :  mvt
No. of Parameters   :  16
[VAR GARCH CC]      : [0+14+2]
No. of Series       :  2
No. of Observations :  1184
Log-Likelihood        :  7613.743
Av.Log-Likelihood :  6.431 

Optimal Parameters
---------------------------------------------------
                        Estimate  Std. Error
[SP500_return].mu       0.000436    0.000071
[SP500_return].ar1      0.967133    0.009397
[SP500_return].ma1     -0.989801    0.000019
[SP500_return].omega   -1.079997    0.010065
[SP500_return].alpha1  -0.297082    0.021932
[SP500_return].beta1    0.889947    0.000205
[SP500_return].gamma1   0.060124    0.006241
[NIKKEI_return].mu      0.000280    0.000456
[NIKKEI_return].ar1    -0.564162    0.016720
[NIKKEI_return].ma1     0.513227    0.016645
[NIKKEI_return].omega  -0.520552    0.192267
[NIKKEI_return].alpha1 -0.100018    0.044924
[NIKKEI_return].beta1   0.938522    0.022040
[NIKKEI_return].gamma1  0.226924    0.039240
[Joint]C1               0.150311    0.029645
[Joint]mshape          21.530978   13.061260
                           t value Pr(>|t|)
[SP500_return].mu       6.1718e+00 0.000000
[SP500_return].ar1      1.0292e+02 0.000000
[SP500_return].ma1     -5.1200e+04 0.000000
[SP500_return].omega   -1.0730e+02 0.000000
[SP500_return].alpha1  -1.3546e+01 0.000000
[SP500_return].beta1    4.3407e+03 0.000000
[SP500_return].gamma1   9.6342e+00 0.000000
[NIKKEI_return].mu      6.1457e-01 0.538840
[NIKKEI_return].ar1    -3.3741e+01 0.000000
[NIKKEI_return].ma1     3.0834e+01 0.000000
[NIKKEI_return].omega  -2.7074e+00 0.006780
[NIKKEI_return].alpha1 -2.2264e+00 0.025989
[NIKKEI_return].beta1   4.2583e+01 0.000000
[NIKKEI_return].gamma1  5.7830e+00 0.000000
[Joint]C1               5.0703e+00 0.000000
[Joint]mshape           1.6485e+00 0.099258

Information Criteria
---------------------
                    
Akaike       -12.834
Bayes        -12.765
Shibata      -12.834
Hannan-Quinn -12.808


Elapsed time : 4.794555

ここでは，[Joint]と書いてある部分に着目してみます．C1が相関パラメーター，mshapeが自由度を表しています．相関パラメーターは先ほどの結果と同じような大きさですが，自由度は大きくなっていることがわかります．先ほどは個々の分散について時間変化を考えていなかったため，同時に大きく動く可能性について自由度を小さくすることで調整していたと考えられます（自由度を小さくすると，裾部分の相関が増えるため）．しかし分散が時間変化することを許すことで，同時に大きく動く可能性を分散によって表現することができるようになり，自由度が大きくなったと考えられます．

Dynamic DCD モデル

DCDモデルだけでもとても柔軟な表現をすることができますが，やはり相関がどのように時間変化するか，といったことを調べることは重要です．DCDモデルのコピュラのパラメーターを変化させることで相関の動きを表すことができたら，DCCモデルよりも広い範囲の相関を見ることができます．先ほどのRのコードを見た時にお気づきかもしれませんが，rmgarchパッケージはその時間変化を考慮したDynamic DCDモデルの推定も行うことができます．先ほどのコードのtime.varying=FALSEの部分をTRUEにするだけで，実際に分析を行うことができます．用意したデータを用いて推定した結果は次のようになりました．

*-------------------------------------------------*
*                  Copula GARCH Fit               *
*-------------------------------------------------*

Distribution        :  mvt
DCC Order           :  1 1
Asymmetric          :  FALSE
No. of Parameters   :  18
[VAR GARCH DCC UncQ]: [0+14+3+1]
No. of Series       :  2
No. of Observations :  1184
Log-Likelihood        :  7613.966
Av.Log-Likelihood :  6.431 

Optimal Parameters
---------------------------------------------------
                        Estimate  Std. Error
[SP500_return].mu       0.000436    0.000072
[SP500_return].ar1      0.967133    0.009710
[SP500_return].ma1     -0.989801    0.000020
[SP500_return].omega   -1.079997    0.009890
[SP500_return].alpha1  -0.297082    0.022440
[SP500_return].beta1    0.889947    0.000213
[SP500_return].gamma1   0.060124    0.006279
[NIKKEI_return].mu      0.000280    0.000452
[NIKKEI_return].ar1    -0.564162    0.016562
[NIKKEI_return].ma1     0.513227    0.016645
[NIKKEI_return].omega  -0.520552    0.188960
[NIKKEI_return].alpha1 -0.100018    0.045090
[NIKKEI_return].beta1   0.938522    0.021669
[NIKKEI_return].gamma1  0.226924    0.039257
[Joint]dcca1            0.013889    0.027853
[Joint]dccb1            0.679697    0.813742
[Joint]mshape          22.527427   31.755845
                           t value Pr(>|t|)
[SP500_return].mu       6.0610e+00 0.000000
[SP500_return].ar1      9.9599e+01 0.000000
[SP500_return].ma1     -5.0386e+04 0.000000
[SP500_return].omega   -1.0920e+02 0.000000
[SP500_return].alpha1  -1.3239e+01 0.000000
[SP500_return].beta1    4.1755e+03 0.000000
[SP500_return].gamma1   9.5759e+00 0.000000
[NIKKEI_return].mu      6.1906e-01 0.535879
[NIKKEI_return].ar1    -3.4064e+01 0.000000
[NIKKEI_return].ma1     3.0833e+01 0.000000
[NIKKEI_return].omega  -2.7548e+00 0.005872
[NIKKEI_return].alpha1 -2.2182e+00 0.026543
[NIKKEI_return].beta1   4.3311e+01 0.000000
[NIKKEI_return].gamma1  5.7805e+00 0.000000
[Joint]dcca1            4.9865e-01 0.618026
[Joint]dccb1            8.3527e-01 0.403564
[Joint]mshape           7.0939e-01 0.478080

Information Criteria
---------------------
                    
Akaike       -12.833
Bayes        -12.760
Shibata      -12.833
Hannan-Quinn -12.805


Elapsed time : 5.495759

AICを見ると，ほんの少しではありますが，時間変化なしのDCDモデルの方が低いという結果になりました．このことから，時間変化を考慮することは大きく精度が向上するわけではないことがわかりました（今回のデータセットの場合は，ということに注意）．しかし，Dynamicにすることで相関の時間変化について図示することができるようになります．その結果が次のグラフです．DCCモデルと比較するために下図に，前回の結果を載せてみました（上図のグラフは相関パラメーターのものなので実際にはDCCモデルと直接比較することは避けた方がいいかもしれません．）．

f:id:GRG:20170318235033j:plain

f:id:GRG:20170312232227j:plain

tコピュラの相関パラメーターの変動がDCCモデルの結果と同じように大きいということは，変数同士が同時に大きく変動する傾向も大きく動いていることも意味しています．なぜなら，先ほど説明したTail Dependenceは相関パラメーターと大きな関係があるからです．tコピュラの場合，相関パラメーターが高くなるとTDも高くなります．

Dynamic DCDモデルはとても柔軟なモデルであり，これ以上追加したら複雑なモデルになりすぎてしまいます（できないことはないが，パラメーターの推定精度が心配）．では，Dynamic DCDモデルは利点だらけのモデルなのかというと，やはりそうではありません．はじめに挙げられる欠点としては，DCDモデルはコピュラの欠点を引き継いでいるという点です．そして，変数が増えるとパラメータ数が急激に増えるので，パラメーター推定に時間がかかったり，推定の精度が悪くなりがちです．このような欠点を加味すると，CCCモデルやDCCモデルの方を用いる方がいい時もありますので注意が必要です．

まとめ

今回ははじめに前回触れたコピュラを実際のデータに当てはめて結果を確認しました．その結果tコピュラが選択した4つのコピュラの中で最適であることがわかりました．その後，DCDモデルやDynamic DCDモデルのパラメーター推定を行い，時間変化を考慮しても今回はAICは大きく向上しなかったことや，DCCモデルの結果と同様に相関が時間によって大きく変化していることが明らかとなりました．DCDモデルは利点も多いですが，欠点も少なからずあるので，注意して使用していきましょう．

次回の予定は考えていませんが，興味があるものがたくさんあるので，それらの紹介をしていけたらなぁと思います．