动态时间扭曲算法提供了两个时间序列之间的距离概念,这两个时间序列的速度可能不同。如果我有N个序列相互比较,我可以通过两两应用算法构造一个NXN对称矩阵和一个nule对角。然而,对于长的二维序列,这是非常缓慢的。因此,我试图对代码进行矢量化,以加速矩阵计算。重要的是,我还想提取定义最佳对齐的索引
到目前为止,我的成对比较代码:
import math
import numpy as np
seq1 = np.random.randint(100, size=(100, 2)) #Two dim sequences
seq2 = np.random.randint(100, size=(100, 2))
def seqdist(seq1, seq2): # dynamic time warping function
ns = len(seq1)
nt = len(seq2)
D = np.zeros((ns+1, nt+1))+math.inf
D[0, 0] = 0
cost = np.zeros((ns,nt))
for i in range(ns):
for j in range(nt):
cost[i,j] = np.linalg.norm(seq1[i,:]-seq2[j,:])
D[i+1, j+1] = cost[i,j]+min([D[i, j+1], D[i+1, j], D[i, j]])
d = D[ns,nt] # distance
matchidx = [[ns-1, nt-1]] # backwards optimal alignment computation
i = ns
j = nt
for k in range(ns+nt+2):
idx = np.argmin([D[i-1, j], D[i, j-1], D[i-1, j-1]])
if idx == 0 and i > 1 and j > 0:
matchidx.append([i-2, j-1])
i -= 1
elif idx == 1 and i > 0 and j > 1:
matchidx.append([i-1, j-2])
j -= 1
elif idx == 2 and i > 1 and j > 1:
matchidx.append([i-2, j-2])
i -= 1
j -= 1
else:
break
matchidx.reverse()
return d, matchidx
[d,matchidx] = seqdist(seq1,seq2) #try it
这里是对代码的一次重写,使其更易于
numba.jit
。这并不是一个完全矢量化的解决方案,但是我看到这个基准测试的加速因子是230以下是一些变化:
cost
使用spatial.distance_matrix
计算李>idx
的定义被一堆丑陋的if语句所取代,这使得编译代码更快李>min([D[i, j+1], D[i+1, j], D[i, j]])
替换为min(D[i, j+1], D[i+1, j], D[i, j])
,也就是说,我们不取列表的min,而是取三个值的min。这导致了jit
下的惊人加速李>matchidx
作为numpy数组预先分配,并在输出之前被截断为正确的大小李>相关问题 更多 >
编程相关推荐