多次元Ornstein-Uhlenbeck過程の確率密度関数の導出メモ
はじめに
この記事では，多次元Ornstein-Uhlenbeck過程と呼ばれる確率過程について，ある初期値からこの確率過程に従って時間発展する状態変数の確率密度関数を導出します．
多次元Ornstein-Uhlenbeck過程
Ornstein-Uhlenbeck (OU) 過程は，ブラウン運動する粒子の振る舞いを記述するために提案された確率過程であり，線形な確率微分方程式によって記述されます*1．多次元の状態変数を持つOU過程を，ここでは多次元OU過程と呼び，以下の (伊藤型) 確率微分方程式で定義します．
但し，はこの系の状態を表す状態変数であり，点は時間
による微分を表します．
は系のダイナミクスを特徴付ける行列であり，異なる
個の固有値を持ち，全ての固有値の実部は負であるとします．
は
\begin{align}
E[\boldsymbol{w}(t)] = 0, \ E[\boldsymbol{w}(t)\boldsymbol{w}(s)^{\top}] = W\delta(t-s)
\end{align}
を満たす白色ガウスノイズとします (はDiracのデルタ関数です)．
状態変数の確率密度関数
状態変数がある初期値
から式(1)に従って時間発展するとします．このとき，時刻
の状態変数
は，ガウス分布の再生性より多次元ガウス分布に従います．この多次元ガウス分布の平均
，共分散行列
を以下のように定義します．
\begin{align}
\boldsymbol{m}(t) &= E[\boldsymbol{x}(t)] \\
V(t) &= E[(\boldsymbol{x}(t) - \boldsymbol{m}(t)) (\boldsymbol{x}(t) - \boldsymbol{m}(t))^{\top}]
\end{align}
このとき，と
の時間発展は以下のように記述されます．
\begin{align}
\boldsymbol{m}(t) &= e^{tA} \boldsymbol{x}(0) \tag{2} \\
V(t) &= {\rm vec}^{-1}\left( \int_0^t e^{s(A \oplus A)} {\rm vec}(W) ds \right) \tag{3}
\end{align}
但し，はある行列
について
\begin{align}
{\rm vec}(M) = [M_{1,1}, \ldots, M_{n,1}, M_{1,2}, \ldots, M_{n,2}, \ldots, M_{1,n}, \ldots, M_{n,n}]^{\top}
\end{align}
と定義される関数であり，はKronecker和*2を表します．また特に，
のとき，
\begin{align}
\boldsymbol{m}(t) &\to \boldsymbol{0} \\
V(t) &\to -{\rm vec}^{-1}\left( (A \oplus A)^{-1} {\rm vec}(W) \right)
\end{align}
となります．
式(2), (3)の導出
以下のように，と
の時間発展を記述する方程式を得ることができます．なお，式(5)の導出には伊藤の補題*3を利用します．
\begin{align}
\dot{\boldsymbol{m}}(t) &= \frac{d}{dt} E[\boldsymbol{x}(t)] = A \boldsymbol{m}(t) \tag{4} \\
\dot{V}(t) &= \frac{d}{dt} E[(\boldsymbol{x}(t) - \boldsymbol{m}(t)) (\boldsymbol{x}(t) - \boldsymbol{m}(t))^{\top}] = AV(t) + V(t)A^{\top} + W \tag{5}
\end{align}
式(4)を初期条件の下で解くことにより，式(2)を導くことができます．また，式(5)はKronecker和を用いることにより，以下のように書き直すことができます．
\begin{align}
\frac{d}{dt} {\rm vec}(V(t)) = (A \oplus A) {\rm vec}(V(t)) + {\rm vec}(W)
\end{align}
この式を初期条件の下で解くことにより，式(3)を導くことができます．
