Como posso alinhar / sincronizar dois sinais?

Estou fazendo alguma pesquisa, mas fiquei paralisado na fase de análise (deveria ter prestado mais atenção às minhas aulas de estatísticas).

Eu coletei dois sinais simultâneos: vazão integrada para volume e alteração na expansão torácica. Eu gostaria de comparar os sinais e, finalmente, esperar obter volume do sinal de expansão torácica. Mas primeiro tenho que alinhar / sincronizar meus dados.

Como a gravação não inicia exatamente ao mesmo tempo e a expansão do tórax é capturada por períodos mais longos, preciso encontrar os dados que correspondem aos meus dados de volume no conjunto de dados de expansão do tórax e ter uma medida de quão bem eles estão alinhados. Não tenho muita certeza de como fazer isso se os dois sinais não iniciarem exatamente ao mesmo tempo ou entre dados em diferentes escalas e em diferentes resoluções.

Anexei um exemplo dos dois sinais ( https://docs.google.com/spreadsheet/ccc?key=0As4oZTKp4RZ3dFRKaktYWEhZLXlFbFVKNmllbGVXNHc ), informe-me se houver mais alguma coisa que eu possa fornecer.

A pergunta pergunta como encontrar a quantidade pela qual uma série temporal ("expansão") fica atrasada para outra ("volume") quando a série é amostrada em intervalos regulares, mas diferentes .

Nesse caso, ambas as séries exibem um comportamento razoavelmente contínuo, como as figuras mostrarão. Isso implica (1) pouca ou nenhuma suavização inicial pode ser necessária e (2) a reamostragem pode ser tão simples quanto a interpolação linear ou quadrática. Quadrático pode ser um pouco melhor devido à suavidade. Após a reamostragem, o atraso é encontrado maximizando a correlação cruzada , conforme mostrado no encadeamento. Para duas séries de dados amostrados por offset, qual é a melhor estimativa do offset entre eles? .

Para ilustrar , podemos usar os dados fornecidos na pergunta, empregando Ro pseudocódigo. Vamos começar com a funcionalidade básica, correlação cruzada e reamostragem:

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

Este é um algoritmo bruto: um cálculo baseado em FFT seria mais rápido. Mas para esses dados (envolvendo cerca de 4000 valores), é bom o suficiente.

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
}

Baixei os dados como um arquivo CSV separado por vírgula e retirei o cabeçalho. (O cabeçalho causou alguns problemas para o R que eu não queria diagnosticar.)

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

NB Esta solução assume que cada série de dados está em ordem temporal, sem lacunas em nenhuma delas. Isso permite usar índices nos valores como proxies para o tempo e escalar esses índices pelas frequências de amostragem temporal para convertê-los em tempos.

Acontece que um ou ambos os instrumentos flutuam um pouco com o tempo. É bom remover essas tendências antes de prosseguir. Além disso, como há um afunilamento do sinal de volume no final, devemos cortá-lo.

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)

Reamostro as séries menos frequentes para obter o máximo de precisão possível no resultado.

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

Agora a correlação cruzada pode ser calculada - por eficiência, buscamos apenas uma janela razoável de defasagens - e o atraso em que o valor máximo é encontrado pode ser identificado.

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

A saída nos diz que a expansão diminui o volume em 1,85 segundos. (Se os últimos 3,5 segundos de dados não foram cortados, a saída seria 1,84 segundos.)

É uma boa ideia verificar tudo de várias maneiras, de preferência visualmente. Primeiro, a função de correlação cruzada :

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)

Em seguida, vamos registrar as duas séries no tempo e plotá-las juntas nos mesmos eixos .

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)

Parece muito bom! No entanto, podemos ter uma noção melhor da qualidade do registro com um gráfico de dispersão . Eu vario as cores pelo tempo para mostrar a progressão.

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

Procuramos os pontos a serem rastreados ao longo de uma linha: variações que refletem não linearidades na resposta com atraso de expansão do volume. Embora existam algumas variações, elas são bem pequenas. No entanto, como essas variações mudam ao longo do tempo pode ter algum interesse fisiológico. O maravilhoso das estatísticas, especialmente seu aspecto exploratório e visual, é como elas tendem a criar boas perguntas e idéias, além de respostas úteis .