用于动态时间规整的矢量化循环

Question

我正在尝试通过矢量化来加速 https://github.com/pierre-rouanet/dtw/blob/master/dtw.py 处的 dtw(x, y, dist=lambda x, y: norm(x - y, ord=1))) 。第一个循环很简单，但我不知道如何向量化第二个循环：

for i in range(r):
    for j in range(c):
        D[i+1, j+1] += min(D[i, j], D[i, j+1], D[i+1, j])

主要问题是，即使我迭代 i，每个 D[i+1, j+1] 都依赖于 D[i+1, j]。

是否可以对其进行矢量化，还是我必须使用 Cython？

对于形状为1000x2的x和y，原来的代码需要15s，我现在的代码需要1.8s，主要是在第二个循环中。

编辑：最小工作示例

np.random.seed(0); A = np.random.randn(4, 3)
r, c = np.array(A.shape)-1
for i in range(r):
    for j in range(c):
        A[i+1, j+1] += min(A[i, j], A[i, j+1], A[i+1, j])
A

应该给：

array([[ 1.76405235,  0.40015721,  0.97873798],
       [ 2.2408932 ,  2.2677152 , -0.57712067],
       [ 0.95008842,  0.79873121, -0.68033952],
       [ 0.4105985 ,  0.55464207,  0.77393398]])

Answer 1

尝试以下操作：

D[1:,1:] += np.min(np.dstack((D[:-1,:-1],D[:-1,1:],D[1:,:-1])),axis=2)

Answer 2

我终于做到了！解决方案是遍历对角线。指数很难得到正确的。谢谢大家！

r, c = np.array(D.shape)-1
for a in range(1, r+c):
    # We have I>=0, I<r, J>0, J<c and J-I+1=a
    I = np.arange(max(0, a-c), min(r, a))
    J = I[::-1] + a - min(r, a) - max(0, a-c)
    # We have to use two np.minimum because np.minimum takes only two args.
    D[I+1, J+1] += np.minimum(np.minimum(D[I, J], D[I, J+1]), D[I+1, J])

用于动态时间规整的矢量化循环

Vectorized loop for dynamic time wrapping

python

numpy