2009-03-13 9 views
69

Ich muss Autokorrelation einer Reihe von Zahlen tun, die, wie ich es verstehe, nur die Korrelation der Menge mit sich selbst ist.Wie kann ich mit numpy.correlate Autokorrelationen durchführen?

Ich habe es versucht mit numpy Korrelat-Funktion, aber ich glaube nicht, das Ergebnis, da es fast immer einen Vektor gibt, wo die erste Zahl ist nicht die größte, wie es sein sollte.

Also, diese Frage ist eigentlich zwei Fragen:

  1. Was genau ist numpy.correlate tun?
  2. Wie kann ich es (oder etwas anderes) für die Autokorrelation verwenden?
+0

Siehe auch: http://stackoverflow.com/questions/12269834/is-there-any-numpy-autocorrellation-function-with-standardized-output für Informationen über normalisierte Autokorrelation. – amcnabb

Antwort

79

Ihre erste Frage zu beantworten, ist numpy.correlate(a, v, mode) die Faltung von a mit der Rückseite v Durchführung und gibt die von dem angegebenen Modus abgeschnitten Ergebnisse. Die definition of convolution, C (t) = Σ -∞ < i < ∞ ein i v t + i wo -∞ < t < ∞, von für die Ergebnisse ermöglicht -∞ bis ∞, aber man kann natürlich nicht Speichern Sie ein unendlich langes Array. . So ist es werden abgeschnitten hat, und das ist, wo der Modus kommt Es gibt 3 verschiedene Modi: voll, gleich, & gültig:

  • „voll“ Modus liefert Ergebnisse für jede t wo beide a und v haben einige Überschneidungen.
  • Der Modus "same" liefert ein Ergebnis mit der gleichen Länge wie der kürzeste Vektor (a oder v).
  • Der Modus "valid" gibt nur Ergebnisse zurück, wenn a und v einander vollständig überlappen. Die documentation für numpy.convolve gibt mehr Details zu den Modi.

Für Ihre zweite Frage, ich denke numpy.correlateist Sie die Autokorrelation zu geben, wird es geben Sie nur ein wenig mehr als gut. Die Autokorrelation wird verwendet, um herauszufinden, wie ähnlich ein Signal oder eine Funktion für sich selbst bei einer bestimmten Zeitdifferenz ist. Bei einer Zeitdifferenz von 0 sollte die Autokorrelation am höchsten sein, da das Signal mit sich selbst identisch ist. Sie haben also erwartet, dass das erste Element im Autokorrelations-Ergebnisfeld das größte ist. Die Korrelation beginnt jedoch nicht bei einer Zeitdifferenz von 0. Sie beginnt bei einer negativen Zeitdifferenz, schließt sich auf 0 und wird dann positiv.Das heißt, Sie erwarten:

Autokorrelation (a) = Σ -∞ < i < ∞ ein i v t + i wobei 0 < = t < ∞

Aber was du hast war:

Autokorrelation (a) = Σ -∞ ∞ < i < a i v t + i wo -∞ < t < ∞

Was Sie tun müssen, um die letzte Hälfte des Korrelationsergebnis wird nehmen, und das sollte die Autokorrelation sein die Sie suchen. Eine einfache Python-Funktion zu tun, das wäre:

def autocorr(x): 
    result = numpy.correlate(x, x, mode='full') 
    return result[result.size/2:] 

werden Sie, natürlich, müssen Fehler sicherstellen Überprüfung, dass x ist eigentlich ein 1-d-Array. Diese Erklärung ist wahrscheinlich auch nicht die mathematisch strengste. Ich habe Unendlichkeiten herumgeworfen, weil die Definition der Faltung sie benutzt, aber das gilt nicht unbedingt für die Autokorrelation. Also, der theoretische Teil dieser Erklärung mag etwas wackelig sein, aber hoffentlich sind die praktischen Ergebnisse hilfreich. Thesepages auf Autokorrelation sind ziemlich hilfreich, und können Ihnen einen viel besseren theoretischen Hintergrund geben, wenn Sie nichts dagegen haben, durch die Notation und schwere Konzepte waten.

