Взаимная корреляция нескольких последовательностей без цикла for

Прочитав это и попробовав np.correlate и cv2.matchTemplate, у меня все еще есть вопрос, который я, кажется, не могу решать.

У меня есть два массива numpy, каждый с формой (6000,50). 6000 последовательностей с каждыми 50 значениями. Теперь я хотел бы сделать взаимную корреляцию двух одномерных последовательностей этого массива, чтобы обнаружить сдвиг во времени. Я кратко попробовал openCV, но для меня это возвращает одно число (я ожидаю наивысшую корреляцию), поэтому теперь я использую numpy.correlate следующим образом:

np.correlate(x[2500], y[2500], mode='same')

(На графике взаимной корреляции я не ищу самый высокий пик, а ищу первый пик, используя это. Пример см. на графике)

введите здесь описание изображения

Как и следовало ожидать, я хотел бы сделать это для всех 6000 последовательностей, но надеясь избежать повторения. Я надеялся, что это сработает:

np.correlate(x, y, mode='same')

Но это дает мне следующую ошибку: ValueError: object too deep for desired array.

Есть ли какие-либо изменения, которые возможны с NumPy или OpenCV. Или мне придется сделать это так :(

for i in range(x.shape[0]):
    np.correlate(x[i], y[i], mode='same')

python numpy opencv cross-correlation
person Mattijn    schedule 18.12.2013    source источник
comment
Почему вы хотите избежать повторения? Даже если бы было что-то вроде того, что вы хотите, он просто упаковывает итерацию внутри. Сложность будет не ниже.   -  person Skyler    schedule 18.12.2013
comment
@Skyler, я подумал, может быть, это связано с каким-то матричным вычислением.   -  person Mattijn    schedule 18.12.2013
comment
Ну я вижу. Но я не знаю такой функции :-(   -  person Skyler    schedule 18.12.2013

Ответы (1)


arrow_upward
3
arrow_downward

scipy.ndimage.correlate1d похоже на то, что вам нужно , но он вещает только на первый массив, второй должен быть строго 1D, так что тут не повезло. И функции в scipy.signal выполняют многомерную корреляцию, а не 1D, как вам нужно. Таким образом, в стеке нет ничего, что могло бы решить вашу проблему.

Ради интереса вы всегда можете сделать это с помощью БПФ и теоремы взаимной корреляции:

def correlate1(a, b):
    c = np.empty_like(a)
    for j in range(len(a)):
        c[j] = np.correlate(a[j], b[j], 'same')
    return c

def correlate2(a, b):
    n = a.shape[-1]
    a_fft = np.fft.fft(a, n=2*n)
    b_fft = np.fft.fft(b, n=2*n)
    cc = np.fft.ifft(a_fft * b_fft.conj()).real
    return np.concatenate((cc[..., -n//2:], cc[..., :(n-1)//2 + 1]), axis=-1)

С вашим вариантом использования это не очень хорошая идея:

In [11]: a = np.random.rand(6000, 50)
    ...: b = np.random.rand(6000, 50)
    ...: 

In [12]: np.allclose(correlate1(a, b), correlate2(a, b))
Out[12]: True

In [13]: %timeit correlate1(a, b)
10 loops, best of 3: 37.5 ms per loop

In [14]: %timeit correlate2(a, b)
10 loops, best of 3: 71.8 ms per loop

Но у этого подхода есть свои достоинства, в основном для больших последовательностей:

In [15]: a = np.random.rand(50, 6000)
    ...: b = np.random.rand(50, 6000)
    ...: 

In [16]: %timeit correlate1(a, b)
1 loops, best of 3: 516 ms per loop

In [17]: %timeit correlate2(a, b)
10 loops, best of 3: 89.2 ms per loop
person Jaime    schedule 18.12.2013
comment
Спасибо за ответ. «Веселая» часть кажется очень полезной для больших последовательностей. Что-то, что может произойти однажды! Как насчет заполнения одномерных массивов нулем или единицей, а затем применения одной из функций scipy.signal? - person Mattijn; 18.12.2013
comment
Нет, функции scipy.signal выполняют корреляцию в обоих направлениях, и вы не можете преодолеть это с помощью заполнения. Правильно было бы иметь правильный ndimage.correlate1d... - person Jaime; 18.12.2013
comment
Ok. Спасибо за объяснение. Позвольте мне написать запрос на numpy github, может быть, корреляция1d будет обновлена ​​​​однажды - person Mattijn; 18.12.2013