Подтвердить что ты не робот

Как сопоставить два временных ряда с разрывами и разными базами времени?

У меня есть два временных ряда данных 3D-акселерометра, которые имеют разные временные базы (часы начинаются в разное время, с очень небольшим ползучести во время выборки), а также содержат много пробелов разного размера (из-за задержек, связанных с запись на отдельные флэш-устройства).

Акселерометры, которые я использую, - это недорогой GCDC X250-2. Я запускаю акселерометры с самым высоким коэффициентом усиления, поэтому данные имеют значительный уровень шума.

Временные ряды имеют около 2 миллионов точек данных (более часа в 512 выборок/сек) и содержат около 500 интересующих событий, где типичное событие охватывает 100-150 выборок (200-300 мс каждый). На многие из этих событий влияют сбои данных во время записи флэш-памяти.

Итак, данные не являются первозданными, и даже не очень хороши. Но мой контроль над глазом показывает, что он явно содержит интересующую меня информацию. (Я могу разместить сюжеты, если это необходимо.)

Акселерометры находятся в аналогичных средах, но только умеренно связаны, что означает, что я могу сказать, какие события соответствуют каждому акселерометру, но пока что я не увенчался успехом в программном обеспечении. Из-за физических ограничений устройства также монтируются в разных направлениях, где оси не совпадают, но они настолько близки к ортогональным, насколько я мог их изготовить. Так, например, для 3-осевых акселерометров A и B, + Ax отображает на -By (вверх-вниз), + Az отображает на -Bx (слева-справа) и + Ay отображает на -Bz (фронт-назад).

Моя первоначальная цель - сопоставить ударные события на вертикальной оси, хотя в конечном итоге мне хотелось бы: a) автоматически обнаружить отображение остов, b) коррелировать активность на отображаемых тузах и c) выявить различия в поведении между двумя акселерометрами ( таких как скручивание или сгибание).

Характер данных временного ряда делает Python numpy.correlate() непригодным. Я также посмотрел пакет R Zoo, но с ним не справился. Я посмотрел на различные области анализа сигналов для справки, но я не добился прогресса.

У кого-нибудь есть какие-то подсказки, что я могу сделать, или подходы, которые я должен исследовать?

Обновление 28 февраля 2011 г. Добавлено несколько графиков здесь, в которых приведены примеры данных.

4b9b3361

Ответ 1

Моя интерпретация вашего вопроса: учитывая два очень длинных, шумных временных ряда, найдите сдвиг одного, который соответствует большим "ударам" в одном сигнале, к большим ударам в другом сигнале.

Мое предложение: интерполировать данные так, чтобы они равномерно распределялись, исправляли и сглаживали данные (предполагая, что фаза быстрых колебаний неинтересна) и выполняют взаимную корреляцию по одной точке в момент времени (предполагая небольшой сдвиг выровнят данные).

import numpy
from scipy.ndimage import gaussian_filter
"""
sig1 and sig 2 are assumed to be large, 1D numpy arrays
sig1 is sampled at times t1, sig2 is sampled at times t2
t_start, t_end, is your desired sampling interval
t_len is your desired number of measurements
"""

t = numpy.linspace(t_start, t_end, t_len)
sig1 = numpy.interp(t, t1, sig1)
sig2 = numpy.interp(t, t2, sig2)
#Now sig1 and sig2 are sampled at the same points.

"""
Rectify and smooth, so 'peaks' will stand out.
This makes big assumptions about your data;
these assumptions seem true-ish based on your plots.
"""
sigma = 10 #Tune this parameter to get the right smoothing
sig1, sig2 = abs(sig1), abs(sig2)
sig1, sig2 = gaussian_filter(sig1, sigma), gaussian_filter(sig2, sigma)

"""
Now sig1 and sig2 should look smoothly varying, with humps at each 'event'.
Hopefully we can search a small range of shifts to find the maximum of the 
cross-correlation. This assumes your data are *nearly* lined up already.
"""
max_xc = 0
best_shift = 0
for shift in range(-10, 10): #Tune this search range
    xc = (numpy.roll(sig1, shift) * sig2).sum()
    if xc > max_xc:
        max_xc = xc
        best_shift = shift
print 'Best shift:', best_shift
"""
If best_shift is at the edges of your search range,
you should expand the search range.
"""

Ответ 2