+4

In aktuellen Builds von numpy kann der Modus 'same' angegeben werden, um genau das zu erreichen, was die A. Levy vorgeschlagen hat. Der Rumpf der Funktion könnte dann 'return numpy.correlate (x, x, mode = 'selbe')' –

+10

@DavidZwicker lesen, aber die Ergebnisse sind anders! 'np.korrelieren (x, x, mode = 'voll') [len (x) // 2:]! = np.korrelieren (x, x, modus = 'gleich')'. Zum Beispiel, x = [1,2,3,1,2]; np.correlate (x, x, mode = 'voll'); '{' >>> array ([2, 5, 11, 13, 19, 13, 11, 5, 2]) '}' np.korrelieren (x, x, mode = 'same'); '{' >>> array ([11, 13, 19, 13, 11])}. Der richtige ist: 'np.correlate (x, x, mode = 'voll') [len (x) -1:];' {'>>> array ([19, 13, 11, 5, 2]) }} see ** Der erste Gegenstand ** ist ** der größte **. – Developer

+7

Beachten Sie, dass diese Antwort die unnormale Autokorrelation angibt. – amcnabb

12

Die Autokorrelation gibt es in zwei Versionen: Statistisch und Faltung. Beide machen dasselbe bis auf ein kleines Detail: Erstere wird auf das Intervall [-1,1] normiert. Hier ist ein Beispiel dafür, wie man den statistisch man tun:

def acf(x, length=20): 
    return numpy.array([1]+[numpy.corrcoef(x[:-i], x[i:]) \ 
     for i in range(1, length)]) 
+5

Sie wollen '' numpy.corrcoef [x: -i], x [i:]) [0,1] '' in der zweiten Zeile als Rückgabewert von '' corrcoef'' ist eine 2x2 Matrix – luispedro

+0

Was ist der Unterschied zwischen den statistischen und den Faltungsautokorrelationen? –

9

