【问题标题】:Speed up distance calculations, sliding window加速距离计算,滑动窗口
【发布时间】:2018-06-05 22:47:00
【问题描述】:

我有两个时间序列 A 和 B。长度为 m 的 A 和长度为 n 的 B。 m << n。两者都有维度d

我通过在 B 上滑动 A 来计算 A 和 B 中所有子序列之间的距离。 在 python 中,代码如下所示。

def sliding_dist(A,B)
    n = len(B)
    dist = np.zeros(n)
    for i in range(n-m):
        subrange = B[i:i+m,:]
        distance = np.linalg.norm(A-subrange)
        dist[i] = distance
    return dist

现在这段代码需要很长时间才能执行,而且我有很多计算要做。 我需要加快计算速度。我的猜测是我可以通过使用卷积和频域乘法(FFT)来做到这一点。但是,我一直无法实现。

有什么想法吗? :) 谢谢

【问题讨论】:

  • 这个问题是关于优化代码的。这些问题属于这里:codereview.stackexchange.com
  • 您确实可以使用 FFT 通过频域加速卷积。但是,这里没有卷积。
  • 嗯,不是在当前的实现中。难道不能将问题重写为等效的卷积吗?
  • 您是否尝试过创建一个超级A,它只是长度为n 的A 副本的串联,这样您只需执行一次减法和一次调用norm()?不确定到底是什么瓶颈。你能指定nm的数量级吗?
  • m ~ 10, n ~ 500

标签: python algorithm fft sliding-window


【解决方案1】:

norm(A - subrange) 本身不是卷积,但可以表示为:

sqrt(dot(A, A) + dot(subrange, subrange) - 2 * dot(A, subrange))

如何快速计算每一项:

  • dot(A, A) - 这只是一个常数。

  • dot(subrange, subrange) - 这可以使用递归方法在 O(1)(每个位置)内计算出来。

  • dot(A, subrange) - 这在此上下文中的卷积。所以这可以通过convolution theorem在频域中计算出来。1

但请注意,如果子范围大小仅为 10,您不太可能看到性能提升。


1。又名fast convolution

【讨论】:

  • 用 dot() 计算 norm() 的想法可能会有所改进。
【解决方案2】:

使用矩阵运算实现,就像我在评论中提到的那样。想法是逐步评估规范。在你的情况下,我的价值是:

d[i] = sqrt((A[0] - B[i])^2 + (A[1] - B[+1])^2 + ... + (A[m-1] - B[i+m-1])^2)

前三行计算平方和,最后一行做sqrt()。

加速约 60 倍。

import numpy
import time

def sliding_dist(A, B):
    m = len(A)
    n = len(B)
    dist = numpy.zeros(n-m)
    for i in range(n-m):
        subrange = B[i:i+m]
        distance = numpy.linalg.norm(A-subrange)
        dist[i] = distance
    return dist

def sd_2(A, B):
    m = len(A)
    dist = numpy.square(A[0] - B[:-m])
    for i in range(1, m):
        dist += numpy.square(A[i] - B[i:-m+i])
    return numpy.sqrt(dist, out=dist)

A = numpy.random.rand(10)
B = numpy.random.rand(500)
x = 1000
t = time.time()
for _ in range(x):
    d1 = sliding_dist(A, B)
t1 = time.time()
for _ in range(x):
    d2 = sd_2(A, B)
t2 = time.time()

print numpy.allclose(d1, d2)
print 'Orig %0.3f ms, second approach %0.3f ms' % ((t1 - t) * 1000., (t2 - t1) * 1000.)
print 'Speedup ', (t1 - t) / (t2 - t1)

更新

这是您在矩阵运算中需要的规范的“重新实现”。如果您想要numpy 提供的其他规范,则它不灵活。由于 norm() 接收参数轴,因此可以使用不同的方法来创建 B 个滑动窗口的矩阵并对整个数组进行规范化。这是该方法的实现,但速度提升约为 40 倍,比以前慢。

def sd_3(A, B):
    m = len(A)
    n = len(B)
    bb = numpy.empty((len(B) - m, m))
    for i in range(m):
        bb[:, i] = B[i:-m+i]
    return numpy.linalg.norm(A - bb, axis=1)

【讨论】:

  • 这太棒了!你能解释一下为什么会这样吗? :) 我没有完全关注。
  • @CarlRynegardh 我添加了一些 cmets。我希望这就足够了:-)
猜你喜欢
  • 2017-03-12
  • 2018-08-27
  • 2017-11-15
  • 2013-07-08
  • 1970-01-01
  • 2015-12-16
  • 1970-01-01
  • 2016-02-13
  • 2011-01-10
相关资源
最近更新 更多