からっぽのしょこ

読んだら書く!書いたら読む!同じ事は二度調べ(たく)ない

3.5.1:粒子フィルタ【白トピックモデルのノート】

はじめに

 『トピックモデルによる統計的潜在意味解析』の学習時のメモです。基本的な内容は、数式の行間を読んで埋めたものになります。本と併せて読んでください。

 この記事では、3.5.1項の粒子フィルタについて書いています。

 数学よく解らない自分が理解できるレベルまで落として数式を書き下していますので、分かる人にはかなりくどいです。

【前節の内容】

www.anarchive-beta.com

【他の節一覧】

www.anarchive-beta.com

【この節の内容】

3.5.1 粒子フィルタ

 粒子フィルタは、時系列データを解析する動的システムにおいて、潜在変数のサンプリングを可能にする技術である。ある時刻において、それまでに与えられた観測データの下での潜在変数の事後分布を重み付きサンプルによって近似する。

 時刻tの観測データ$x_t$は、潜在変数$z_t$に依存する確率分布$p(x_t | z_t)$によって生成されたと仮定する。また、潜在変数$z_t$は1つ前の時刻t-1における潜在変数$z_{t-1}$に依存した確率分布$p(z_t | z_{t-1})$によって生成されたと仮定する。(図3.7(a))

$$ \begin{aligned} x_t &\sim p(x_t | z_t) \\ z_t &\sim p(z_t | z_{t-1}) \end{aligned} $$

・重みの導入

 時刻t+1において、それまでに観測されたデータ$x_{1:t}$が与えられた下で、観測データ$x_{t+1}$を生成する潜在変数$z_{t+1}$の条件付き確率分布$p(z_{t+1} | x_{t:1})$は、提案分布$q(z_{1:t} | x_{1:t})$を$\frac{q(z_{1:t} | x_{1:t})}{q(z_{1:t} | x_{1:t})} = 1$を掛ける形で導入して

$$ \begin{align} p(z_{t+1} | x_{1:t}) &= \int p(z_{t+1}, z_{1:t} | x_{1:t}) dz_{1:t} \\ &= \int p(z_{t+1} | z_{1:t}, x_{1:t}) p(z_{1:t} | x_{1:t}) dz_{1:t} \\ &= \int p(z_{t+1} | z_{1:t}) \frac{ p(z_{1:t} | x_{1:t}) }{ q(z_{1:t} | x_{1:t}) } q(z_{1:t} | x_{1:t}) dz_{1:t} \\ &= \int p(z_{t+1} | z_{1:t}) \omega(z_{1:t}) q(z_{1:t} | x_{1:t}) dz_{1:t} \tag{3.160} \end{align} $$

となる。ここで

$$ \omega(z_{1:t}) = \frac{ p(z_{1:t} | x_{1:t}) }{ q(z_{1:t} | x_{1:t}) } \tag{3.161} $$

とおく。

 提案分布$q(z_{1:t} | x_{1:t})$からのサンプル$z_{1:t}^{(s)} \sim q(z_{1:t} | x_{1:t}) \ (s = 1, 2, \cdots, S)$を用いて、$p(z_{1:t} | x_{1:t})$を

$$ \begin{aligned} p(z_{1:t} | x_{1:t}) &= \frac{ p(z_{1:t} | x_{1:t}) }{ q(z_{1:t} | x_{1:t}) } q(z_{1:t} | x_{1:t}) \\ &= \omega(z_{1:t}) q(z_{1:t} | x_{1:t}) \approx \tilde{p}_S (z_{1:t} | x_{1:t}) \end{aligned} $$

と近似する。$\tilde{p}_S (z_{1:t} | x_{1:t})$は

$$ \begin{align} \tilde{p}_S (z_{1:t} | x_{1:t}) &= \sum_{s=1}^S \frac{ \omega(z_{1:t}^{(s)}) }{ \sum_{s=1}^S \omega(z_{1:t}^{(s)}) } \delta(z_{1:t} = z_{1:t}^{(s)}) \\ &= \sum_{s=1}^S \bar{\omega}(z_{1:t}^{(s)}) \delta(z_{1:t} = z_{1:t}^{(s)}) \tag{3.162} \end{align} $$

である。ここで、$\bar{\omega}(z_{1:t}^{(s)}) = \frac{ \omega(z_{1:t}^{(s)}) }{ \sum_{s=1}^S \omega(z_{1:t}^{(s)}) },\ 0 \leq \bar{\omega}(z_{1:t}^{(s)}) \leq 1,\ \sum_{s=1}^S \bar{\omega}(z_{1:t}^{(s)}) = 1$とおく。(数学的に高度なので読み飛ばし可とのことなので今回はこの理由を読み飛ばした…無念)
 次からは、この重みの逐次更新式を求めていく。

