マルコフ連鎖解析と定常分布

この例では、固有値分解を計算して普通のマルコフ連鎖からシンボリック定常分布を導出する方法を示します。

定常分布は、ステップ数または遷移が増加する場合の、制限的で時間非依存のマルコフ過程の分布を表します。

前のイメージで示すように、状態 A から F の間の (正の) 遷移確率を定義します。

syms a b c d e f cCA cCB positive;

さらに、遷移確率を抑制する仮定を追加します。これは後で適切な定常分布を選択する際に役立ちます。

assumeAlso([a, b, c, e, f, cCA, cCB] < 1 & d == 1);

遷移行列を定義します。状態 A から F が列と行 1 から 6 にマッピングされます。各行の値の和は 1 になります。

P = sym(zeros(6,6));
P(1,1:2) = [a 1-a];
P(2,1:2) = [1-b b];
P(3,1:4) = [cCA cCB c (1-cCA-cCB-c)];
P(4,4) = d;
P(5,5:6) = [e 1-e];
P(6,5:6) = [1-f f];
P

P = 
$$\left(\begin{array}{cccccc}a& 1-a& 0& 0& 0& 0\\ 1-b& b& 0& 0& 0& 0\\ \mathrm{cCA}& \mathrm{cCB}& c& 1-\mathrm{cCA}-\mathrm{cCB}-c& 0& 0\\ 0& 0& 0& d& 0& 0\\ 0& 0& 0& 0& e& 1-e\\ 0& 0& 0& 0& 1-f& f\end{array}\right)$$

マルコフ連鎖の各状態について可能なすべての解析的定常分布を計算します。これはeigを抽出する問題で、遷移確率の値によっては、対応する固有値が 1 に等しくなる可能性があります。

[V,D] = eig(P');

解析的な固有ベクトル

V

V = 
$$\begin{array}{l}\left(\begin{array}{ccccccccc}\frac{b-1}{a-1}& 0& 0& -\frac{\left(c-d\right)\hspace{0.17em}\left(\mathrm{cCB}-b\hspace{0.17em}\mathrm{cCA}-b\hspace{0.17em}\mathrm{cCB}+c\hspace{0.17em}\mathrm{cCA}\right)}{{\sigma }_1}& \frac{b-1}{a-d}& 0& 0& -1& 0\\ 1& 0& 0& -\frac{\left(c-d\right)\hspace{0.17em}\left(\mathrm{cCA}-a\hspace{0.17em}\mathrm{cCA}-a\hspace{0.17em}\mathrm{cCB}+c\hspace{0.17em}\mathrm{cCB}\right)}{{\sigma }_1}& 1& 0& 0& 1& 0\\ 0& 0& 0& -\frac{c-d}{c+\mathrm{cCA}+\mathrm{cCB}-1}& 0& 0& 0& 0& 0\\ 0& 1& 0& 1& 0& 1& 0& 0& 0\\ 0& 0& \frac{f-1}{e-1}& 0& 0& 0& -\frac{f-1}{d-e}& 0& -1\\ 0& 1-d& 1& 0& 0& 0& 1& 0& 1\end{array}\right)\\ \\ \mathrm{where}\\ \\ \mathrm{ }{\sigma }_1=\left(c+\mathrm{cCA}+\mathrm{cCB}-1\right)\hspace{0.17em}\left(a+b-a\hspace{0.17em}c-b\hspace{0.17em}c+c^2-1\right)\end{array}$$

解析的な固有値

diag(D)

ans = 
$$\left(\begin{array}{c}1\\ 1\\ c\\ d\\ a+b-1\\ e+f-1\end{array}\right)$$

厳密に 1 に等しい固有値を求めます。固有値の条件の決定にあいまいさがある場合は、エラーで停止します。このように、このステップが正常に終了すれば、以下のインデックスのリストは信頼できるものです。

ix = find(isAlways(diag(D) == 1,'Unknown','error'));
diag(D(ix,ix))

ans = 
$$\left(\begin{array}{c}1\\ 1\\ d\end{array}\right)$$

解析的定常分布を抽出します。固有ベクトルは 1 ノルムまたは sum(abs(X)) で正規化されてから表示されます。.

for k = ix'
    V(:,k) = simplify(V(:,k)/norm(V(:,k)),1);
end
Probability = V(:,ix)

Probability = 
$$\begin{array}{l}\left(\begin{array}{ccc}\frac{b-1}{\left(a-1\right)\hspace{0.17em}{\sigma }_1}& 0& \frac{{\sigma }_5}{{\sigma }_2}\\ \frac{1}{{\sigma }_1}& 0& \frac{{\sigma }_6}{{\sigma }_2}\\ 0& 0& -\frac{c-1}{{\sigma }_3\hspace{0.17em}\left(c+\mathrm{cCA}+\mathrm{cCB}-1\right)}\\ 0& 1& \frac{1}{{\sigma }_3}\\ 0& 0& 0\\ 0& 0& 0\end{array}\right)\\ \\ \mathrm{where}\\ \\ \mathrm{ }{\sigma }_1=\sqrt{\frac{{\left(b-1\right)}^2}{{\left(a-1\right)}^2}+1}\\ \\ \mathrm{ }{\sigma }_2={\sigma }_3\hspace{0.17em}\left(c+\mathrm{cCA}+\mathrm{cCB}-1\right)\hspace{0.17em}\left(a+b-c-1\right)\\ \\ \mathrm{ }{\sigma }_3=\sqrt{\frac{{\left(c-1\right)}^2}{{\left(c+\mathrm{cCA}+\mathrm{cCB}-1\right)}^2}+\frac{{{\sigma }_6}^2}{{\sigma }_4}+\frac{{{\sigma }_5}^2}{{\sigma }_4}+1}\\ \\ \mathrm{ }{\sigma }_4={\left(c+\mathrm{cCA}+\mathrm{cCB}-1\right)}^2\hspace{0.17em}{\left(a+b-c-1\right)}^2\\ \\ \mathrm{ }{\sigma }_5=\mathrm{cCB}-b\hspace{0.17em}\mathrm{cCA}-b\hspace{0.17em}\mathrm{cCB}+c\hspace{0.17em}\mathrm{cCA}\\ \\ \mathrm{ }{\sigma }_6=\mathrm{cCA}-a\hspace{0.17em}\mathrm{cCA}-a\hspace{0.17em}\mathrm{cCB}+c\hspace{0.17em}\mathrm{cCB}\end{array}$$

最初の固有ベクトル ケースでは、A または B の定常状態の確率は遷移確率 a と b の関数です。この依存関係を可視化します。

fsurf(Probability(1), [0 1 0 1]);
xlabel a
ylabel b
title('Probability of A');

figure(2);
fsurf(Probability(2), [0 1 0 1]);
xlabel a
ylabel b
title('Probability of B');

定常分布は以下を裏付けています (状態 A から F は行インデックス 1 から 6 に対応することを思い出してください)。