Как я могу выровнять / синхронизировать два сигнала?

21

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

Я собрал два одновременных сигнала: скорость потока, интегрированная для объема, и изменение объема груди. Я хотел бы сравнить сигналы и в конечном итоге надеяться получить объем из сигнала расширения груди. Но сначала я должен выровнять / синхронизировать свои данные.

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

Я приложил пример двух сигналов ( https://docs.google.com/spreadsheet/ccc?key=0As4oZTKp4RZ3dFRKaktYWEhZLXlFbFVKNmllbGVXNHc ), пожалуйста, дайте мне знать, если есть что-то еще, что я мог бы предоставить.

person157
источник
Я не знаю этого достаточно хорошо, чтобы дать ответ, и не уверен, что это решает вопрос, но один из подходов синхронизации сигналов называется «регистрация», которая является подмножеством функционального анализа данных. Эта тема обсуждается в книге FDA Рамси и Сильвермана. Основная идея заключается в том, что наблюдаемые сигналы могут быть «искажены» (например, если мы интересовались механизмами жевания людей, но у нас есть данные о людях, жующих с разными скоростями - в этом случае ось времени «искажается») и регистрация пытается определить базовый сигнал по общему, «не искривленному» масштабу.
Макро
1
Вы уже собрали все свои данные? Это один экспериментальный предмет? Если вы только начинаете, я бы рассмотрел возможность отделения сигнала от вашего пояса и использования его в качестве триггера (или даже просто для отметки временной метки) записи вашего потока. Обычно системы сбора данных имеют эту возможность со вспомогательным или триггерным портом. Я уверен, что есть способы отличить его, просто используя ваши данные, как предложил Макро, но добавление этого дополнительного шага избавит вас от многих догадок.
Йонска
1
Я думаю, вы хотите оценить только фиксированную задержку. Вы можете использовать взаимную
это
1
Возможно, вы захотите задать этот вопрос на dsp.SE, где они также думают о синхронизации сигналов.
Дилип Сарват
1
@Thias верен, но похоже, что первые серии должны быть пересчитаны, чтобы они имели общие интервалы.
whuber

Ответы:

24

Вопрос состоит в том, как найти величину, на которую один временной ряд («расширение») отстает от другого («объем»), когда ряды отбираются с регулярными, но разными интервалами.

В этом случае обе серии демонстрируют достаточно непрерывное поведение, как показано на рисунках. Это подразумевает, что (1) может потребоваться небольшое или нулевое начальное сглаживание, и (2) повторная выборка может быть такой же простой, как линейная или квадратичная интерполяция. Квадратичный может быть немного лучше из-за гладкости. После повторной выборки задержка определяется путем максимизации взаимной корреляции , как показано в потоке. Для двух рядов выборочных данных со смещением, какова наилучшая оценка смещения между ними? ,


Чтобы проиллюстрировать это , мы можем использовать данные, представленные в вопросе, используя Rпсевдокод. Начнем с базовой функциональности, взаимной корреляции и повторной выборки:

cor.cross <- function(x0, y0, i=0) {
  #
  # Sample autocorrelation at (integral) lag `i`:
  # Positive `i` compares future values of `x` to present values of `y`';
  # negative `i` compares past values of `x` to present values of `y`.
  #
  if (i < 0) {x<-y0; y<-x0; i<- -i}
  else {x<-x0; y<-y0}
  n <- length(x)
  cor(x[(i+1):n], y[1:(n-i)], use="complete.obs")
}

Это грубый алгоритм: расчет на основе БПФ будет быстрее. Но для этих данных (включая около 4000 значений) этого достаточно.

resample <- function(x,t) {
  #
  # Resample time series `x`, assumed to have unit time intervals, at time `t`.
  # Uses quadratic interpolation.
  #
  n <- length(x)
  if (n < 3) stop("First argument to resample is too short; need 3 elements.")
  i <- median(c(2, floor(t+1/2), n-1)) # Clamp `i` to the range 2..n-1
  u <- t-i
  x[i-1]*u*(u-1)/2 - x[i]*(u+1)*(u-1) + x[i+1]*u*(u+1)/2
}

Я загрузил данные в виде файла CSV, разделенного запятыми, и удалил его заголовок. (Заголовок вызвал некоторые проблемы для R, которые я не хотел диагностировать.)

data <- read.table("f:/temp/a.csv", header=FALSE, sep=",", 
                    col.names=c("Sample","Time32Hz","Expansion","Time100Hz","Volume"))

NB. Это решение предполагает, что каждая серия данных находится во временном порядке без пропусков ни в одной из них. Это позволяет ему использовать индексы в значениях в качестве прокси для времени и масштабировать эти индексы по частотам временной выборки, чтобы преобразовать их во времена.

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

n.clip <- 350      # Number of terminal volume values to eliminate
n <- length(data$Volume) - n.clip
indexes <- 1:n
v <- residuals(lm(data$Volume[indexes] ~ indexes))
expansion <- residuals(lm(data$Expansion[indexes] ~ indexes)

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

e.frequency <- 32  # Herz
v.frequency <- 100 # Herz
e <- sapply(1:length(v), function(t) resample(expansion, e.frequency*t/v.frequency))

Теперь можно вычислить взаимную корреляцию - для эффективности мы ищем только разумное окно лагов - и можно определить лаг, в котором найдено максимальное значение.

lag.max <- 5       # Seconds
lag.min <- -2      # Seconds (use 0 if expansion must lag volume)
time.range <- (lag.min*v.frequency):(lag.max*v.frequency)
data.cor <- sapply(time.range, function(i) cor.cross(e, v, i))
i <- time.range[which.max(data.cor)]
print(paste("Expansion lags volume by", i / v.frequency, "seconds."))

Вывод говорит нам, что расширение отстает от объема на 1,85 секунды. (Если последние 3,5 секунды данных не были обрезаны, результат будет 1,84 секунды.)

Это хорошая идея, чтобы проверить все несколькими способами, желательно визуально. Во-первых, функция взаимной корреляции :

plot(time.range * (1/v.frequency), data.cor, type="l", lwd=2,
     xlab="Lag (seconds)", ylab="Correlation")
points(i * (1/v.frequency), max(data.cor), col="Red", cex=2.5)

кросс-корреляционный график

Далее, давайте зарегистрируем две серии во времени и нанесем их вместе на одной оси .

normalize <- function(x) {
  #
  # Normalize vector `x` to the range 0..1.
  #
  x.max <- max(x); x.min <- min(x); dx <- x.max - x.min
  if (dx==0) dx <- 1
  (x-x.min) / dx
}
times <- (1:(n-i))* (1/v.frequency)
plot(times, normalize(e)[(i+1):n], type="l", lwd=2, 
     xlab="Time of volume measurement, seconds", ylab="Normalized values (volume is red)")
lines(times, normalize(v)[1:(n-i)], col="Red", lwd=2)

Зарегистрированные участки

Это выглядит довольно хорошо! Тем не менее, мы можем лучше понять качество регистрации на графике рассеяния . Я меняю цвета по времени, чтобы показать прогрессию.

colors <- hsv(1:(n-i)/(n-i+1), .8, .8)
plot(e[(i+1):n], v[1:(n-i)], col=colors, cex = 0.7,
     xlab="Expansion (lagged)", ylab="Volume")

разброс точек

Мы ищем точки для отслеживания вперед-назад вдоль линии: отклонения от этого отражают нелинейности в запаздывающем отклике расширения на объем. Хотя есть некоторые варианты, они довольно маленькие. Тем не менее, как эти изменения меняются со временем, может представлять некоторый физиологический интерес. Замечательная вещь в статистике, особенно в ее исследовательском и визуальном аспектах, заключается в том, что она имеет тенденцию создавать хорошие вопросы и идеи наряду с полезными ответами.

Whuber
источник
1
Черт возьми, ты потрясающий. Кросс-корреляция - это именно то, что я представлял (я знал, что для этого должно быть название), но ваш ответ / объяснение выходило за рамки. Большое спасибо!
person157
Сейчас у меня нет времени на полное объяснение, но в книгах «Числовые рецепты» появляется отличный отчет. Например, посмотрите на главу 13.2, «корреляции и автокорреляции с помощью БПФ,» в Numerical Recipes в C . Вы также можете посмотреть на acfфункцию Р.
whuber
Новичок в 'r', пожалуйста, будьте любезны: функция нормализации, использованная в комбинированном графике (2-й последний график), не будет работать для меня, есть ли обновление этой функции, так как этот ответ был опубликован?
CmKndy
1
@CmKndy Я тоже был новичком R, когда опубликовал этот ответ и просто забыл дать определение этой функции. Вот оригинал:normalize <- function(x) { x.max <- max(x); x.min <- min(x); dx <- x.max - x.min; if (dx==0) dx <- 1; (x-x.min) / dx }
whuber
Отлично, спасибо @whuber. Если бы вы могли опубликовать ответ, подобный этому, когда вы были
новичком