Numpy：用索引向量化数学公式

Question

我正在使用以下方程式模拟马尔可夫生与死过程：（抱歉，MathJax 不可用）。

$ b_i = \textbf{P}(X_{n+1} = i +1 | X_n = i)$

$ d_i = \textbf{P}(X_{n+1} = i -1 | X_n = i)$

$\pi_i b_i = \pi_{i+1} d_{i+1}$

$\pi_i = \pi_0 \frac{b_0 b_1 b_2 \dots b_{i-1}}{d_1 d_2 \dots d_i},$ $i = 1, \dots , m$

$\sum \pi_i = 1$

我使用这段代码来计算 pi_i 的：

import numpy as np

m = 10
# b_0 = 0.5, b_1, ..., b_{m-1} = 0.25
b = 0.25 * np.ones(m-1)
b[0] = 0.5

# d1, ..., d_{m-1} = 0.25, d_m = 0.5
d = 0.25 * np.ones(m-1) # d[0] corresponds to d_1 
d[-1] = 0.5

pi = np.ones(m)

for i in range(1, len(pi)):
    #print(f'i={i}, b[0:i] = {b[0:i]}, d[0:i] =  {d[0:i]}')
    pi[i] = pi[0] * np.prod(b[0:i]) / np.prod(d[0:i])

pi_normalized = pi / pi.sum()
print(pi_normalized)

我对矢量化这个 for 循环很感兴趣：

for i in range(1, len(pi)):
    #print(f'i={i}, b[0:i] = {b[0:i]}, d[0:i] =  {d[0:i]}')
    pi[i] = pi[0] * np.prod(b[0:i]) / np.prod(d[0:i])

问题一：如何"vectorize"上面的for循环？ "vectorizing" 我的意思是不使用 numpy.vectorize。相反，我想通过使用数组操作来加速代码，比如 c_i = a_i b_i --> c = a * b where a,b and c are arrays.

问题2：是否有通用的method/algorithm到"vectorize"带索引的数学公式？

Answer 1

您可以为此使用 cumprod，这是可能有效的 sudo 代码

result = pi[0] * np.cumprod(b) / np.cumprod(d)

我认为问题 2 很模糊，这取决于您拥有的数学公式的类型。通常，要向量化涉及索引的数学运算，您可以使用 np.cumsum、np.cumprod、np.arange 等函数，具体取决于可以对整个数组应用的运算。您可以了解如何使用矢量化和广播（矢量化的重要部分）here，在广播和矢量化部分

Numpy：用索引向量化数学公式

Numpy: Vectorize mathematical formula with indices

python

numpy

vectorization