функция для вычисления бикогерентности

Question

Дорогие все, что я ищу функцию numpy / scipy для вычисления бикогерентности и автобикогеренции перед изучением 3-волнового взаимодействия.Спасибо за всевозможную помощь Никола

Answer 1

Лучший пакет для этого в Python Land - http://pypi.python.org/pypi/nitime

Он имеет несколько оценок когерентности, но я не очень внимательно к ним относился.Это пакет для нейровизуализации, но алгоритмы намеренно используют только nupy и scipy, поэтому могут использоваться другими приложениями.

Answer 2

Вот функция, которая опирается на функцию scipy.spectrogram (версия scipy> 0,17) и вычисляет бикогерентность между двумя сигналами.

Определение от Hagihira 2001 и Hayashi 2007. См. Wikipedia-bicoherence

Надеюсь, это поможет.

С уважением,

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

Answer 3

Возможно, это Matlab toolbox поможет;в общем, довольно просто перевести Matlab на Python.

Answer 4

Если вы ссылаетесь на нормализованную перекрестную спектральную плотность (, как определено в википедии ), тогда matplotlib.mlab.cohere сделает это.