矢量化 numpy 多条件嵌套循环

Question

在尝试在 Python 中由 Felix Scholkmann、Jens Boss 和 Martin Wolf 在 Python 中生成 Noisy 周期和准周期信号中的自动峰值检测时，我在实现中遇到了绊脚石。

在尝试优化时，我注意到嵌套的 for 循环正在造成处理时间瓶颈（平均需要 115394 毫秒才能完成）。有没有更有效的方法来构造嵌套的 for 循环？

注：参数 signal 是算法将处理的坐标列表，其形式为

-48701.0 -20914.0 -1757.0 -49278.0 -106781.0 -88139.0 -13587.0 28071.0 11880.0 -13375.0 -18056.0 -15248.0 -12476.0 -9832.0 -26365.0 -65734.0 -81657.0 -41566.0 6382.0 872.0 -30666.0 -20261.0 17543.0 6278.0 ...

该列表有 32768 行长。

该函数返回在另一个函数中处理的检测到的峰的索引。

def ampd(signal):
    s_time = range(1, len(signal)+1)
    [fitPolynomial, fitError] = np.polyfit(s_time, signal, 1)
    fitSignal = np.polyval([fitPolynomial, fitError], s_time)
    dtrSignal = signal - fitSignal
    N = len(dtrSignal)
    L = math.ceil(N/2.0)-1
    creation_start = time.time()
    np.random.seed(1969)
    LSM = np.random.uniform(0, 2, size=(L, N))
    creation_elapsedTime = time.time() - creation_start
    print('LSM created in %s ms' % int(creation_elapsedTime * 1000))
    loop_start = time.time()
    for k in range(1, L):
        for i in range(k+2, N-k+1):
            if signal[i-1]>signal[i-k-1] and signal[i-1]>signal[i+k-1]:
                    LSM[k,i] = 0
    loop_elapsedTime = time.time() - loop_start
    print('Loop completed in %s ms' % int(loop_elapsedTime * 1000))
    G = np.sum(LSM, axis=1)
    l = min(enumerate(G), key=itemgetter(1))[0]
    MLSM = LSM[0:l]
    S = np.std(MLSM, ddof=1)
    found_indices = np.where(MLSM == ((S-1) == 0))
    del LSM
    del MLSM
    return found_indices[1]

Answer 1

这是一个仅使用一个循环的解决方案

for k in range(1, L):
    mat=1-((signal[k+1:N-k]>signal[1:N-2*k]) & (signal[k+1:N-k]>signal[2*k+1:N]))
    LSM[k,k+2:N-k+1]*=mat

它更快，似乎确实给出了相同的解决方案。您比较切片（如 Ami Tavory 建议的那样），将比较与 & 相结合，给出一个True/False数组;使用1-操作，您可以将其转换为零和一，零对应于满足条件的位置。最后，您将行乘以结果。