⏳
時系列解析初心者が始める時系列解析(4)
統計的モデリングと情報量基準AIC
参考：https://elf-c.he.u-tokyo.ac.jp/courses/381
統計的モデリングとモデル評価
統計的モデリングとは
-
統計的モデル①：情報抽出の「道具」である。
- モデルを作成するために必要な要素（入力）
- 理論
- 経験的知識（ドメイン知識）
- データ
- 作成したモデルから得られる要素（出力）
- 情報抽出 / 知識発見
- 予測シミュレーション
- 制御 / 意思決定
- モデルを作成するために必要な要素（入力）
-
統計的モデル②：「確率分布」で表現する。
-
基本的な形
連続型確率分布 離散型確率分布 正規分布 コーシー分布 ピアソン分布 二項分布 指数分布 二重指数分布 \chi{^2} 多項分布 一様分布 ポアソン分布 -
プラスアルファの様々な情報を取り入れる必要がある。
- 他の変数の情報
- 過去の変動の情報
- 時間経過による変化
-
-
①、②共に、いかに有効活用するかで得られるモデルは異なる。特に、それぞれの課題設定にあった出力となっているかは重要である。
モデル評価
- モデリング時に様々な情報から色々なモデルを構築することになる。
- そのためにもどのモデルが適しているか評価し選択することが重要。
- その尺度として用いられるのが情報量基準である。
Kullback-Leibler情報量
情報量基準とは
- ”モデルの良さ”を予測能力で評価する。
- 予測は点推定ではなく”予測分布”で行う。（信頼区間とか？）
- 分布の近さを”K-L情報量”で評価する。
K-L情報量によるモデルの評価
-
Kullback-Leibler情報量：真の分布とモデルの分布の乖離具合を測る尺度
\begin{align}I(g;f) &= E_Y\log{[{\frac{g(Y)}{f(Y)}}]}\\ &=\int{\log{[\frac{g(x)}{f(x)}}}]dG(x)\\ &= \begin{cases}\int{\log{[\frac{g(x)}{f(x)}]}g(x)dx}~~連続分布モデル\\ ~~ \\ \sum_{i \in J}{g_i \log{[\frac{g_i}{f_i}]}}~~離散分布モデル\end{cases} \end{align} -
K-L情報量の性質
\begin{align}&(i)~I(g;f) ≥ 0\\ &(ii)~I(g;f)=0 \iff g(x)=f(x)\end{align} -
K-L情報量は距離ではないことに注意。そのため、「対称性がある」という距離の性質は成り立たないということになる。
d(g,f)≠d(f,g)
KL情報量の推定と最尤法
K-L情報量の推定
- そもそもK-L情報量に用いられる
g(y) - そのため、
I(g(y);f(y)) - その推定をするにあたり、K-L情報量の便利な性質がある。それが
\log g f
I(g;f) = E_Y\log{\frac{g}{f}} = E_Y\log{g}-E_Y\log{f} - この等式の中で注目すべき点が
E_Y\log{f(Y)} - また、
E_Y\log{g} f - そういった点から、モデリング時には、K-L情報量を求める代わりに、平均対数尤度を使えるらしい。数式で関係を表すと下記の通り。
E_Y\log{f(Y)}:大 \iff I(g;f):小 - ちなみに、K-L情報量は０に近くづくほど良いです。
平均対数尤度E\log{f}
- 平均対数尤度を求めることで、モデリングを評価し、選択につながるということであったが、平均対数尤度においても未知の分布を含んだ数式になっている。
E_Y{\log{f(Y)}}=\int{f(y)dG(y)} - 上記の式における
G(y) dG(y)=g(y)dy
G(y) \rightarrow \hat{G}_n(y) = \frac{1}{n}\sum_{i=1}^{n}{I(y,X_i)} - この式を対数尤度として求めると、、、(
n
- 上記式を改めて平均化し、大数の法則を当てはめると、平均対数尤度に収束していく。
\frac{1}{n}\sum_{i=1}^{n}\log{f(X_i)} \rightarrow E_Y{\log{f(Y)}}
最尤法
- 今までの
f(y) \theta
f(y~|~\theta), ~~ \theta \equiv (\theta_1, \cdots, \theta_k)' - この式を対象に最尤法を適応する。最尤法は近似的にK-L情報量を最小化することになる。
対数尤度：l(\theta) = \sum_{i=1}^n{\log{f(X_i~|~\theta)}} \equiv \log{f(X~|~\theta)}
最尤法：\max_\theta{l(\theta)} \rightarrow \hat\theta = \hat\theta(X) - 主な用語
用語 数式 最尤推定量 \hat\theta 最大対数尤度 l(\hat\theta) 最尤モデル f(y~\lvert ~\hat\theta)
- 主な用語
最尤推定値の求め方
- 尤度方程式を解く
- 数値的最適化を用いる --> 時系列解析においては、データの性質上、主にこちら
最尤推定量の性質
- 尤度方程式
\frac{\partial l(\theta)}{\partial\theta}=0 \theta_0 -
\hat\theta_n n \rightarrow+\infin \theta_0 -
\sqrt{n}(\hat\theta_n - \theta_0) \rightarrow N(0, I(\theta_0)^{-1}) I
- ちなみに、最尤推定量は解を求めた時、その解に確率変数が含まれている状態であり、最尤推定値はその解に試行結果を当てはめ、実際の数値として算出したものである。こちらのURLが参考になるかも。https://detail.chiebukuro.yahoo.co.jp/qa/question_detail/q10215490355
バイアス補正とAICの導出
複数のモデルの比較
- 出てくる用語
| 名称 | 数式 |
|---|---|
| モデル | |
| パラメータ | |
| 最大対数尤度 | |
| Fisher情報行列 | |
| Hessianの期待値 |
-
k - 残念ながらそういうわけにはいかない。なぜなら、低次数モデルと高次数モデルでは、高次数モデルの方が良い結果を出してしまうため。
- 通常の機械学習同様、高次数モデルが予測精度が高いのかというと、そういうわけではない。
-
理由・原因・対策
| 項目 | 内容 |
|---|---|
| 理由 | 最大対数尤度 |
| 原因 | パラメータ推定と平均対数尤度の推定に同じデータを２回用いたため。 |
| 対策 | バイアスを評価し補正する。 |
-
どのようにバイアスを補正するのか？それは、バイアス補正をするための項
b(G) -
b(G) n
\log{f(X~\lvert~\hat\theta(X))}-nb(G)
-
-
b(G) D D
b(G) = E_X[D] = E_X[D_1 + D_2 + D_3]
| それぞれの点 | |
|---|---|
| 対数尤度グラフ上の |
|
| 対数尤度グラフ上の |
|
| 平均対数尤度グラフ上の |
|
| 平均対数尤度グラフ上の |
|
|
||
|---|---|---|
-
D1・D2・D3
b(G)=\frac{1}{n}tr[I(G)J(G)^{-1}] -
この値を情報量基準
IC TIC
IC = -2\log{f(x~\lvert~\hat\theta)}+2nb(G)
TIC = -2\log{f(x~\lvert~\hat\theta)}+2tr[I(G)J(G)^{-1}] -
ただし、ここから使うのは
AIC -
k
AIC = -2\log{f(x~\lvert~\hat\theta)}+2k
-
AICとTICの関係
nb(G)=tr[I(G)J(G)^{-1}] k
-
複数作成したモデル（モデル族）内に真の分布を含む場合、下記の関係が成り立つらしい。
I(\theta)=J(\theta) -
上記の性質が満たされるのであれば、下記の式が得られる。
- この時、
I_k
\begin{align}nb(G) &= tr[I(G)J(G)^{-1}]\\ &= tr[I_k]\\ &= k\end{align}
- この時、
-
ちなみに、
TIC AIC I(\theta)=J(\theta) - ただ、これは誤解で、
I(\theta),J(\theta) - また、高次モーメントを含んでいるため、分散が大きくなってしまう。
- ただ、これは誤解で、
AICによるモデルの選択例
- どういったふうに選択するのかは、動画を見て欲しい。基準としては、AICの値がいかに小さいか。その時の次数はいくつかというところである。
モデルの予測誤差分散
-
予測誤差については、以下の等式が成り立つ。
- バイアス：モデルの悪さ
- 分散：モデルの不安定さ
予測誤差=バイアス+分散
-
つまり、AICが最小のモデルというのは、、、
- バイアスと分散を適度に小さくしたモデル
- それは、期待予測誤差を最小にするモデルを選択したことになる。
その他の情報量基準
- 情報量基準というのは、K-L情報量のバイアスを補正して導くもの。
AIC_c
-
AIC - それへの対策として、有限修正を施す。具体的にはバイアス補正項を置き換える。そうすることで、データ少ない場合でもある程度の精度を出せるようにする。
nb_c(G) = \frac{n(p+1)}{n-p-2}
AIC_c = -2\log{f(X~\lvert~\hat\theta(X))}+2\frac{n(p+1)}{n-p-2}
GIC
- 特徴として、統計的汎関数として定義できる任意の推定量に適用可能。つまり、最尤推定量に限らず、他の推定量を用いることができる。
- 例）ロバスト推定量、正則化法、ベイズ推定の一部
- また、高次モーメントを含むようなものでもうまく補正してくれる。
- ただし、計算が複雑で面倒。
GIC = -2\log{f(X~\lvert~\hat\theta(X))}+2nb(\hat{G}_{GIC})
EIC
-
バイアス補正量を解析的にではなく、ブートストラップによって数値的に求める。
-
特長としては以下3つ。
- 解析的近似が不要（今までのような計算）
- 計算実装も比較的容易
- 最尤推定量以外の様々な推定量やモデルに適用可能
-
弱点としては、ブートストラップ法のため、データ生成・推定を繰り返す。すると、計算量が多くなってしまう。
ABIC
-
良いモデルを求める方法には、モデル選択だけでなく、「パラメータに制約を課す」という方法が存在する。
-
そこで登場するのがベイズモデル。
-
f(x~\lvert~\theta) \theta \pi{(\theta~|~\lambda)} \lambda x
p(x~|~\lambda)=\int{f(x~|~\theta)\pi(\theta~|~\lambda)d\theta} -
この式をパラメータを
\lambda p \lambda
\begin{align} ABIC &= -2\log[\max_\lambda{p(x|\lambda)}]+2q\\ &= -2\max_{\lambda}\log{\int{f(x~|~\theta)\pi(\theta~|~\lambda)d\theta}} +2q \end{align} -
注意事項としては、
\lambda \theta
BIC
- 周辺尤度を用いて、それを近似的に最大化するようなモデルを求めるというもの。これまでの情報量基準とは少しロジックが異なるらしい。
- 周辺尤度については、事前分布
\pi(\theta_m)
p(x) = \int{f(x|\theta_m)\pi(\theta_m)d\theta_m} - 以上を用いて、
BIC \thickapprox
\begin{align}BICm = -2\log{p(x)} &=-2\log{\int{f(x|\theta_m)\pi(\theta_m)d\theta_m}}\\ &\thickapprox -2\log{f(x|\hat\theta)}+m\log{N} \end{align}
交差検証法（Cross Validation）
- 全データを推定用データと評価用データに分ける
- 推定用データでモデルを推定
- 評価用データでモデルを評価（予測２乗誤差など）
- 1の分割の仕方を変えて、すべての場合について、2,3を繰り返し、評価量の平均を求める。
- 分割の仕方
- Leave-one-out：1個のデータだけ評価に用いる
- k分割法：全体のデータをk分割し、そのうちの1つを評価に用いる
AICに関する批判
- よくある批判として「AICには真の次数への一致性がない」ということが挙げられる。
- 一致性：サンプル数を増やすことで、本来のものに近似する
- ただ、それは、AICにとって最も重要な問題ではない。以下がその理由。
- モデリングの目的は、「良い」モデルを求めることで、「真の」モデルを求めることではない。
- 次数の一致性は良いモデルを求めるための「必要条件」でも「十分条件」でもない。
- 「真」の次数は一般に存在しない。存在したとしても、真の次数が推定されたモデルが予測にも良いとは限らない。
- 「真」の次数より高くてもパラメータが一致性を持てば、モデルは一致する。
Discussion