• @y_yammt
• 某Web系エンジニア

• 入門ベイズ統計 第7章「ベイズ更新とカルマン・フィルター」を見ていきます。
• 7.1 リアル・タイム推定
• 7.2 カルマン・フィルターのモデル
• 7.3 カルマン・フィルターの漸化式
• 7.4 シミュレーション例 ← 省略します

• 簡単な例として平均値のリアル・タイム推定を次のスライドで紹介します。

$t-1$時点での平均値:

$$\mathbf{\mu}_{t-1} = \frac{1}{t-1}\sum_{i=1}^{t-1}\mathbf{x}_i$$

• 時系列に対する逐次推定法の一種。

• 状態は観測できず、推定しなければならないことに注意。
• 状態と観測値の例:

状態:

製品の真の不良率

観測値:

サンプリングによって得られた不良率

• 二段構成
• システム方程式に従い、1つ前の状態から今の状態を予測する。
• 今の観測値から、今の状態を補正する。

• システム方程式
$$\theta_t = G_t \theta_{t-1} + w_t$$

$G_t$: (科学的に)すでに知られた量 (例: 物理の運動方程式)。

$w_t$: 平均$0$、分散$\tau^2_t$の正規分布に従う。

• 観測方程式
$$y_t = F_t \theta_{t-1} + v_t$$

$F_t$: (科学的に)すでに知られた量 (例: 物理の運動方程式)。

$v_t$: 平均$0$、分散$\sigma^2_t$の正規分布に従う。

• 以下の公式は証明なしで使います (PRML 2.3を見てね)

• 1変数の正規分布ですと次のようになります。

• とりあえず以下の導出をします。
• システム方程式からの過去の状態から現在の状態への予測

• $(t-1)$時点での状態$\theta_{t-1}$の推定値が$\hat{\theta}_{t-1}$であるとします。また、ベイズ的に$\theta_{t-1}$の真の値はこれを中心とした分散$\delta^2_{t-1}$の正規分布に従うとします。つまり、 $$\mathcal{N}(\theta_{t-1} | \hat{\theta}_{t-1}, \delta^2_{t-1})$$
                                  ... $p(\theta_{t-1})$
• $\hat{\theta}_{t-1}$にしても$\delta^2_{t-1}$にしても初期値が決まれば、(後々の導出で)連鎖的に決めることができることがわかるので、特に気にせず進めていきます。

• システム方程式の定義に戻ると $\theta_t = G_t \theta_{t-1} + w_t$
• $w_t$: 平均$0$、分散$\tau^2_t$の正規分布に従うので、別の書き方をすると、 $$\mathcal{N}(\theta_t | G_t \theta_{t-1}, \tau^2_t)$$

• $(t-1)$時点での状態$\theta_{t-1}$は、 $$\mathcal{N}(\theta_{t-1} | \hat{\theta}_{t-1}, \delta^2_{t-1})$$
                                  ... $p(\theta_{t-1})$
• $\theta_{t-1}$が与えられた場合の、$y_t$を観測する前の$\theta_t$は、 $$\mathcal{N}(\theta_t | G_t \theta_{t-1}, \tau^2_t)$$
                                  ... $p(\theta_t | \theta_{t-1})$
• 以上より、以下の関係を得ます。 $$\mathcal{N}(\theta_t | G_t \hat{\theta}_{t-1}, G^2_t \delta^2_{t-1} + \tau^2_t)$$
                                ... $p(\theta_t)$

• システム方程式から得られる代表値(平均値)で一点先予測できるのでは?
• $\tilde{\theta}_t \equiv G_t \hat{\theta}_{t-1}$ のように。
• 推定はできるが、いまだ$(t-1)$時点であって、本格的な推定には$t$時点の観測値$y_t$を待たねばならない。

• $t$時点で観測値$y_t$が得られたとします。
• この値は一つ前の状態からシステム方程式と観測方程式から得られる代表値(平均値):$$\tilde{y}_t \equiv F_t \tilde{\theta}_t = F_t G_t \tilde{\theta}_{t-1}$$とは異なります(が固定された推定値としては使っている)。

• $y_t$と$\tilde{y}_t$の差は、$$e_t \equiv y_t - \tilde{y}_t = F_t \theta_t + v_t - F_t G_t \tilde{\theta}_{t-1} = F_t (\theta_t - G_t \tilde{\theta}_{t - 1}) + v_t$$
• $\tilde{y}_t$を固定値と見なしていることに注意。
• 上の式は以下のように書くこともできます。$$\mathcal{N}(e_t | F_t (\theta_t - G_t \hat{\theta}_{t-1}), \sigma^2_t)$$
                                  ... $p(e_t | \theta_{t})$
• $e_t$の値はすでにわかっているので、知りたいのは $p(\theta_{t} | e_t)$ → 正規分布のベイズ更新則

• $p(\theta_{t} | e_t) \propto p(e_t | \theta_t) p(\theta_t)$
• → $w(\theta_{t} | e_t) \propto p(e_t | \theta_t) w(\theta_t)$
• $p(\theta_t)$ のところを、 $y_t$を見る前の$\theta_t$の確率に置き換えました。
• $w(\theta_{t} | e_t)$の求め方については第1章の(1.8.11)を参照してください。

• [注意点]
  $p(e_t | \theta_t)= \mathcal{N}(e_t | F_t (\theta_t - G_t \hat{\theta}_{t-1}), \sigma^2_t)$ については、$\theta_t$ に係数 $F_t$ がついていますので、ベイズ更新を適用する場合は、変数変換 $e'_t = \frac{e_t}{F_t} + G_t \hat{\theta}_{t-1}$ を施してから計算する必要があります。

結果は以下のとおり。

• もう少し簡単に書くこともできます(p.126)が省略します。

• 入門ベイズ統計 第7章「ベイズ更新とカルマン・フィルター」を見ました。
• ちょっと導出が煩雑でしたね。
• なお、平均値のリアルタイム推定はオンライン学習(パーセプトロンとか)の更新分の平均化にも使えますので押さえておくとよいかもしれません。