Cher tous Je recherche une fonction numpy/scipy pour calculer la bicohérence et l'auto-bicohérence avant l'étude de l'interaction 3-ondes. Merci pour toute l'aide possible nicolafonction pour le calcul de la bicohérence
Répondre
Peut-être que cette Matlab toolbox aidera; il est assez facile de traduire Matlab en Python, en général.
Le meilleur package pour cela dans la terre python est http://pypi.python.org/pypi/nitime
Il a plusieurs estimateurs de cohérence, mais je ne regarde pas très attentivement ceux-ci. C'est un paquetage pour la neuro-imagerie, mais les algorithmes n'utilisent que numpy et scipy, intentionnellement, pour pouvoir être utilisés par d'autres applications.
Si vous vous référez à la densité spectrale croisée normalisée (as defined in wikipedia) alors matplotlib.mlab.cohere ferait l'affaire.
Voici une fonction qui repose sur la fonction scipy.spectrogram (version scipy> 0.17) et calcule la bicohérence entre deux signaux.
Définition Hagihira 2001 et 2007. Voir Hayashi Wikipedia-bicoherence
Hope this helps.
Cordialement,
def compute_bicoherence(s1, s2, rate, nperseg=1024, noverlap=512):
""" Compute the bicoherence between two signals of the same lengths s1 and s2
using the function scipy.signal.spectrogram
"""
from scipy import signal
import numpy
# compute the stft
f1, t1, spec_s1 = signal.spectrogram(s1, fs = rate, nperseg = nperseg, noverlap = noverlap, mode = 'complex',)
f2, t2, spec_s2 = signal.spectrogram(s2, fs = rate, nperseg = nperseg, noverlap = noverlap, mode = 'complex')
# transpose (f, t) -> (t, f)
spec_s1 = numpy.transpose(spec_s1, [1, 0])
spec_s2 = numpy.transpose(spec_s2, [1, 0])
# compute the bicoherence
arg = numpy.arange(f1.size/2)
sumarg = arg[:, None] + arg[None, :]
num = numpy.abs(
numpy.mean(spec_s1[:, arg, None] * spec_s1[:, None, arg] * numpy.conjugate(spec_s2[:, sumarg]),
axis = 0)
) ** 2
denum = numpy.mean(
numpy.abs(spec_s1[:, arg, None] * spec_s1[:, None, arg]) ** 2, axis = 0) * numpy.mean(
numpy.abs(numpy.conjugate(spec_s2[:, sumarg])) ** 2,
axis = 0)
bicoh = num/denum
return f1[arg], bicoh
# exemple of use and display
freqs, bicoh = compute_bicoherence(s1, s2, rate)
f = plt.figure(figsize = (9, 9))
plt.pcolormesh(freqs, freqs, bicoh,
# cmap = 'inferno'
)
plt.colorbar()
plt.clim(0, 0.5)
plt.show()