是否可以使用矢量化在 numpy 中加速这种平均自相关计算？

Question

我希望加快以下在 numpy 中使用标准自相关函数的平均自相关代码。是否可以对其进行矢量化？

输入：X -> shape(100,3000) = 3000 项的 100 天时间序列

输出：float -> 找到每个项目的 lag5 自相关的平均值（跨越 3000 个项目）

import numpy as np
import time
A = np.random.rand(100,3000)
start_time = time.time()
np.nanmean(np.array([np.corrcoef(np.array([A[:-5,i],A[5:,i]]))[1,0] for i in range(A.shape[1])]),axis=0)
print("--- %s seconds ---" % (time.time() - start_time))

--- 0.35528588 秒 ---

Answer 1

我假设您主要关心更小的运行时间，而不是如何实现它。以下方法通过改进相关系数的计算实现了超过 3 倍的加速。

更具体地说，它将对 np.corrcoef 的调用替换为对自定义函数 fast_corrcoef 的调用。本质上，fast_corrcoef 遵循 numpy 计算相关系数的方法，只是它不执行大量冗余计算。主要是，它利用了这样一个事实，即对于 np.corrcoef 的每次调用，调用 np.cov 来计算协方差矩阵，这需要计算两个序列 S[:-lag] 和 S[lag:] 的平均值，做了很多不必要的计算：不是像 np.cov 那样单独计算平均值，fast_corrcoef 通过首先计算 S[lag:-lag] 的总和然后计算两个系列的平均值来计算它们中间金额对于大小为 n 的时间序列，这仅产生 (n - 2 * lag) + 2 * lag + 2 = n + 2 加法，而不是 2n 加法，即本质上节省了所有加法的 50% 用于计算平均值，如果 A.shape = (100, 3000) 和 lag = 5 总计 (2 * (100 - 5) - (100 - 5 + 2)) * 3000 = 279000 保存的添加。

此外，由于协方差矩阵是对称的，因此不需要在计算时进行全矩阵乘法运算，因此对 np.dot 的调用也已替换为表达式 np.cov它只计算实际需要的树向量点积，即只计算协方差一次而不是两次。

import numpy as np
import time


def fast_corrcoef(X, lag):
    n = X.shape[1]
    mid_sum = X[0, lag:].sum()
    X -= np.array([(mid_sum + X[0, :lag].sum()) / n, (mid_sum + X[1, -lag:].sum()) / n])[:, None]
    return (X[0,:] * X[1,:]).sum() / (np.sqrt((X[0,:] ** 2).sum()) * np.sqrt((X[1,:] ** 2).sum()))
    
A = np.random.rand(100,3000)

lag = 5

def orig(A):
    return np.nanmean(np.array([np.corrcoef(np.array([A[:-5,i],A[5:,i]]))[1,0] for i in range(A.shape[1])]),axis=0)

def faster(A, lag):
    corrcoefs = np.empty([A.shape[1]])
    for i in range(A.shape[1]):
        corrcoefs[i] = fast_corrcoef(np.array([A[:-lag, i], A[lag:, i]]), lag)
    return np.nanmean(corrcoefs, axis=0)

start_time = time.time()
res_orig = orig(A)
runtime_orig = time.time() - start_time
print(f'orig: runtime: {runtime_orig:.4f}s; result: {res_orig}')


start_time = time.time()
res_faster = faster(A, lag)
runtime_faster = time.time() - start_time
print(f'fast: runtime: {runtime_faster:.4f}s; result: {res_faster}')

assert abs(res_orig - res_faster) < 1e-16

speedup = runtime_orig / runtime_faster

print(f'speedup: {speedup:.4f}x')

此代码打印了以下内容（当我在 Google Colab 中运行它时）：

orig: runtime: 0.4214s; result: -0.00961872402216293
fast: runtime: 0.1138s; result: -0.009618724022162932
speedup: 3.7030x

是否可以使用矢量化在 numpy 中加速这种平均自相关计算？

Is it possible to speed up this mean autocorrelation calculation in numpy using vectorization?

python

numpy

vectorization

autocorrelation