Если данные содержат пробелы неизвестных размеров, которые различаются в каждом временном ряду, я бы отказался от попытки сопоставить целые последовательности и вместо этого попытался перекрестно сопоставить пары коротких окон в каждом временном ряду, например, перекрывая окна в два раза длина типичного события (300 проб длиной). Найдите потенциальные совпадения с высокой взаимной корреляцией по всем возможностям, а затем наложите последовательное ограничение порядка на возможные совпадения для получения последовательностей совпадающих окон.

Оттуда у вас есть более мелкие проблемы, которые легче анализировать.

Ответ 3

Похоже, вы хотите свести к минимуму функцию (Ax '+ By) + (Az' + Bx) + (Ay '+ Bz) для пары значений: А именно, смещение по времени: t 0 и временной масштаб: t r. где Ax '= t r * (Ax + t 0) и т.д.

Я бы посмотрел на двунаправленную SciPy optimize. И я использовал бы маску или временно нулевую информацию (как Ax ', так и By), например, через "пробелы" (при условии, что пробелы могут программно определяется).

Чтобы сделать процесс более эффективным, начните с грубой выборки A и B, но установите точность в fmin (или любом выбранном вами оптимизаторе), который соизмерим с вашей выборкой. Затем приступайте к получению более тонких окон полного набора данных до тех пор, пока ваши окна не станут узкими и не будут сбрасываться.

Изменить - совпадающие оси

Что касается вопроса о попытке определить, какая ось является колинейной с данной осью, и не зная о характеристиках ваших данных, я могу указать на аналогичный вопрос. Посмотрите pHash или любой другой метод, описанный в этой, чтобы помогают идентифицировать похожие формы сигналов.

Ответ 4

Это не технический ответ, но он может помочь вам придумать один:

  • Преобразуйте сюжет в изображение и вставьте его в программу достойного изображения, например gimp или photoshop.
  • разбивать графики на дискретные изображения всякий раз, когда есть пробел
  • помещает первую серию графиков в горизонтальную линию
  • поместите вторую серию в горизонтальную линию под ней
  • визуально идентифицировать первое коррелированное событие
  • Если два события не выстроены вертикально:
    • выберите любой экземпляр слева и все справа от него в этой строке
    • перетащите эти вещи вправо до тех пор, пока они не выстроились в линию

Это в значительной степени работает с аудиоредактором, поэтому если вы преобразуете его в простой аудиоформат, такой как несжатый WAV файл, вы можете манипулировать им напрямую в чем-то вроде Audacity. (Разумеется, это будет ужасно, но вы сможете легко перемещать участки данных.)

На самом деле, у дерзости есть язык сценариев, называемый nyquist, поэтому, если вам не нужна программа для обнаружения корреляций (или вы хотя бы готовы отложить этот шаг на данный момент), вы, вероятно, можете использовать некоторые сочетание маркеров громкости и nyquist для автоматизации выравнивания и экспорта чистых данных в выбранный вами формат после того, как вы отметите точки корреляции.

Ответ 5

Я предполагаю, что вам придется вручную создать таблицу смещений, которая выравнивает "совпадения" между рядами. Ниже приведен пример способа получения этих совпадений. Идея состоит в том, чтобы сдвинуть данные влево-вправо до тех пор, пока она не выровняется, а затем скорректируйте масштаб до тех пор, пока он не "соответствует". Попробуйте.

library(rpanel)

#Generate the x1 and x2 data
n1 <- rnorm(500)
n2 <- rnorm(200)
x1 <- c(n1, rep(0,100), n2, rep(0,150))
x2 <- c(rep(0,50), 2*n1, rep(0,150), 3*n2, rep(0,50))

#Build the panel function that will draw/update the graph
lvm.draw <- function(panel) {
       plot(x=(1:length(panel$dat3))+panel$off, y=panel$dat3, ylim=panel$dat1, xlab="", ylab="y", main=paste("Alignment Graph   Offset = ", panel$off, "   Scale = ", panel$sca, sep=""), typ="l")
       lines(x=1:length(panel$dat3), y=panel$sca*panel$dat4, col="red")
       grid()
       panel
}

#Build the panel
xlimdat <- c(1, length(x1))
ylimdat <- c(-5, 5)
panel <- rp.control(title = "Eye-Ball-It", dat1=ylimdat, dat2=xlimdat, dat3=x1, dat4=x2, off=100, sca=1.0, size=c(300, 160))
rp.slider(panel, var=off, from=-500, to=500, action=lvm.draw, title="Offset", pos=c(5, 5, 290, 70), showvalue=TRUE)
rp.slider(panel, var=sca, from=0, to=2, action=lvm.draw, title="Scale", pos=c(5, 70, 290, 90), showvalue=TRUE)