ファクターの因果性調査 ～グレンジャー因果検定をPythonでやってみよう～

モチベーション

トレーダーは往々にしてマーケットを説明するための指標を作ったり探したりしますが、ここでは、その指標が対象としている銘柄に与えるインパクトを古典的なアプローチであるグレンジャー因果性を用いて検証したいと思います。
仮想通貨市場でもよくビットコインのドミナンスが～とか言われていて、そういった指標が実際にビットコインの将来価格に対して関係があるのかを見るのに役立ちます。
下で色々詳しく解説しておりますが、要するにベクトル自己回帰モデルの統計検定になります。

グレンジャー因果性の解説

グレンジャー因果性とは

(Wikipedia参照: こちらからアクセス)
グレンジャー因果検定は、ある時系列が別の時系列の予測に役立つかどうかを判断するための統計的仮説検定で、1969年に初めて提案された。
通常、回帰は「単なる」相関関係を反映するものだが、クライヴ・グレンジャーは、ある時系列の過去の値が別の時系列の将来の値を予測する能力を測定することにより、経済学における因果関係を検証できると主張した。
「真の因果関係」は非常に哲学的な問題であり、前後に続く事象に因果関係があると仮定する前後即因果の誤謬が起こり得るため、計量経済学者はグレンジャー検定が「予測的因果関係」しか見つけられないものだと主張する。
「因果関係」という用語を単独で使用することは誤用で、グレンジャー因果性は「前後関係（precedence）、またはグレンジャー自身が1977年に主張したように「時間的な関連（temporally related）と説明される方が適切である。
グレンジャー因果性は、XがYを引き起こすかどうかを検定するのではなく、XによってYを予測できるかどうかを検定するものである。

通常、Xの過去の値（およびYの過去の値）に対する一連のt-検定およびF-検定を通じて、Xの値が将来のYの値に関して統計的に有意な情報を提供できると証明された場合に、「時系列Xから時系列Yへのグレンジャー因果がある」と言われる。

方法

$y$と$x$を定常性を持つ時系列とする。$x$から$y$へのグレンジャー因果がないという帰無仮説を検定するために、最初に、$y$の単変量自己回帰モデルに含める適切な（異なる時点の）$y$の値の数を決める：

$$
\displaystyle y_{t}=a_{0}+a_{1}y_{t-1}+a_{2}y_{t-2}+\cdots +a_{m}y_{t-m}+\epsilon_{t}
$$

解説 : 定常性を持つ時系列とはホワイトノイズのような時系列のことを指します。例えば、日経平均株価の終値の時系列データは定常的ではなく、日経平均株価の日次価格変化率の時系列データ（ただし重なりはない）等が定常性を持つ時系列データと言えます。
自己回帰性に関しては、もし自己回帰しなさそうであれば、無難に$y_t = \epsilon_t$とするのが無難かもしれません。

次に、$x$の過去の値を含めることで、自己回帰モデルを拡張する。

$$\displaystyle y_{t}=a_{0}+a_{1}y_{t-1}+a_{2}y_{t-2}+\cdots +a_{m}y_{t-m}+b_{p}x_{t-p}+\cdots +b_{q}x_{t-q}+\epsilon_{t} = a_0 + \sum_{k=1}^{m}a_k y_{t-k} + \sum_{k=p}^{q}b_k x_{t-k}$$

この回帰モデルには、$t$統計量に従って個別に有意であることが示された$x$の過去の値をすべて含める。ただし、$x$の過去の値が集合的に回帰モデルの説明力を強めることが、$F$検定（帰無仮説は「$x$によって集合的に説明力は追加されない」）で示された場合に限る。上記の拡張回帰モデルの表記では、$x$が有意となる最小の時間差が$p$、最大の時間差が$q$となる。

$x$から$y$へのグレンジャー因果がないという帰無仮説が採択されるのは、回帰モデルの中に有意な$x$の過去の値が残らなかった場合のみである。

解説 : ここまで見ると、実際には定常性を持つ時系列データに対して、シフトさせながら関係がありそうな変数を探しているだけと捉えることもできますね。もっと多いファクターがある場合にも、同じように追加していくだけです。

ノンパラメトリック検定

上記の線形性に基づく方法は、平均でグレンジャー因果性を検定するのに適している。ただし、分散などのより高次なモーメントでグレンジャー因果性を検出することはできない。グレンジャー因果性のノンパラメトリック検定は、この問題に対処するように設計されている。これらの検定におけるグレンジャー因果性の定義は一般的なもので、線形自己回帰モデルなどのモデルを仮定しない。グレンジャー因果性のノンパラメトリック検定は、高次モーメントや非線形性などを含む、より優れたパラメトリックモデルを構築するための診断ツールとして使用することができる。

