2013-11-21 14 views
9

Voglio calcolare l'auto-covarianza di 3 array X1, X2 e Y che sono tutti processi stazionari casuali. C'è qualche funzione in sciPy o in un'altra libreria in grado di risolvere questo problema?Come calcolare l'auto-covarianza in Python

risposta

4

Secondo la stima standard del coefficiente autocovarianza per segnali discreti, che può essere espresso per equazione:

enter image description here

... dove x(i) è un dato segnale (cioè vettore 1D specifico), k significa lo spostamento di x(i) segnale da k campioni, N è la lunghezza di x(i) segnale, e:

enter image description here

... che è semplice media, possiamo scrivere:

''' 
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)) 

Se volete normalizzare il coefficiente autocovarianza, che diventerà il coefficiente di autocorrelazione espresso come:

enter image description here

... quanto basta aggiungere al codice di cui sopra solo due linee aggiuntive:

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

Ecco piena script:

''' 
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()) 
+0

Numpy ha già tutto il necessario per calcolare la [correlazione] (https://docs.scipy.org/doc/numpy/reference/generated/numpy.correlate.html). (Che può anche essere accelerato con [scipy.signal.fftconvolve] (http://docs.scipy.org/doc/scipy-0.15.1/reference/generated/scipy.signal.fftconvolve.html). Che poi solo deve essere moltiplicato per la [varianza] (http://docs.scipy.org/doc/numpy/reference/generated/numpy.var.html) per ottenere l'autocovarianza. – Celelibi

1

Get campione auto covarianza:

# 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 
0

Un piccolo tweak alle risposte precedenti, che evita python for loop e utilizza invece operazioni di array numpy. Questo sarà più veloce se hai molti dati.

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 

Confrontando questo contro codice @bluevoxel s', utilizzando una serie temporale di 50.000 punti dati e calcolando l'autocorrelazione per un singolo valore fisso di lag, il codice di ciclo pitone for media di circa 30 millisecondi e l'utilizzo di array numpy ha una media superiore a 0,3 milli-secondi (in esecuzione sul mio laptop).

Problemi correlati