时间戳

首页 > 解决方案 > 概率矩阵 Numpy 的终端概率

python - 概率矩阵 Numpy 的终端概率

问题描述

我有一个矩阵 m 表示从状态转换到状态的概率。

例如，对于下面的示例，我将始终卡在状态 1、3、4 和状态 2 中，我将随机转换到 4 个状态之一。

import numpy as np
m = np.eye(4)
m[1] = 0.25
print(m)

[[1.   0.   0.   0.  ]

 [0.25 0.25 0.25 0.25]

 [0.   0.   1.   0.  ]

 [0.   0.   0.   1.  ]]

如何找到表示无限转换后最终状态的矩阵？

例如，如果我这样做，我会得到状态 1,3,4 --> 100% 坚持在 1,3,4 但状态 2 --> 1/3 机会结束于所有其他状态的直观结果。由于来自状态 2 的所有案例最终通过多次转换在 1、3、4 之间平均分配。

t = m
for _ in range(100_000):
    t = t @ t
print(t)

[[1.         0.         0.         0.        ]

 [0.33333333 0.         0.33333333 0.33333333]

 [0.         0.         1.         0.        ]

 [0.         0.         0.         1.        ]]

如何在不使用重复乘法的情况下计算这个？我认为它对应于矩阵的特征向量/特征值，但是当我计算它时我得到了一些非常不同的东西。

np.linalg.eig(m)

[[0.        , 0.9486833 , 0.        , 0.        ],

[1.        , 0.31622777, 0.31622777, 0.31622777],

[0.        , 0.        , 0.9486833 , 0.        ],

[0.        , 0.        , 0.        , 0.9486833 ]]

有没有一种方法可以使用 numpy 进行计算？我需要它来处理任意矩阵，但是会有一个已知的终端状态列表和从所有其他状态到达这些状态的正概率。

目前我正在考虑使用重复乘法方法，但感觉不是最理想的，应该有一个函数/技巧可以在不循环的情况下计算。

我正在阅读这篇文章，但没有完全理解该方法是什么以及如何实施它。

https://math.dartmouth.edu/archive/m20x06/public_html/Lecture14.pdf

我也看过这个问题。人们似乎给出了一些手工求解的技巧，但没有给出通用算法：

https://math.stackexchange.com/questions/2003258/calculating-the-probability-of-reaching-each-charging-state-in-markov-chain

标签： pythonnumpymatrixprobabilitymarkov-chains

解决方案

我的朋友指出了以下技巧。

特征分解意味着我们可以将原始矩阵写为

V x D x V^-1

其中 D 是具有特征值的对角矩阵，V 是特征向量。

如果我们将它自己无限次相乘，它就是

V x D^inf x V^-1

我们可以使用下面的 numpy 计算。

d, v = np.linalg.eig(m)
v @ np.diag(d >= 1).astype(int) @ np.linalg.inv(v)

因为如果对角线值 < 1，那么当我们相乘时，它们将趋于 0（假设我们有一个具有有效概率的矩阵，并且所有状态都可以达到终端状态）。

推荐阅读