Comment puis-je utiliser numpy.correlate pour réaliser une autocorrélation ?

Question

Je dois faire l'auto-corrélation d'un ensemble de nombres, ce qui, selon ce que je comprends, est simplement la corrélation de l'ensemble avec lui-même.

J'ai essayé de le faire en utilisant la fonction de corrélation de numpy, mais je ne crois pas en le résultat, car il donne presque toujours un vecteur où le premier nombre n'est pas le plus grand, comme il devrait l'être.

Par conséquent, cette question se divise en deux questions :

1. Que fait précisément numpy.correlate?
2. Comment puis-je l'utiliser (ou autre chose) pour faire de l'auto-corrélation?

Answer 1

Pour répondre à votre première question, numpy.correlate(a, v, mode) effectue la convolution de a avec l'inverse de v et donne les résultats tronqués par le mode spécifié. La définition de la convolution, C(t)=∑ _{-∞ < i < ∞} a_iv_t+i où -∞ < t < ∞, permet d'obtenir des résultats de -∞ à ∞, mais vous ne pouvez évidemment pas stocker un tableau de taille infinie. Il doit donc être tronqué, et c'est là que le mode intervient. Il existe 3 modes différents : complet (full), même (same) et valide (valid) :

• Le mode complet (full) renvoie des résultats pour chaque t où a et v se chevauchent.
• Le mode même (same) renvoie un résultat de même longueur que le vecteur le plus court (a ou v).
• Le mode valide (valid) ne renvoie des résultats que lorsque a et v se chevauchent complètement. La documentation pour numpy.convolve fournit plus de détails sur les modes.

Pour votre deuxième question, je pense que numpy.correlate vous donne en fait l'autocorrélation, mais vous donne un peu plus aussi. L'autocorrélation est utilisée pour déterminer à quel point un signal, ou une fonction, est similaire à lui-même à un certain retard. Pour un retard de 0, l'auto-corrélation devrait être la plus élevée car le signal est identique à lui-même, donc vous vous attendiez à ce que le premier élément du tableau de résultat d'autocorrélation soit le plus grand. Cependant, la corrélation ne commence pas à un retard de 0. Elle commence à un retard négatif, proche de 0, puis devient positive. Vous vous attendiez donc à :

autocorrélation(a) = ∑ _{-∞ < i < ∞} a_iv_t+i où 0 <= t < ∞

Mais ce que vous avez obtenu était :

autocorrélation(a) = ∑ _{-∞ < i < ∞} a_iv_t+i où -∞ < t < ∞

Vous devez prendre la dernière moitié de votre résultat de corrélation, et cela devrait être l'autocorrélation que vous recherchez. Une simple fonction Python pour cela serait :

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

Vous aurez bien sûr besoin de vérifier les erreurs pour vous assurer que x est effectivement un tableau à une dimension. De plus, cette explication n'est probablement pas la plus rigoureuse sur le plan mathématique. J'ai parlé d'infinis car la définition de la convolution les utilise, mais cela ne s'applique pas nécessairement à l'autocorrélation. Ainsi, la partie théorique de cette explication peut être légèrement confuse, mais j'espère que les résultats pratiques vous seront utiles. Ces pages sur l'autocorrélation sont plutôt utiles et peuvent vous donner un arrière-plan théorique bien meilleur si vous ne vous souciez pas de vous plonger dans la notation et les concepts complexes.

Answer 2

L'auto-corrélation existe en deux versions : statistique et convolution. Elles font toutes les deux la même chose, à l'exception d'un petit détail : la version statistique est normalisée pour être dans l'intervalle [-1,1]. Voici un exemple de comment vous réalisez la version statistique :

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

Answer 3

Utilisez la numpy.corrcoef fonction au lieu de numpy.correlate pour calculer la corrélation statistique pour un décalage de t :

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

Answer 4

Comme je viens de rencontrer le même problème, je voudrais partager quelques lignes de code avec vous. En fait, il y a plusieurs publications assez similaires sur l'autocorrélation dans stackoverflow maintenant. Si vous définissez l'autocorrélation comme a(x, L) = sum(k=0,N-L-1)((xk-xbar)*(x(k+L)-xbar))/sum(k=0,N-1)((xk-xbar)**2) [c'est la définition donnée dans la fonction a_correlate de l'IDL et cela correspond à ce que je vois dans la réponse 2 de la question #12269834], alors ce qui suit semble donner les résultats corrects:

import numpy as np
import matplotlib.pyplot as plt

# générer des données
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
# utiliser seulement la seconde moitié
acor = acor[len(acor)/2:]

plt.plot(acor)
plt.show()

Comme vous pouvez le constater, j'ai testé cela avec une courbe sinus et une distribution aléatoire uniforme, et les deux résultats semblent correspondre à ce que je m'attendais. Notez que j'ai utilisé mode="same" au lieu de mode="full" comme les autres l'ont fait.

Answer 5

1) Voici la documentation pour numpy.correlate. Le code à l'intérieur du fichier ressemble à ceci :

mode = _mode_from_name(mode)
return multiarray.correlate(a,v,mode)

multiarray.correlate pointe vers un fichier .pyd (c'est-à-dire un fichier DLL), donc pour obtenir le fonctionnement interne, vous devriez probablement demander aux développeurs de numpy.

2) Si vous ne croyez pas aux résultats de numpy, vous pouvez essayer la fonction correlate de SciPy.