Spectral Coherence in Julia

You are essentially trying to perform clustering of the signals, group similar signals together. Since you have a relatively small number of signals, you can use a wide range of tools, including computing all pairwise distances between signals etc.

Without having much insight into the actual properties of the signals, it’s hard to judge whether or not DTW is a suitable metric in this case, but you could attempt to compute the DTW cost matrix using something like (naive)

d = DTW(radius=radius, dist=SqEuclidean())
D = [d(a,b) for a in signals, b in signals]

For a slightly more efficient version, see SpectralDistances.distmat (may not be needed in your case).

If you want to further process the distance matrix D into a clustering, you could pass it into any of the methods from Clustering.jl that accept a distance matrix, such as

using Clustering
cr = hclust(Symmetric(sqrt.(D)))
assignments = cutree(cr,k=3) # k is desired number of clusters

You could also try the built-in clustering dbaclust in DynamicAxisWarping

I do not think that MatrixProfile is a suitable tool in this case, at least not the standard techniques associated with the matrix profile. You are certainly looking for something more towards computing distances/similarity and clustering.

2 Likes

@baggepinnen

Thanks a lot for the recomendation. I will check all of this out and give it a try.

Just in case is matters. My signals happened randomly in time (they are little earthquakes). There should be between 2 and 5 sources located in different locations, and the signals arrives to my receivers located a few meters away from the sources. In theory each source should produce exactly the same signal, that travels a constant path to the receivers, hence the waveform should be exactly the same (except for randomly occurring noise). This is not entirely true as the source-time function differs slightly for every earthquake, but our hope is to be able to cluster the signals and say, “well all this events look similar so they should be from source A”.

I’m by no means a seismologist, but I’m not sure if this is correct. Each signal will be colored by the transfer function from the source to the receiver, and unless all signals propagate through exactly the same path, they will differ even without the presence of any random noise.

Since you were originally considering coherence-based measures, I assume that the propagation path is reasonably linear (an LTI system). If that’s the case, you could potentially try to fit pairwise linear models between each pair of signals, with one signal as input and one as output. If you account for the delay between the signals, there should be a linear system that relates each pair of signals that originates from the same event. In spirit, this would be similar to a coherence analysis, but could potentially perform better in this low-data regime where you only have a single transient event to compare. The metric you would compare would then be, e.g., the prediction error that the estimated model makes on the data (normalized to account for length of data record etc.)

1 Like

Hello again @baggepinnen.

Ok this is quite deeper that I thought hehehe.

Yes, you are right, the path is rather linear and uncomplicated (no weird stuff that would make scattered signals), the source time function is similar, and signals propagate through a path that is very similar.

Here a quick example of some signals selected by coherence. Some signals are very similar to each other, but I “see” at least two families here.

I thought maybe something out of DynamicAxisWarping or SpectralDistances could give me a measure of how similar they are, and help me differentiate them. I am not really familiar with pairwise models, can you suggest any tool in julia I can use to make some test?

My suggestion is still to do something like this

A good explanation of coherence and how to calculate it is given about 2/3 of the way down this page. I assume this is what DSP.jl does.

It looks like the cross-correlation function should also be considered for what you are doing. It looks like DSP.jl is your friend here, though I have not used this function.

A friend did research on syncro-squeezing for the medical field. This was described as a time-frequency approach. I am not at all sure if it would be helpful, but may be worth taking a look at.

@baggepinnen
Dear Frederik thanks a lot for all your input :slight_smile:
I will review all of this information and I am sure I will come up with some nice results.
Thanks again and cheers!

Hi Jake
This is nice input, thanks a lot :slight_smile: