如何計算Python中的自協方差

Question

我想計算3個數組X1，X2和Y的自協方差，它們都是平穩的隨機過程。 sciPy或其他庫中是否有任何函數可以解決這個問題？

Answer 1

Statsmodels具有自相關和互協方差函數

http://statsmodels.sourceforge.net/devel/generated/statsmodels.tsa.stattools.acovf.html http://statsmodels.sourceforge.net/devel/generated/statsmodels.tsa.stattools.ccovf.html

加的相關函數和部分自相關 http://statsmodels.sourceforge.net/devel/tsa.html#descriptive-statistics-and-tests

Answer 2

根據自協方差係數爲離散信號的標準估計，這可以由等式表示：

...其中x(i)是一個給定的信號（即特定1D向量），k代表x(i)信號由k樣本移位，N是x(i)信號的長度，並且：

...這是簡單的平均，我們可以這樣寫：

''' 
Calculate the autocovarriance coefficient. 
''' 

import numpy as np 

Xi = np.array([1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5]) 
N = np.size(Xi) 
k = 5 
Xs = np.average(Xi) 

def autocovariance(Xi, N, k, Xs): 
    autoCov = 0 
    for i in np.arange(0, N-k): 
     autoCov += ((Xi[i+k])-Xs)*(Xi[i]-Xs) 
    return (1/(N-1))*autoCov 

print("Autocovariance:", autocovariance(Xi, N, k, Xs)) 

如果你想標準化自協方差係數，這將成爲自相關係數表示爲：

...比你只需要添加到上面的代碼只是兩個額外的線路：

def autocorrelation(): 
    return autocovariance(Xi, N, k, Xs)/autocovariance(Xi, N, 0, Xs) 

這裏完整的腳本：

''' 
Calculate the autocovarriance and autocorrelation coefficients. 
''' 

import numpy as np 

Xi = np.array([1, 2, 3, 4, 5, 1, 2, 3, 4, 5, 1, 2, 3, 4, 5]) 
N = np.size(Xi) 
k = 5 
Xs = np.average(Xi) 

def autocovariance(Xi, N, k, Xs): 
    autoCov = 0 
    for i in np.arange(0, N-k): 
     autoCov += ((Xi[i+k])-Xs)*(Xi[i]-Xs) 
    return (1/(N-1))*autoCov 

def autocorrelation(): 
    return autocovariance(Xi, N, k, Xs)/autocovariance(Xi, N, 0, Xs) 

print("Autocovariance:", autocovariance(Xi, N, k, Xs)) 
print("Autocorrelation:", autocorrelation()) 

Answer 3

獲取樣本自協方差：

# cov_auto_samp(X,delta)/cov_auto_samp(X,0) = auto correlation 
def cov_auto_samp(X,delta): 
    N = len(X) 
    Xs = np.average(X) 
    autoCov = 0.0 
    times = 0.0 
    for i in np.arange(0, N-delta): 
     autoCov += (X[i+delta]-Xs)*(X[i]-Xs) 
     times +=1 
    return autoCov/times 

Answer 4

對以前的答案進行了小小的調整，避免了python for循環，而是使用numpy數組操作。如果你有很多數據，這會更快。

def lagged_auto_cov(Xi,t): 
    """ 
    for series of values x_i, length N, compute empirical auto-cov with lag t 
    defined: 1/(N-1) * \sum_{i=0}^{N-t} (x_i - x_s) * (x_{i+t} - x_s) 
    """ 
    N = len(time_series) 

    # use sample mean estimate from whole series 
    Xs = np.mean(Xi) 

    # construct copies of series shifted relative to each other, 
    # with mean subtracted from values 
    end_padded_series = np.zeros(N+t) 
    end_padded_series[:N] = Xi - Xs 
    start_padded_series = np.zeros(N+t) 
    start_padded_series[t:] = Xi - Xs 

    auto_cov = 1./(N-1) * np.sum(start_padded_series*end_padded_series) 
    return auto_cov 

此針對@bluevoxel的代碼，用一個時間序列的50000個數據點，計算用於滯後的單一固定值的自相關相比，蟒for循環代碼平均約爲30毫秒和使用numpy陣列的平均速度超過0.3毫秒（運行在我的筆記本電腦上）。