・重みの更新式の導出

 式(3.60)の積分を$\tilde{p}_S (z_{1:t} | x_{1:t})$によって近似する。

$$ \begin{align} p(z_{t+1} | x_{1:t}) &= \int p(z_{t+1} | z_{1:t}) p(z_{1:t} | x_{1:t}) dz_{1:t} \\ &\approx \int p(z_{t+1} | z_{1:t}) \tilde{p}_S (z_{1:t} | x_{1:t}) dz_{1:t} \\ &= \int p(z_{t+1} | z_{1:t}) \sum_{s=1}^S \bar{\omega}(z_{1:t}^{(s)}) \delta(z_{1:t} = z_{1:t}^{(s)}) dz_{1:t} \\ &= \sum_{s=1}^S \bar{\omega} (z_{1:t}^{(s)}) p(z_{t+1} | z_{1:t}^{(s)}) \tag{3.165} \end{align} $$

デルタ関数の性質$\int f(x) \delta(x = a) dx = f(x = a)$を用いている。
 ここで、提案分布を

$$ \begin{align} q(z_{1:t} | x_{1:t}) &= q(z_t, z_{1:t-1} | x_{1:t}) \\ &= q(z_t | x_{1:t}, z_{1:t-1}) q(z_{1:t-1} | x_{1:t}) \\ &= q(z_t | x_{1:t}, z_{1:t-1}) q(z_{t-1}, z_{1:t-2} | x_{1:t}) \\ &= q(z_t | x_{1:t}, z_{1:t-1}) q(z_{t-1} | x_{1:t}, z_{1:t-2}) q(z_{1:t-2} | x_{1:t}) \\ &= q(z_t | x_{1:t}, z_{1:t-1}) q(z_{t-1} | x_{1:t}, z_{1:t-2}) \cdots q(z_1 | x_{1:t}, z_0) q(z_0 | x_{1:t}) \\ &= q(z_t | x_{1:t}, z_{1:t-1}) q(z_{t-1} | x_{1:t-1}, z_{1:t-2}) \cdots q(z_1 | x_1, z_0) q(z_0) \\ &= q(z_0) \prod_{l=1}^t q(z_l | x_{1:l}, z_{1:l-1}) \tag{3.166} \end{align} $$

とおく。ただし図 3.7(a)より、時刻tの観測値$x_t$が1つ前の時刻の潜在変数$z_{t-1}$に影響しないことを考慮して、$q(z_{1:t-1} | x_{1:t}) = q(z_{1:t-1} | x_{1:t-1})$と設定し、置き換えている。

 これを用いて、重みの更新式を導出する。更新式は$\omega(z_{1:t}^{(s)}) = \omega(z_{1:t-1}^{(s)}) F(\cdot)$の形になることから、$\frac{\omega(z_{1:t}^{(s)})}{\omega(z_{1:t-1}^{(s)})}$を求めればよい。よって、式(3.161)より

$$ \begin{align} \frac{\omega(z_{1:t}^{(s)})}{\omega(z_{1:t-1}^{(s)})} &= \frac{p(z_{1:t}^{(s)} | x_{1:t})}{q(z_{1:t}^{(s)} | x_{1:t})} \frac{q(z_{1:t-1}^{(s)} | x_{1:t-1})}{p(z_{1:t-1}^{(s)} | x_{1:t-1})} \\ &= \frac{p(z_{1:t}^{(s)} | x_{1:t})}{p(z_{1:t-1}^{(s)} | x_{1:t-1})} \frac{q(z_{1:t-1}^{(s)} | x_{1:t-1})}{q(z_{1:t}^{(s)} | x_{1:t})} \\ &= \frac{ p(z_{1:t}^{(s)}, x_{1:t}) p(x_{1:t-1}) }{ p(x_{1:t}) p(z_{1:t-1}^{(s)}, x_{1:t-1}) } \frac{ q(z_{1:t-1}^{(s)} | x_{1:t-1}) }{ q(z_t^{(s)} | x_{1:t}, z_{1:t-1}^{(s)}) q(z_{1:t-1}^{(s)} | x_{1:t}) } \\ &\propto \frac{p(z_{1:t}^{(s)}, x_{1:t})}{p(z_{1:t-1}^{(s)}, x_{1:t-1})} \frac{1}{q(z_t^{(s)} | x_{1:t}, z_{1:t-1}^{(s)})} \\ &= \frac{ p(z_t^{(s)}, x_t, z_{1:t-1}^{(s)}, x_{1:t-1}) }{ p(z_{1:t-1}^{(s)}, x_{1:t-1}) } \frac{1}{q(z_t^{(s)} | x_{1:t}, z_{1:t-1}^{(s)})} \\ &= \frac{ p(z_t^{(s)}, x_t | z_{1:t-1}^{(s)}, x_{1:t-1}) }{ q(z_t^{(s)} | x_{1:t}, z_{1:t-1}^{(s)}) } \\ &= \frac{ p(x_t | z_t^{(s)}, z_{1:t-1}^{(s)}, x_{1:t-1}) p(z_t^{(s)} | z_{1:t-1}^{(s)}, x_{1:t-1}) }{ q(z_t^{(s)} | x_{1:t}, z_{1:t-1}^{(s)}) } \tag{3.168} \end{align} $$

