マルコフ連鎖 (Markov chain)

マルコフ過程のなかで、取りうる状態が離散的（有限的）であるものを特にマルコフ連鎖と言います。前回、お話したランダムウォークはマルコフ連鎖であったということです。もし、粒子の動き$X_n$が$1$か$-1$ではなく、$X_n \in [-1, 1]$のように連続的に移動することが可能であれば、$S_n \in \mathbb{R}$となり、マルコフ連鎖ではなくなります。 この場合をブラウン運動ということも前回お話しました。なので、ランダムウォークもブラウン運動もマルコフ過程ですが、マルコフ連鎖なのはランダムウォークだけです。

さて、まずは天気の状態を晴と雨の2つの状態で定義します。このことを $$ S=\{晴, 雨\} $$ と表します。$S$のような、状態の集合を状態空間といいます。

次に状態の変化を考えます。ある日の天気を考慮した次の日の天気の確率が次のようになるとします。 $$ \begin{array}{c|c||c} \hline ある日(t=n) & 次の日(t=n+1)&確率\\ \hline 晴&晴&0.8\\ 晴&雨&0.2\\ 雨&晴&0.4\\ 雨&雨&0.6\\ \hline \end{array} $$ これを状態遷移図で表すと、次図のようになります。

時刻$t$のとき、晴である確率を$q_1^{(t)}$, 雨である確率を$q_2^{(t)}$とし、確率行ベクトルを $$ \boldsymbol{\pi}^{(t)}:=\begin{bmatrix} q_1^{(t)},\ q_2^{(t)} \end{bmatrix} $$ とします（分かりやすいようにカンマで区切りました）。条件として \begin{align*} q_i \geq 0\ (i&=1,2)\\ q_1^{(t)}+q_2^{(t)}&=1 \end{align*} が成り立ちます。1番目の式は確率は0以上の値を取るということで、2番目の式は晴と雨以外の状態を取らないということを表します。 このとき、$t=n+1$のときの確率行ベクトルは、高校で習う漸化式と同じように考えれば、 \begin{align*} \boldsymbol{\pi}^{(n+1)}= \begin{bmatrix} q_1^{(n+1)},\ q_2^{(n+1)} \end{bmatrix} &=\begin{bmatrix} 0.8q_1^{(n)}+0.4q_2^{(n)},\ 0.2q_1^{(n)}+0.6q_2^{(n)} \end{bmatrix}\\ &=\begin{bmatrix} q_1^{(n)},\ q_2^{(n)} \end{bmatrix} \underbrace{\begin{bmatrix} 0.8&0.2\\ 0.4&0.6 \end{bmatrix}}_{:=P}\\ &=\boldsymbol{\pi}^{(n)}P \end{align*} となります。この$P$を遷移行列(transition matrix)といいます。$P$は同じ行の要素をすべて足し合わせると$1$になります。つまり、$P=[p_{ij}]$とすると、 $$ \sum_j p_{ij}=1 $$ が成り立ちます。要は、次の状態は状態空間内のどれかになるということです。なお、マルコフ連鎖はこの遷移行列か状態遷移図のいずれかで表されます。遷移行列はグラフを描くときに用いる隣接行列と同じです。

少し具体的な例を見てみましょう。1日目($t=0$のとき)晴である、つまり$\boldsymbol{\pi}^{(0)}=\begin{bmatrix} 1&0 \end{bmatrix}$であるとします。この$\boldsymbol{\pi}^{(0)}$を初期分布と呼びます。2日目($t=1$のとき)の天気は、 $$ \boldsymbol{\pi}^{(1)}=\boldsymbol{\pi}^{(0)}P =\begin{bmatrix} 1&0 \end{bmatrix} \begin{bmatrix} 0.8&0.2\\ 0.4&0.6 \end{bmatrix} =\begin{bmatrix} 0.8&0.2 \end{bmatrix} $$ と予測できます。よって2日目は晴の確率が$0.8$で、雨の確率が$0.2$です。3日目($t=2$のとき)も同様に、 $$ \boldsymbol{\pi}^{(2)}=\boldsymbol{\pi}^{(1)}P =\begin{bmatrix} 0.8&0.2 \end{bmatrix} \begin{bmatrix} 0.8&0.2\\ 0.4&0.6 \end{bmatrix} =\begin{bmatrix} 0.72&0.28 \end{bmatrix} $$ と計算できます。ところで、上の2つの計算をまとめて、 \begin{align*} \boldsymbol{\pi}^{(2)}&=\boldsymbol{\pi}^{(1)}P\\ &=\left(\boldsymbol{\pi}^{(0)}P\right)P\\ &=\boldsymbol{\pi}^{(0)}P^2\\ &=\begin{bmatrix} 1&0 \end{bmatrix} \begin{bmatrix} 0.8&0.2\\ 0.4&0.6 \end{bmatrix}^2 =\begin{bmatrix} 0.72&0.28 \end{bmatrix} \end{align*} と一気に計算することもできます。これを一般的に考えれば、 \begin{align*} \boldsymbol{\pi}^{(n)}&=\boldsymbol{\pi}^{(n-1)}P\\ &=\left(\boldsymbol{\pi}^{(n-2)}P\right)P\\ &=\cdots\\ &=\boldsymbol{\pi}^{(0)}P^n \end{align*} が成り立つことが分かります。

