Wie kann ich zwei Signale ausrichten / synchronisieren?

Ich mache einige Nachforschungen, stecke aber in der Analysephase fest (hätte meinen Statistikvorträgen mehr Aufmerksamkeit schenken sollen).

Ich habe zwei gleichzeitige Signale gesammelt: die integrierte Durchflussrate für das Volumen und die Änderung der Brustausdehnung. Ich möchte die Signale vergleichen und letztendlich hoffen, die Lautstärke aus dem Brustexpansionssignal abzuleiten. Aber zuerst muss ich meine Daten ausrichten / synchronisieren.

Da die Aufnahme nicht genau zur gleichen Zeit beginnt und die Brustexpansion über einen längeren Zeitraum erfasst wird, muss ich die Daten, die meinen Volumendaten entsprechen, im Brustexpansionsdatensatz finden und ein Maß dafür haben, wie gut sie ausgerichtet sind. Ich bin mir nicht ganz sicher, wie ich vorgehen soll, wenn die beiden Signale nicht genau zur gleichen Zeit oder zwischen Daten mit unterschiedlichen Maßstäben und Auflösungen starten.

Ich habe ein Beispiel der beiden Signale angehängt ( https://docs.google.com/spreadsheet/ccc?key=0As4oZTKp4RZ3dFRKaktYWEhZLXlFbFVKNmllbGVXNHc ). Bitte teilen Sie mir mit, ob ich weitere Informationen bereitstellen kann.

Bei der Frage wird gefragt, um wie viel Zeitreihen ("Expansion") einem anderen Zeitreihen ("Volume") nacheilen, wenn die Reihen in regelmäßigen, aber unterschiedlichen Intervallen abgetastet werden.

In diesem Fall zeigen beide Serien ein relativ kontinuierliches Verhalten, wie die Abbildungen zeigen werden. Dies impliziert, dass (1) möglicherweise keine oder nur eine geringe anfängliche Glättung erforderlich ist und (2) die Neuabtastung so einfach wie die lineare oder quadratische Interpolation sein kann. Quadratisch kann aufgrund der Glätte etwas besser sein. Nach dem erneuten Abtasten wird die Verzögerung durch Maximieren der Kreuzkorrelation ermittelt , wie im Thread gezeigt. Was ist die beste Schätzung des Versatzes zwischen zwei versetzten abgetasteten Datenreihen? .

Zur Veranschaulichung können wir die in der Frage angegebenen Daten unter Verwendung Rdes Pseudocodes verwenden. Beginnen wir mit der Grundfunktionalität, der Kreuzkorrelation und dem Resampling:

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")
}

Dies ist ein grober Algorithmus: Eine FFT-basierte Berechnung wäre schneller. Aber für diese Daten (mit ungefähr 4000 Werten) ist es gut genug.

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
}

Ich habe die Daten als kommagetrennte CSV-Datei heruntergeladen und den Header entfernt. (Der Header verursachte einige Probleme für R, die ich nicht diagnostizieren wollte.)

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

Hinweis: Bei dieser Lösung wird davon ausgegangen , dass jede Datenreihe in zeitlicher Reihenfolge vorliegt und keine Lücken aufweist. Dies ermöglicht es, Indizes in die Werte als Proxy für die Zeit zu verwenden und diese Indizes durch die zeitlichen Abtastfrequenzen zu skalieren, um sie in Zeiten umzuwandeln.

Es stellt sich heraus, dass eines oder beide dieser Instrumente im Laufe der Zeit etwas abweichen. Es ist gut, solche Trends zu entfernen, bevor Sie fortfahren. Da sich das Lautstärkesignal am Ende verjüngt, sollte es abgeschnitten werden.

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)

Ich nehme die weniger häufigen Serien erneut auf, um das Ergebnis so genau wie möglich zu machen.

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

Jetzt kann die Kreuzkorrelation berechnet werden - aus Effizienzgründen wird nur ein vernünftiges Zeitfenster durchsucht - und die Zeitverzögerung, in der der Maximalwert gefunden wird, kann identifiziert werden.

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."))

Die Ausgabe gibt an, dass die Erweiterung der Lautstärke um 1,85 Sekunden nacheilt. (Wenn die letzten 3,5 Sekunden der Daten nicht abgeschnitten würden, würde die Ausgabe 1,84 Sekunden betragen.)

Es ist eine gute Idee, alles auf verschiedene Arten zu überprüfen, am besten visuell. Erstens die Kreuzkorrelationsfunktion :

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)

Als nächstes registrieren wir die beiden Serien in der Zeit und zeichnen sie zusammen auf den gleichen Achsen .

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)

Es sieht ziemlich gut aus! Mit einem Streudiagramm können wir uns jedoch ein besseres Bild von der Registrierungsqualität machen . Ich ändere die Farben von Zeit zu Zeit, um den Fortschritt zu zeigen.

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

Wir suchen nach Punkten, die entlang einer Linie vor- und zurückverfolgt werden können: Abweichungen davon spiegeln Nichtlinearitäten in der zeitverzögerten Reaktion der Expansion auf die Lautstärke wider. Obwohl es einige Variationen gibt, sind sie ziemlich klein. Es kann jedoch von physiologischem Interesse sein , wie sich diese Schwankungen im Laufe der Zeit ändern. Das Wunderbare an Statistiken, insbesondere an ihrem explorativen und visuellen Aspekt, ist, wie sie dazu neigen, gute Fragen und Ideen zusammen mit nützlichen Antworten hervorzubringen.