Wie ich gerade in das gleiche Problem lief, Ich mag würde mit Ihnen ein paar Zeilen Code teilen. In der Tat gibt es einige ziemlich ähnliche Posts über Autokorrelation im Stackoverflow. Wenn Sie die Autokorrelation als a(x, L) = sum(k=0,N-L-1)((xk-xbar)*(x(k+L)-xbar))/sum(k=0,N-1)((xk-xbar)**2) definieren [das ist die Definition in IDL a_correlate Funktion gegeben, und es stimmt, was ich in der Antwort 2 von Frage sehen #12269834], dann scheint die folgenden die richtigen Ergebnisse zu erhalten:

import numpy as np 
import matplotlib.pyplot as plt 

# generate some data 
x = np.arange(0.,6.12,0.01) 
y = np.sin(x) 
# y = np.random.uniform(size=300) 
yunbiased = y-np.mean(y) 
ynorm = np.sum(yunbiased**2) 
acor = np.correlate(yunbiased, yunbiased, "same")/ynorm 
# use only second half 
acor = acor[len(acor)/2:] 

plt.plot(acor) 
plt.show() 

Wie Sie sehen, ich habe das mit einer Sinuskurve und einer gleichmäßigen Zufallsverteilung getestet, und beide Ergebnisse sehen so aus, als ob ich sie erwarten würde. Beachten Sie, dass ich mode="same" anstelle von mode="full" als die anderen verwendet habe.

13

Mit der numpy.corrcoef Funktion anstelle von numpy.correlate die statistischen Korrelation für eine Verzögerung von t zu berechnen:

def autocorr(x, t=1): 
    numpy.corrcoef(numpy.array([x[0:len(x)-t], x[t:len(x)]])) 
+0

Beziehen sich "Korrelationskoeffizienten" nicht auf die bei der Signalverarbeitung verwendete Autokorrelation und nicht auf die in der Statistik verwendete Autokorrelation? https://en.wikipedia.org/wiki/Autocorrelation#Signal_processing –

+0

@DanielPendergast Ich bin nicht so vertraut mit Signalverarbeitung. Aus den numpigen Dokumenten: "Pearson-Produkt-Moment-Korrelationskoeffizienten zurückgeben.". Ist das die Signalverarbeitungsversion? –

0

Ich denke, die wirkliche Antwort auf die Frage des OP kurz und bündig in diesem Auszug aus der Numpy.correlate Dokumentation enthalten ist :

mode : {'valid', 'same', 'full'}, optional 
    Refer to the `convolve` docstring. Note that the default 
    is `valid`, unlike `convolve`, which uses `full`. 

Dies bedeutet, dass, wenn ohne ‚mode‘ Definition verwendet, kehrt die Numpy.correlate Funktion einen Skalar, den gleichen Vektor für seine zwei Eingabeargumente, wenn angegeben (dh - wenn verwendet, um Autokorrelation durchzuführen).

1

Ich benutze Talib.CORREL für Autokorrelation wie diese, ich vermute, dass Sie mit anderen Paketen das gleiche tun könnten:

def autocorrelate(x, period): 

    # x is a deep indicator array 
    # period of sample and slices of comparison 

    # oldest data (period of input array) may be nan; remove it 
    x = x[-np.count_nonzero(~np.isnan(x)):] 
    # subtract mean to normalize indicator 
    x -= np.mean(x) 
    # isolate the recent sample to be autocorrelated 
    sample = x[-period:] 
    # create slices of indicator data 
    correls = [] 
    for n in range((len(x)-1), period, -1): 
     alpha = period + n 
     slices = (x[-alpha:])[:period] 
     # compare each slice to the recent sample 
     correls.append(ta.CORREL(slices, sample, period)[-1]) 
    # fill in zeros for sample overlap period of recent correlations  
    for n in range(period,0,-1): 
     correls.append(0) 
    # oldest data (autocorrelation period) will be nan; remove it 
    correls = np.array(correls[-np.count_nonzero(~np.isnan(correls)):])  

    return correls 

# CORRELATION OF BEST FIT 
# the highest value correlation  
max_value = np.max(correls) 
# index of the best correlation 
max_index = np.argmax(correls) 
5

Ihre Frage 1 wurde in mehreren ausgezeichneten Antworten hier bereits ausführlich diskutiert.

Ich dachte, mit Ihnen ein paar Zeilen Code zu teilen, mit denen Sie die Autokorrelation eines Signals nur auf der Grundlage der mathematischen Eigenschaften der Autokorrelation berechnen können. Das heißt, dass die Autokorrelation in der folgenden Weise berechnet werden:

  1. den Mittelwert aus dem Signal subtrahieren und erhält ein unverzerrte Signal

  2. Berechne die Fourier-Transformation des unvoreingenommenen Signals

  3. compute die spektrale Leistungsdichte des Signals, indem die quadratische Norm von jedem Wert der Fourier - Transformation des unverzerrten Signals

    die inverse Fourier - Transformation von die spektrale Leistungsdichte

  4. Normalisieren der inversen Fourier-Transformation der spektralen Leistungsdichte durch die Summe der Quadrate der das nicht vorgespannten Signals und nehmen nur die Hälfte des resultierenden Vektors

Den Code umwandeln, dies zu tun ist folgende:

def autocorrelation (x) : 
    """ 
    Compute the autocorrelation of the signal, based on the properties of the 
    power spectral density of the signal. 
    """ 
    xp = x-np.mean(x) 
    f = np.fft.fft(xp) 
    p = np.array([np.real(v)**2+np.imag(v)**2 for v in f]) 
    pi = np.fft.ifft(p) 
    return np.real(pi)[:x.size/2]/np.sum(xp**2) 
+0

ist es möglich, dass etwas nicht stimmt? Ich kann die Ergebnisse nicht mit anderen Autokorrelationsfunktionen vergleichen. Die Funktion sieht ähnlich aus, wirkt aber etwas gequetscht. – pindakaas

+0

@pindakaas könnten Sie genauer sein? Bitte geben Sie an, welche Diskrepanzen Sie mit welchen Funktionen finden. – Ruggero

+0

Warum nicht 'p = np.abs (f)'? – dylnan