時系列分析というと、SARIMAモデルや状態空間モデルなどがよく使われているかと思います。 私自身、これらのモデルについて一つの系列データについてモデルを適用したことはありますが、複数の系列データが影響するようなモデルについては扱ってきませんでした。
今回はある系列データが予測対象の系列データに影響を与えている状況を考え、これをベイズ構造時系列モデルが適用して考えてみたいと思い、実際にやってみたので、今回はそのメモです。
ベイズ構造時系列モデル
構造時系列モデル
まずは構造時系列モデルについて簡単に確認します。
構造時系列モデルは、注目している時系列を、直接解釈可能な複数の時系列の和で表すモデルで、任意の時系列の足し合わせであるため解釈性が高いというメリットがあります。 時系列データに対する状態空間モデルの例(ローカルレベルモデル・構造時系列モデル) | AVILEN AI Trend
参考にした論文によれば、これを式で表すと下記のようになります。
観測方程式は、
状態方程式は、
となっています。この辺は普通の状態方程式・観測方程式を想像してもらって問題ないと思います。
また、Googleのブログでは下記のように書かれていました。
この場合は、潜在変数を仮定せず、状態方程式と観測方程式が一つになっている場合を仮定しているようです。 この場合、注目する系列は複数の時系列の和で表現され、個別の時系列が周期性やトレンドを表現することで、全体の系列の動きの解釈が容易になります。
ベイズ構造時系列モデルの概要
今回のテーマになった、ベイズ構造時系列モデルについて確認します。 こちらは参考文献によれば、下記の論文が元になっているようです。
理屈についてざっくりとは下記にまとめを書いてみますが、理解が怪しい部分があるので詳細が気になる方は上の論文を直接ご確認ください。
問題意識
構造時系列モデルは、理論上多数の系列データを考慮して着目する系列データを表現することができる一方、適切でない系列データを説明変数に採択しまうことによるオーバーフィッティングの恐れがあります。 そのため、変数選択を適切に行うことが非常に重要になってきます。
アプローチ
上記の問題意識を解決するために、spike-and-slab事前分布を仮定します。 この事前分布をMCMCによるサンプリングが可能で、ベイズ平均化法を適用してシミュレーションを行うことで変数選択を実現することを考えます。
spike-and-slab事前分布
変数選択のアプローチにもいくつかと思いますが、この手法ではspike-and-slab事前分布を使用したベイジアン変数選択を行うようです。
それぞれの説明変数について、係数が1 or 0となると考えることで変数選択を実現する事を考えるようです。*1
このとき、係数が1になる確率を使用することで、説明変数の係数はベルヌーイ分布を使用することで表現することが出来るようになり、
のように表現することができるようになります。
ベイズ平均化法
spike-and-slab事前分布を用いることで、どの説明変数が使用されるかを確率で表現できるようになりました。 次は、これを使って実際にどの説明変数が使用されるかをシミュレーションすることで最終的な説明変数の係数を得ることを考えます。
この場合、MCMC等を使用することでどの説明変数が採択されるかを確率的にシミュレーションすることが可能です。 これを使って、モデルを多数シミュレーションによって作成し、その平均を最終的なモデルとして使用するベイズモデル平均化法（BMA）を使用することで、説明変数の係数を求めることで最終的なモデルを得ます。
BMAはランダムフォレストにも使用されているようで、ランダムフォレストの学習過程のようにモデルが多数作られていくと考えて問題なさそうです。
使えそうなライブラリ
上で紹介した難しい理屈については私はあまり良く分かっていませんが、ベイズ構造時系列モデルを構築する場合には、下記のライブラリが使用できるようです。
- Python
- TensorFlow Probability
- Prophet
- R
- bsts
簡単な例で良ければ、Python・Rともにサンプルコードに習って書けばできそうです。
サンプルコード
実際にベイズ構造時系列モデルのサンプルコードを見ていきます。
TensorFlow Probability
Googleのブログにあったやつをそのままはこんな感じ。
CO2の予測は状態空間モデルでもやったことがありますが、電力需要予測に関しては気温の時系列データを説明変数として使用しています。 こういったことも表現できるんですね。
自分でもやってみる
上のサンプルコードを参考に、自分でもやってみたいと思います。 今回は、このブログに関する検索履歴について、Google search consoleからエクスポートしたデータ使ってモデリングをやってみたいと思います。
データはざっくりこんな感じになってます。
このときの、最後28日間におけるページリンクのクリック数をモデリングすることを考えます。
基本的にGoogleのブログと同じ流れでやるとこんな感じにモデリングされました。
信頼区間内に大体観測値が収まっているので、予測としては問題なさそうですね。 予測を個々の要素に分解した結果はこんな感じです。
何れにせよ、信頼区間内に収まっているように見受けられるので、ある程度妥当なモデリングにはなっているのではないかと思います。
書いたコードの残骸
使ったコードやらデータやらはこちらです。誰得なのかはよくわかりませんが、一応。
参考文献
感想
久しぶりにベイズっぽいものをやったので、非常に難しかったです。
イメージは、景気指標とか、webサービスの分析は、説明変数の候補が多くなりそうなんで、もしかしたら使えそうですね。 実際に自分のブログの検索流入についてモデリングしてみましたが、そこそこ合ってそうで良かったです。 ただし、実際には表示回数や検索順位は既知ではないので、本当の意味で予測する場合にはそこの予測も含めて考える必要があるのかなとは思いました。
機会があったら検討してみたいと思います。
*1:この辺あまり良くわかっていないので、間違ってたらごめんなさい