【途中式の途中式】

  1. 式(3.161)より置き換える。
  2. $p(\cdot)$と$q(\cdot)$を揃えるために、分母を入れ替える。
  3. 項を変形する。
    • $p(A_1 | A_2, B) = \frac{p(A_1, A_2 | B)}{p(A_2 | B)}$より、前の因子の分母分子をそれぞれ変形する。
    • $p(A_{1:2} | B) = p(A_1, A_2 | B) = p(A_1 | A_2, B) p(A_2 | B)$より、後ろの因子の分母を変形する。
  4. 式を整理する。
    • $z_{1:t}^{(s)}, z_{1:t-1}^{(s)}$に関係する項のみ取り出す。
    • 提案分布を$q(z_{1:t-1}^{(s)} | x_{1:t}) = q(z_{1:t-1}^{(s)} | x_{1:t-1})$と設定したので、約分する。
  5. 式を整理するために、分子の項からtについての変数を取り出しておく。
  6. $\frac{p(A_1, A_2)}{p(A_2)} = p(A_1 | A_2)$より、変形する。


 更に、提案分布を$q(z_t^{(s)} | x_{1:t}, z_{1:t-1}^{(s)}) = p(z_t^{(s)} | z_{1:t-1}^{(s)}, x_{1:t-1})$とすると

$$ \begin{align} \frac{\omega(z_{1:t}^{(s)})}{\omega(z_{1:t-1}^{(s)})} &\propto \frac{ p(x_t | z_t^{(s)}, z_{1:t-1}^{(s)}, x_{1:t-1}) p(z_t^{(s)} | z_{1:t-1}^{(s)}, x_{1:t-1}) }{ q(z_t^{(s)} | x_{1:t}, z_{1:t-1}^{(s)}) } \tag{3.168}\\ &= \frac{ p(x_t | z_t^{(s)}, z_{1:t-1}^{(s)}, x_{1:t-1}) p(z_t^{(s)} | z_{1:t-1}^{(s)}, x_{1:t-1}) }{ p(z_t^{(s)} | z_{1:t-1}^{(s)}, x_{1:t-1}) } \\ &= p(x_t | z_t^{(s)}, z_{1:t-1}^{(s)}, x_{1:t-1}) \\ \omega(z_{1:t}^{(s)}) &= \omega(z_{1:t-1}^{(s)}) p(x_t | z_t^{(s)}, z_{1:t-1}^{(s)}, x_{1:t-1}) \tag{3.169} \end{align} $$

重みの逐次更新式が得られる。
 正規化は次の式で行う。

$$ \bar{\omega}(z_{1:t}^{(s)}) = \frac{\omega(z_{1:t}^{(s)})}{\sum_{s=1}^S\omega(z_{1:t}^{(s)})} $$


参考文献

  • 佐藤一誠『トピックモデルによる統計的潜在意味解析』(自然言語処理シリーズ 8)奥村学監修,コロナ社,2015年.

おわりに

 この本の中で紹介されている『予測に生かす統計モデリングの基本』を読んで、粒子フィルタについて理解して戻ってきました。初学者向けに、数式を端折らず載せた上で言葉でも式変形の理由を説明していて、とても良い本でした。
 が、参考書で説明された粒子フィルタとこの本の粒子フィルタが少し違う(汗)。リサンプリングした後に誤差項を加えないと、偏ったままになってしまうのでは??図3.8を見ると、「提案分布からのサンプリング」がシステムモデル(あるいはシステムノイズ)に対応してるっぽいのだが、よく解らん。次の項では、閾値以下になると過去に戻ってリサンプリングして活性化(若返り)を行うとのことだが、その過去時点にリサンプリングするのとの違いが思い浮かばぬ…。
 という疑問を残しつつ次節へ。

【次節の内容】

www.anarchive-beta.com