Pythonでグレンジャー因果性をやる方法

ここでは実験として、ファーストリテイリングとソフトバンクグループの因果性について検証しながら、グレンジャー因果性の勉強が出来ればと思います。

今回の変数は対数差分$y_t = \log{(x_t/x_{t-1})}$で計算したいと思います。(価格変化率でも良いですが、価格価格変化率の分布は厳密には左右対称ではなく、対数差分であれば左右対称になるため、対数差分の方が使い勝手が良いケースもあります。実務上における分かりやすさや結果の解釈性等は価格変化率の方がおすすめであるため、ケースバイケースです。)

ライブラリのインポート

import numpy as np
import pandas as pd
from datetime import datetime
import pandas_datareader.data as web
import statsmodels.api as sm
from statsmodels.tsa.api import VAR
from statsmodels.tsa.base.datetools import dates_from_str

データのダウンロード及び整形

_data = dict()

symbols = ["9983.T","9984.T"]
start = datetime(2000, 1, 1)
end = datetime(2021, 12, 3)

for symbol in symbols:
    _df = web.DataReader(symbol, 'yahoo', start, end)
    _data[symbol] = _df["Close"]

data = pd.concat(_data, axis=1)

データのフォーマットは次の通りです。

	9983.T	9984.T
Date
2021-05-31	89170.0	8256.0
2021-06-01	88900.0	8162.0
2021-06-02	88770.0	8109.0
2021-06-03	85100.0	8208.0
2021-06-04	85100.0	8102.0

データを統計検定できるようにする。

# 対数差分
logdiff = np.log(data).diff().replace([np.inf, -np.inf], np.nan).dropna()

# VARモデルのインスタンス化
model = VAR(logdiff, missing='drop')

さて、ここまで準備出来ましたら、あとは検定していきましょう。

results = model.fit(4)

# 9983 => 9984の因果性
test_results = results.test_causality(causing='9983.T', caused='9984.T')
print(test_results.pvalue) # 0.014741228373033606

# 9984 => 9983の因果性
test_results = results.test_causality(causing='9984.T', caused='9983.T')
print(test_results.pvalue) # 0.3556412686636041

結果に着目しますと、9983から9984の因果性のp値はかなり低いですね。これは、5%を有意水準としますと、9984=>9983は棄却出来そうです。つまり、9984の予測に9983をファクターに使うことが出来る一方で逆は使えないことを意味しています。

インパクト調査(インパルス応答)

グレンジャー因果性でファクター変数が変化した時の、目的変数に与える影響を見てみましょう。

インパルス応答に関しては、公式を参考までにお願いします。

period = 10
irf = results.irf(period)
irf.plot()

対角にあるグラフは自己相関で、その他は、$x_i \Rightarrow x_j$へのインパルス応答を表しています。点線は95％信頼区間です。

reference: https://www.statsmodels.org/dev/vector_ar.html#impulse-response-analysis

予測

さて、予測してみます。statsmodelsの予測では基本的に時刻$t$で将来のT期間を推定するため、全ての時刻で推論する必要があります。
そのため、一行で書くのは中々難しいです。一ステップのウォークフォワード法っぽいですね。

from tqdm.auto import tqdm

lag_order = results.k_ar
T = len(logdiff)
original_index = logdiff.index
predictions = []
index = []

for t in tqdm(range(lag_order, T)):
    res = results.forecast(logdiff.values[t-lag_order:t], 1)[0]
    predictions.append(res)
    index.append(original_index[t])

_preds = pd.DataFrame(predictions, index=index, columns=["9983_pred", "9984_pred"])

この予測シグナルの使い方は読者各自にお任せいたします。

結論

ここでは、グレンジャー因果性を見てきました。
ここで紹介した統計検定手法によれば、ファーストリテイリングの直近4期間の変化を見れば、ソフトバンクグループの株価の将来予測が出来る可能性があることを表しています。
モデルがシンプルなだけに非常に使い勝手が良いというメリットがあります。同時に、こういった統計はやり過ぎると、ランダムな結果にも関わらず、良い結果のみを見てしまうというバイアスにもかかりやすいので、結果の使用には注意が必要です。
では、皆さん、良い分析ライフを！