Estimate changing time lag between signals Suppose I have two timeseries $x_1(t)$ and $x_2(t)$, shown in the image that I drew below. They look almost identical in general form, except some features are shifted slightly in time from one series to the other by an about $\delta t$. Sometimes they are measured first in $x_1$, and sometimes in $x_2$, so sometimes $\delta t$ is positive and sometimes negative. The change in lag between features is random, but constrained to be of a small range so that shifts never overlap and cancel.

I want to know, what is the mean time lag between the two signals?

Even better, would be to know the distribution of time lags so I could estimate the mean, variance and range.

Does anyone have ideas or suggestions for how to do this?



My own thoughts so far are: I can imagine calculating the cross correlation between signals $x_1$ and $x_2$. However, I would expect the result to be strong correlation at zero time lag, because the time lags can be both positive and negative, so that the average lag is zero. However, the width/stretching of the zero peak, must be related to the average time lag, however, I can't figure out how exactly at the moment.

Another thought I have had is to create a moving correlation estimator. In a similar way to how wavelet analysis can be used to estimate how the power of certain frequency components changes as a function of time. I could develop a function that measures the lag correlation as a function of time by breaking one of the signals into smaller chunks and convolving it across the second signal while stepping along.

All suggestions appreciated.
If the signals are as you've drawn them (flat, with abrupt changes), then perhaps this approach might work:


Take the absolute value of the differences between successive time samples.
The non-zero entries are the transition times.
Find these times for both signals.
Subtract the two.


I've made an attempt to do this in the R code below.

Example output:

 &gt; print(differences_in_delay)
  [1]  2 -5  2 -4 -1  5 -4  3 -3  0
 &gt; print(t[which(abs(diff(x2))&gt;0)] - t[which(abs(diff(x1))&gt;0)])
  [1]  2 -5  2 -4 -1  5 -4  3 -3  0


It breaks down when the transitions are within a sample or two of each other, but otherwise seems to work OK.

Will also probably break down if there's noise on your signals.



R Code Below

T &lt;- 1024    
t &lt;- seq(1,T,1)    
Nchanges &lt;- 10    
times_for_changes_for_x1 &lt;- sort(runif(Nchanges,T/10,9*T/10))
differences_in_delay &lt;- round(runif(Nchanges,-5,5))
times_for_changes_for_x2 &lt;- times_for_changes_for_x1 + differences_in_delay
changes_for_x1 &lt;- runif(Nchanges,-5,5)

x1 &lt;- 0*t
x2 &lt;- 0*t
change_no_x1 &lt;- 1
change_no_x2 &lt;- 1

for (tm in t)
{
  if ( (change_no_x1 &lt;= Nchanges) &amp;&amp; (tm &gt; times_for_changes_for_x1[change_no_x1]))
  {
    x1[tm] &lt;- x1[tm] + changes_for_x1[change_no_x1]
    change_no_x1 &lt;- change_no_x1 + 1
  }
  else
  {
    if (tm &gt; 1)
    {
      x1[tm] = x1[tm-1]
    }
  }
  if ( (change_no_x2 &lt;= Nchanges) &amp;&amp; (tm &gt; times_for_changes_for_x2[change_no_x2]))
  {
    x2[tm] &lt;- x2[tm] + changes_for_x1[change_no_x2]
    change_no_x2 &lt;- change_no_x2 + 1
  }
  else
  {
    if (tm &gt; 1)
    {
      x2[tm] = x2[tm-1]
    }
  }
}

plot(x1,type='l')
lines(x2,col="red")

plot(abs(diff(x1)&gt;0)*1, type='l')
lines(abs(diff(x2)&gt;0)*1, col="red")


print(differences_in_delay)
print(t[which(abs(diff(x2))&gt;0)] - t[which(abs(diff(x1))&gt;0)])