粒子フィルターの連続極限

どんな過程が現れるか?

Particles
Process
Author

司馬 博文

Published

1/23/2024

Modified

2/12/2024

概要
粒子フィルターを拡散過程に対して適用することを考える.拡散過程の Euler-Maruyama 離散化に対して構成された粒子フィルターの,タイムステップを \(0\) にする極限 \(\Delta\searrow0\) での振る舞いを議論する.

リサンプリングと粒子フィルターは \((\Omega,\mathcal{F},\operatorname{P})\) 上に定義されているとし,この確率測度 \(\operatorname{P}\) に関する期待値を \(\operatorname{E}\) で表す.

一方で,観測対象は \((\Omega',\mathcal{F}',\operatorname{P}')\) 上の過程 \(z:\mathbb{R}_+\times\Omega'\to\mathbb{R}^d\) とする.

2つの分布の間の関係を導く.

1 リサンプリング

1.1 形式化

積空間 \([N]^{[N]}\)附番の空間 といい,元を \(a\in[N]^{[N]}\) で表す.

積空間 \(\mathbb{R}_+^N\setminus\{0\}^N\)荷重の空間 といい,元を \(w\in\mathbb{R}_+^N\setminus\{0\}^N\) で表す.簡単のため,以降は \(\mathbb{R}_+^N\) とも表してしまう.ノルム \[ \|w\|_1:=\sum_{i=1}^Nw_i \] を考える.

定義(リサンプリング)
  1. 粒子数 \(N\) に関する リサンプリング法 \(r\) とは,荷重の空間 \(\mathbb{R}_+^N\) から附番の空間 \([N]^{[N]}\) への確率核 をいう.すなわち,荷重に応じて,附番の空間上の確率分布を対応させる関数 \[ r:\mathbb{R}_+^N\to\mathcal{P}([N]^{[N]}) \] をいう.

  2. リサンプリング法 \(r\)不偏 であるとは,任意の番号 \(j\in[N]\) に対して,その番号を得る粒子数の期待値が,荷重の通りになることをいう:\(A\sim r(w)\) ならば, \[ \begin{align*} \operatorname{E}[\#A^{-1}(j)]&=\operatorname{E}\left[\sum_{i=1}^N1_{\left\{A(i)=j\right\}}\right]\\ &=N\frac{w_j}{\|w\|_1}\quad(j\in[N]). \end{align*} \]

定義は (Chopin et al., 2022, p. 3199) を参考にした.リサンプリング法の定義は,確率核 \(r\) の終域をどうするかで議論があると思う.初めは \(r:\mathbb{R}_+^N\to\mathcal{P}([N])\) と定めたくなる.しかしこれでは \(N\) 回の非独立的なサンプリング を表現できない.\(r(w^{1:N})^{\otimes N}\) とするわけには行かないし,ベクトル \(N\cdot r(w^{1:N})^\top\) が個数の分布になるようにサンプリングせよ,ということだったら確率核になっていない.結局確率核は \(r:\mathbb{R}_+^N\to\mathcal{P}([N]^{[N]})\) と取るしかない.すると,不偏性の表現が難しくなる.

問題構成がすごく PAC-Bayes っぽい.というより,Gibbs classifier っぽい.

離散空間上の測度 \(\mu\in\mathcal{M}^1([N]^{[N]})\) は,その一点集合上での値 \[ f_\mu(a):=\mu(\{a\}) \] で一意に定まるため,\(\mathcal{M}^1([N]^{[N]})\)\(\mathrm{Map}([N]^{[N]},\mathbb{R}_+)\) と,\(\mathcal{P}([N]^{[N]})\)\([N]^{[N]}\) 上の確率質量関数全体の集合と同一視できる.

さらに,関数空間 \(\mathrm{Map}([N]^{[N]},\mathbb{R}_+)\) に各点収束の位相を入れたものと,弱位相を備えた \(\mathcal{M}^1([N]^{[N]})\) とは同相である

これを,\([N]^{[N]}\) に限らず,任意の可算な離散空間 \(S\) において示す.すると \(S\) 上では,弱収束は全変動収束に等しい.1

一般に (Scheffé, 1947) の定理より,密度が各点収束するならば,分布は全変動収束,特に弱収束する.逆に,\(S\) 上の分布が収束するとき,任意の一点集合 \(\{x\}\subset S\) は開かつ閉であるから境界は空集合であり,(Alexandroff, 1940) の特徴付けから,密度(というより確率質量関数)も収束する.

よって,\(\mathcal{M}^1([N]^{[N]})\) 上の弱収束と,\(\mathrm{Map}([N]^{[N]},\mathbb{R}_+)\) 上の各点収束とは同値.いずれも第2可算だから,点列の収束が等しいことは位相も等しいことを含意する.

以降,\(\mathcal{M}^1([N]^{[N]})\) は全変動ノルムが備わったノルム空間と見る.これは,写像の空間 \(\mathrm{Map}([N]^{[N]},\mathbb{R}_+)\) に各点収束の位相が備わった空間と同一視できる.

1.2 荷重が一様になる極限でも安定なリサンプリング法

1.2.1 モチベーション

リサンプリング法 \(r:\mathbb{R}_+^N\to\mathcal{P}([N]^{[N]})\) に対して,荷重が一様に近づく(=殆どリサンプリングをしなくなっていく)極限における性質を,連続関数 \[ \begin{align*} e:\mathbb{R}_+^N\times(0,1)&\to \mathbb{R}_+^N\\ (v,\Delta)\quad&\mapsto e^{-\Delta v}:=(e^{-\Delta v^1},\cdots,e^{-\Delta v^N}) \end{align*} \] を通じて調べる: \[ \widetilde{r}(v,\Delta):=r\circ e(v,\Delta)=r(e^{-\Delta v}). \]

\(v\in\mathbb{R}_+^N\) が大雑把に方向を意味し,\(\Delta>0\)\(0\) に近づくほど荷重が一様に近づくようなつまみの役割を果たす.

\(\Delta\) をタイムステップとも考えると,粒子フィルターに \(\Delta\searrow0\) の極限が存在するためには,リサンプリング回数が爆発しないために, \[ \widetilde{r}(v,\Delta)(a)=O(\Delta)\quad(v\in\mathbb{R}_+^N,a\ne\mathrm{id}_{[N]}) \tag{1}\] が成り立つ必要がある.

1.2.2 リサンプリング安定性の定式化

条件 (1) が成り立つことは, \[ \iota(v,\Delta):=\frac{r(e^{-\Delta v})}{\Delta}\quad\in\mathcal{M}^1([N]^{[N]}). \] により定まる関数 \[ \iota:\mathbb{R}_+^N\times(0,1)\to\mathcal{M}^1([N]^{[N]}) \] が,任意の \(a\in[N]^{[N]}\setminus\{\mathrm{id}_{[N]}\}\) に関して \(\Delta\searrow0\) に関する極限 \(\lim_{\Delta\searrow0}\iota(v,\Delta)\) を持つことに同値.

ここで,このリサンプリング強度関数 \[ \iota:\mathbb{R}^N_+\times(0,1)\to\mathcal{M}^1([N]^{[N]}) \] は連続になるとする.これはリサンプリング法 \(r:\mathbb{R}_+^N\to\mathcal{P}([N]^{[N]})\) が連続ならば成り立つ. 2

すると結局,条件 (1) は,\(\iota\) が連続な延長 \[ \overline{\iota}:\mathbb{R}_+^N\times[0,1)\to\mathcal{M}^1([N]^{[N]}\setminus\{\mathrm{id}_{[N]}\}) \] を持てば十分.3

これは \(r\) が連続であることに加えて,上の収束が \(v\in\mathbb{R}_+^N\) 上一様に成り立つならば,成り立つ.

\(\overline{\iota}(v):=\overline{\iota}(v,0)\)極限リサンプリング強度 であり,総和を取ったもの \[ \begin{align*} \overline{\iota}^*(w)&:=(\overline{\iota}(w,0)|1)\\ &=\sum_{a\in[N]^{[N]}\setminus\{\mathrm{id}_{[N]}\}}\overline{\iota}(w,0)(a)\quad\in\mathbb{R}_+ \end{align*} \]全リサンプリング率 と呼べる.4

補題(ここまでの議論のまとめ)
  1. リサンプリング法 \(r:\mathbb{R}_+^N\to\mathcal{P}([N]^{[N]})\) が連続で,
  2. \([N]^{[N]}\setminus\{\mathrm{id}_{[N]}\}\) 上の関数の空間 \(\mathcal{M}^1([N]^{[N]}\setminus\{\mathrm{id}_{[N]}\})\) での収束 \[ \overline{\iota}(v)=\lim_{\Delta\searrow0}\frac{r(e^{-\Delta v})}{\Delta} \]\(v\in\mathbb{R}^N_+\) 上一様に成り立つ

とする.このとき,極限リサンプリング強度 \[ \overline{\iota}:\mathbb{R}_+^N\to\mathcal{M}^1([N]^{[N]}\setminus\{\mathrm{id}_{[N]}\}) \]\(\mathcal{M}^1([N]^{[N]}\setminus\{\mathrm{id}_{[N]}\})\) の弱位相について連続である.

\[ f(x,y):=x^{\frac{1}{y}} \] によって \(f:(0,1]\times(0,1)\to[0,1]\) を定める.任意の \(x\in(0,1]\) について,極限 \[ \lim_{y\searrow0}x^{\frac{1}{y}}=\begin{cases} 0&x<1,\\ 1&x=1. \end{cases} \] は存在するが,延長 \[ \overline{f}:(0,1]\times[0,1]\to[0,1] \] は連続ではない.\(\overline{f}(x,0)=\delta_1(x)\) である.

しかし,距離空間 \(X,Y\) に対して,\(f:X\times Y\to\mathbb{R}\) が,境界点 \(\overline{y}\in\partial Y\) に対して,\(f(x,y)\to f(x,\overline{y})\)\(x\in X\) に関して一様に成り立つならば,\(f:X\times\overline{Y}\to\mathbb{R}\) は連続になる.

一様収束するとは,任意の \(\epsilon>0\) に対して,ある \(\delta>0\) が存在して, \[ \lvert y-\overline{y}\rvert<\delta\Rightarrow\|f(-,y)-f(-,\overline{y})\|_\infty<\epsilon \] が成り立つことをいう.このとき,\(f(-,\overline{y}):X\to\mathbb{R}\) は連続になる.5

すると,\(f:X\times\overline{Y}\to\mathbb{R}\) も連続であることが示せる.任意の収束点列 \(\{(x_n,y_n)\}\subset X\times Y\) を取る.収束先が \((x,\overline{y})\in X\times\partial Y\) の形である場合について,\(f(x_n,y_n)\to f(x,\overline{y})\) が成り立つことを示せばよい.

実際このとき,ある \(N\in\mathbb{N}\) が存在して,任意の \(n\ge N\) について \[ \|f(-,y_n)-f(-,\overline{y})\|_\infty<\frac{\epsilon}{2} \] かつ \[ \lvert f(x_n,\overline{y})-f(x,\overline{y})\rvert<\frac{\epsilon}{2} \] であるから, \[ \begin{align*} &\lvert f(x_n,y_n)-f(x,\overline{y})\rvert\\ \le&\lvert f(x_n,y_n)-f(x_n,\overline{y})\rvert+\lvert f(x_n,\overline{y})-f(x,\overline{y})\rvert\\ <&\frac{\epsilon}{2}+\frac{\epsilon}{2}=\epsilon. \end{align*} \]

2 粒子フィルターの構成

目標となる Feynman-Kac 測度

\(\mathbb{R}^d\) 上の伊藤拡散 \(\{z_t\}\subset\mathcal{L}(\Omega';\mathbb{R}^d)\) は, \[ b_i,\sigma_{ij}\in\mathrm{Lip}_b(\mathbb{R}^d)\quad(i,j\in[d]) \] が定める確率微分方程式 \[ z_t=b(z_t)dt+\sigma(z_t)dB_t \] で定まるものとする.

これを参照過程とし,ポテンシャル \(V\in C_b(\mathbb{R}^d)_+\) が定める Feynman-Kac 測度を \(\Pi\in C_c(D_{\mathbb{R}^d}(\mathbb{R}_+))^*\) と表す: \[ \Pi(f):=\frac{1}{\mathcal{Z}}\operatorname{E}\left[f(z_{[0,\tau]})\exp\left(-\int^\tau_0V(z_u)\,du\right)\right] \] \[ \mathcal{Z}:=\operatorname{E}\left[\exp\left(-\int^\tau_0V(z_u)\,du\right)\right],\qquad f\in C_c(D_{\mathbb{R}^d}(\mathbb{R}_+);\mathbb{R}^d). \]

この \(D_{\mathbb{R}^d}(\mathbb{R}_+)\) 上の Feynman-Kac 測度を,離散時間粒子フィルターがどこまで近似できるかを考える.

粒子フィルターの構成

粒子フィルター \(\{X^\Delta_k\}\subset\mathcal{L}(\Omega;M_{Nd}(\mathbb{R}))\) は最も自然に構成する.

粒子数 \(N\) のリサンプリング法 \(r\) に関して,サンプリング間隔を \(\Delta>0\) として, \[ \begin{align*} (X^\Delta_k)^i&:=(X_{k-1}^\Delta)^{A_k(i)}+\Delta\cdot b\left((X_{k-1}^\Delta)^{A_k(i)}\right)\\ &\qquad+\sigma\left((X_{k-1}^\Delta)^{A_k(i)}\right)(B_{k\Delta}^i-B_{(k-1)\Delta}^i)\\ &\qquad\qquad(i\in[N],k\in\mathbb{N})\\ A_k&\overset{\text{iid}}{\sim}r(e^{-\Delta V(X^\Delta_k)}) \end{align*} \] と再帰的に定める.

これを \(D_{\mathbb{R}^d}(\mathbb{R}_+)\)-過程と見たものを \[ Z_t^\Delta:=X_{\left\lfloor\frac{t}{\Delta}\right\rfloor}^\Delta \] と表す.フィルトレーション \(\mathcal{F}_t\)\((B_s)_{s\in[0,t]}\)\(\left(A_{\left\lfloor\frac{s}{\Delta}\right\rfloor}\right)_{s\in[0,t]}\) とが生成するものの,完備右連続化とする.

  1. 任意の単調減少列 \(\{\Delta_n\}\subset\mathbb{R}^+\) に対して,粒子フィルターの列 \[ \{X^{\Delta_n}\}\subset\mathcal{L}(\Omega;D_{M_{Nd}(\mathbb{R})}(\mathbb{R}_+)) \] は一様に緊密である.
定理

Feynman-Kac 測度を定める参照過程 \(\{z_t\}\) とポテンシャル \(V:\mathbb{R}^d\to\mathbb{R}_+\) について,次を仮定する:

  1. \[ b_i,\sigma_{ij}\in\mathrm{Lip}_b(\mathbb{R}^d)\quad(i,j\in[d]) \]
  2. \[ \inf_{x\in\mathbb{R}^d}\inf_{\substack{\theta\in\mathbb{R}^d\\\|\theta\|_2=1}}\|\sigma(x)\theta\|_2>0 \]
  3. \[V\in C_b(\mathbb{R}^d)_+\]

リサンプリング法 \(r\) について,極限リサンプリング測度 \[ \overline{\iota}:\mathbb{R}_+^{N}\to\mathcal{M}^1([N]^{[N]}\setminus\{\mathrm{id}_{[N]}\}) \] が存在し,連続であるとする.6

このとき,これが定める粒子フィルター \(\{X_k^\Delta\}_{k=0}^\infty\) の càdlàg 延長 \(\{Z_t^\Delta\}_{t\in\mathbb{R}_+}\) は,任意の単調減少列 \(\{\Delta_n\}\subset\mathbb{R}^+\) に対して,次を満たす:

  1. ある Lévy 過程 \(\{Z_t\}\) に分布収束する: \[ (Z_t^{\Delta_n})\overset{\text{d}}{\to}(Z_t). \]
  2. 分布収束極限の生成作用素 \(\mathcal{L}\) は,伊藤拡散 \(\{z_t\}\) の生成作用素 \[ \begin{align*} Lf(x)&=\frac{1}{2}\sum_{i,j=1}^d(\sigma\sigma^\top)_{ij}(x)\frac{\partial ^2f}{\partial x_i\partial x_j}(x)\\ &\qquad+\sum_{i=1}^db_i(x)\frac{\partial f}{\partial x_i}(x)\\ &\qquad\qquad(f\in C_c^2(\mathbb{R}^d),x\in\mathbb{R}^d) \end{align*} \] を用いて, \[ \begin{align*} \mathcal{L}f(x)&:=\sum_{n=1}^N\sum_{i=1}^db_i(x^n)\frac{\partial f}{\partial x^n_i}(x)\\ &\qquad+\sum_{n=1}^N\frac{1}{2}\sum_{i,j=1}^d(\sigma\sigma^\top)_{ij}(x^n)\frac{\partial ^2f}{\partial x^n_i\partial x^n_j}(x)\\ &\qquad+\sum_{a\ne1:N}\overline{\iota}(V(x),a)\biggr(f(x^{a(1:N)})-f(x^{1:N})\biggl)\\ &\qquad\qquad(f\in C_c^2(\mathbb{R}^{dN}),x\in\mathbb{R}^{dN},x^n\in\mathbb{R}^d) \end{align*} \] と表せる.
  3. \(\mathcal{V}\in C_b(\mathbb{R}_+\times\mathbb{R}^{dN})_+\) を有界連続なポテンシャルとすると,任意の \(f\in C_b(\mathbb{R}^{dN})\) について, \[ \begin{align*} \lim_{n\to\infty}&\operatorname{E}\left[f(X^{\Delta_n}_{\lfloor\tau/\Delta_n\rfloor})\exp\left(-\sum_{k=0}^{\lfloor\tau/\Delta_n\rfloor-1}\Delta_n\mathcal{V}(k\Delta_n,X^{\Delta_n}_k)\right)\right]\\ =&\operatorname{E}\left[f(Z_\tau)\exp\left(-\int^\tau_0\mathcal{V}(u,Z_u)\,du\right)\right]. \end{align*} \]

References

Alexandroff, A. D. (1940). Additive set-functions in abstract spaces. Matematiceskij Sbornik, 50(8), 307–348.
Chopin, N., Singh, S. S., Soto, T., and Vihola, M. (2022). On Resampling Schemes for Particle Filter with Weakly Informative Observations. The Annals of Statistics, 50(6), 3197–3222.
Dudley, R. M. (2002). Real analysis and probability,Vol. 74. Cambridge University Press.
Scheffé, H. (1947). A useful convergence theorem for probability discributions. The Annals of Mathematical Statistics, 18(3), 434–438.
杉浦光夫. (1980). 解析入門I. 東京大学出版会.

Footnotes

  1. (Dudley, 2002, p. 389) 演習11.1.2.↩︎

  2. \(\mathcal{P}([N]^{[N]})\) 上の弱位相は,質量関数の各点収束に同値.よって,任意の附番 \(a\in[N]^{[N]}\) に対して,\(r(e^{-\Delta w})(a)\in[0,1]\) が連続であることに同値.つまり,\(a\in[N]^{[N]}\) の通りにリサンプリングされる確率が,荷重 \(w\in\mathbb{R}_+^N\) の変化に対して連続的に変化することを要請している.↩︎

  3. ただしもちろん,\(r\) が連続,または \(\iota\) が連続であるという仮定の下で.↩︎

  4. ただし,\((\overline{\iota}(w,0)|1)\)\(1\)\([N]^{[N]}\setminus\{\mathrm{id}_{[N]}\}\) 上の定値関数とした.↩︎

  5. (杉浦光夫, 1980, p. 305) 定理13.3 など.↩︎

  6. その結果,任意の \(a\in[N]^{[N]}\setminus\{\mathrm{id}_{[N]}\}\) に関して,\(\mathbb{R}^{dN}\ni x\mapsto\overline{\iota}(V(x^1,\cdots,x^N))(a)\) が有界連続になる.↩︎