この $$\color{red}{\boldsymbol{\pi}^{(n)}=\boldsymbol{\pi}^{(0)}P^n}$$ は、$n$回目の遷移確率を求めるには、遷移行列を$n$乗すれば良いということを表しています。数学的な正確性を求めると、この式はチャップマン・コルモゴロフの等式(Chapman-Kolmogorov equation): $$ p_{ij}^{(t+s)}=\sum _{k\in S}p_{ik}^{(t)}p_{kj}^{(s)} $$ から導かれます($S$は状態空間)。しかし、違和感なくこの事実は理解できると思うので、特にこの式を深く掘り下げることはしません。

マルコフ連鎖は$t \to \infty$にすると、一定の確率分布に収束する場合があります。収束した確率分布を定常分布と呼びます。上記の天気の例を含め多くの場合、マルコフ連鎖は定常分布を持ちます。 もちろん定常分布を持たない場合も存在し、例えば次の2つの図のような場合は収束しません。

それでは実際に定常分布を求めてみましょう。以後の確率過程は全て定常分布が存在する（エルゴード性を持つ）とします。定常分布$\boldsymbol{\pi}$は$t\to\infty$にしたときの分布なので、素直に考えれば $$ \boldsymbol{\pi}=\boldsymbol{\pi}^{(\infty)}=\lim_{n\to\infty} \boldsymbol{\pi}^{(0)}P^n $$ により求められます。$P^n$を計算するには対角化するのが良い方法です。 $ Q:=\begin{bmatrix} 1&1\\ 1&-2 \end{bmatrix} $ とすると、 $$Q^{-1}PQ=\begin{bmatrix} 1&0\\ 0&0.4 \end{bmatrix},\ Q^{-1}=\frac{1}{3}\begin{bmatrix} 2&1\\ 1&-1 \end{bmatrix} $$ と対角化できるので、 \begin{align*} P^n&=Q\left(Q^{-1}PQ\right)^nQ^{-1}\\ &=\frac{1}{3}\begin{bmatrix} 2+(0.4)^n&1-(0.4)^n\\ 2-2(0.4)^n&1+2(0.4)^n \end{bmatrix} \end{align*} となります。ここで、$\boldsymbol{\pi}^{(0)}=\begin{bmatrix} 1&0 \end{bmatrix}$であるとすると、 $$ \lim_{n\to\infty} P^n=\lim_{n\to\infty} \frac{1}{3}\begin{bmatrix} 2+(0.4)^n&1-(0.4)^n\\ 2-2(0.4)^n&1+2(0.4)^n \end{bmatrix}=\frac{1}{3}\begin{bmatrix} 2&1\\ 2&1 \end{bmatrix} $$より、定常分布は \begin{align*} \boldsymbol{\pi}&=\lim_{n\to\infty} \boldsymbol{\pi}^{(0)}P^n\\ &=\begin{bmatrix} 1&0 \end{bmatrix}\cdot\frac{1}{3}\begin{bmatrix} 2&1\\ 2&1\\ \end{bmatrix}\\ &=\begin{bmatrix} \dfrac{2}{3}&\dfrac{1}{3} \end{bmatrix} \end{align*} となります。これは$\boldsymbol{\pi}^{(0)}=\begin{bmatrix} 0&1 \end{bmatrix}$とした場合でも同じです。よって、定常分布では晴れの確率が$\dfrac{2}{3}$, 雨の確率が$\dfrac{1}{3}$になることが分かります。

定常分布が$\boldsymbol{\pi}=\begin{bmatrix} q_1&q_2 \end{bmatrix}$であると仮定すると、 \begin{align*} &\boldsymbol{\pi}=\boldsymbol{\pi} P\\ \iff&\boldsymbol{\pi}\left(I-P\right)=0\ \ (Iは単位行列)\\ \iff&\begin{bmatrix} q_1&q_2 \end{bmatrix}\begin{bmatrix} 0.2&-0.2\\ -0.4&0.4 \end{bmatrix}=0\\ \Longrightarrow\ &0.2q_1-0.4q_2=0\\ \iff& q_1=2q_2 \end{align*} となります。これと$q_1+q_2=1$より、 $$ q_1=\frac{2}{3},\ q_2=\frac{1}{3} $$ となり、前節で求めた答えと一致します。

$X_t　(t=1,\cdots,n)$の同時確率関数は $$ P(X_1=x_1,\cdots,X_n=x_n)=P(X_1=x_1)\prod_{k=1}^{n} P(X_{k+1}=x_{k+1}|X_k=x_k) $$ です。まず$(X_1=x_1)$が起こり、以後の事象$(X_{k+1}=x_{k+1})$の起こる確率は一つ前の事象$(X_k=x_k)$が起こった下での条件付き確率となります。$P(X_1=x_1)$と$P(X_{k+1}=x_{k+1}|X_k=x_k)$の式中のパラメータ$\theta$が未知の場合、上の式は尤度となり、この式を最大とするようなパラメータ$\theta$を求めることが最尤推定となります。

尤度って何、というのに簡単に答えると、尤度とは観測結果$x_1,\cdots,x_n$が仮説の下で観測される確率のことです。確率なのですが、変数$\theta$で積分しても$1$にはなりません。なので、確率だけどちょっと違うということで、「尤度(likelihood)」つまり「観測結果の尤もらしさ」という名前が付けられています。

　この尤度を最大にすることがパラメータ$\theta$を決めることにどうしてなるのかと言うと、