在 NumPy 中计算自协方差函数向量而不使用 np.correlate

Posted 2023-02-23

【中文标题】在 NumPy 中计算自协方差函数向量而不使用 np.correlate

【英文标题】：Computing autocovariance function vector in NumPy without using np.correlate

【发布时间】：2019-03-15 16:36:43

【问题描述】：

我正在尝试创建一个程序，该程序使用 Hannan-Rissanen 算法来计算 ARMA(p, q) 自回归移动平均随机过程的样本参数。

我遇到困难的主要步骤是计算时间序列的自协方差函数。

程序应该接受一个 n×1 维的列向量 Y 并计算一个 k×1 维的列向量 γ^hat，由下式给出：

acvf equation image

其中 Ybar 是 Y 元素的平均值。

如何有效地计算上述总和？ （显然 for 循环会起作用，但我正试图在矢量化 numpy 操作方面做得更好）因为我将其用作学习经验，所以我不希望使用任何 numpy 函数，而不是像 np.sum 或np.mean.

已提出以下类似问题，但并未完全回答我的问题：

Computing autocorrelation of vectors with numpy（使用np.correlate）

（其他一些人在使用更高级的numpy 函数时遇到了同样的问题，或者没有像我在这里想要的那样吐出向量。）

【参考方案1】：

这是替换np.correlate 的一种方法（我认为这是主要困难；我还假设您不想手动编写fft 代码）：

def autocorr_direct(a, debug=False):
    n, _ = a.shape
    out = np.zeros((n+1, 2*n-1), a.dtype)
    out.reshape(-1)[:2*n*n].reshape(n, 2*n)[::-1, :n] = a*a.T
    if debug:
        print(out.reshape(-1)[:2*n*n].reshape(n, 2*n))
        print(out)
    return out.sum(0)

例如：

>>> a = np.array([[1, 1, 2, -1]]).T
>>> autocorr_direct(a, True)
[[-1 -1 -2  1  0  0  0  0]
 [ 2  2  4 -2  0  0  0  0]
 [ 1  1  2 -1  0  0  0  0]
 [ 1  1  2 -1  0  0  0  0]]
[[-1 -1 -2  1  0  0  0]
 [ 0  2  2  4 -2  0  0]
 [ 0  0  1  1  2 -1  0]
 [ 0  0  0  1  1  2 -1]
 [ 0  0  0  0  0  0  0]]
array([-1,  1,  1,  7,  1,  1, -1])
>>> np.correlate(a[:, 0], a[:, 0], 'full')
array([-1,  1,  1,  7,  1,  1, -1])

注意剪切方形数组a[::-1]*a.T的重塑技巧。

注2；要从一维向量X 获取列向量，请使用X[:, None]。