A "seismic section" shows amplitude for m discrete x values along its horizontal axis times n discrete time values along its vertical axis: Seismic section

Peaks in amplitude (black) are centered on horizons; interfaces between different geological layers. The horizons are mostly continous but not necessarily horizontal nor flat lines. In the image above somebody has "tracked" the green and yellow horizons (using a crayon on a piece of paper?) that form the top and bottom of the "Ordovian" geological layer. These two lines are continous everywhere except for where they cross one of the vertical black lines: faults.

It is easy to implement a computer program that "draws" the colored lines by following the peaks, except for at the faults.

The problem here is that the peak "jump" in time and by how much differs from case to case.

One possible way to automate this is to assume that the shape of the signal around the peak is almost constant along a horizon. We can then reconstruct the amplitude around each peak sample by e.g. a Cosine series expansion with say 5 terms:

Cosine approx

The result is a feature vector for each peak sample. In the case above the feature vector would be:

F1 = [94.2496  192.7706 -211.4520  -82.8782   11.2105]

The feature vectors can then be used as input to a classification/clustering algorithm that figures out which peak samples are on the same horizon.

The problem with the feature vector above is that it does not force the peak samples to be linked together; in fact it disregards their positions in space altogether. This would be fine if each horizon had a very different shape around the peak. Unfortunately this is often not the case.

So we want to include the x postion into the feature vector. However many clustering algorithms are not able to deal with this since they would require that all observations within a cluster have almost the same x position.

Perhaps hierarchical clustering may be a solution to this problem?

As I understand "single-linkage clustering" would only require that each peak had a peak within the same cluster that was its neighbour (adjacent x positions).

This is almost what we want. Ideally I would want single-linkage clustering on the x-position and average linkage clustering on the Cosine coefficients.

Is this possible?

Thanks in advance for any answers!


1 Answer 1


Can't you do this just by choosing an appropriate distance metric? You then should be able to use any clustering algorithm that uses distances - single-linkage, average-linkage, DBSCAN, OPTICS, even k-medoids.

For example, you could use the distance function $$ d(a,b) := \alpha |x(a) - x(b)|^2 + (1-\alpha) \sqrt{\sum_{i=1}^{5} |c_i(a)-c_i(b)|^2}$$

Where $\alpha$ controls the balance between locality (which is deliberately squared, to punish larger differences more) and a similar shape using euclidean distance on attributes $c_1$ to $c_5$.

If you want an even more flexible approach, you could give GDBSCAN a try (Generalized DBSCAN). For this you need to define when two objects are similar (e.g. distance in x) along with a predicate when it starts/continues a cluster (e.g. similarity in the coefficients to the neighbours).

  • $\begingroup$ Thanks Anony-Mousse! I was thinking about this, but was not sure if the clustering algorithms accepted custom distance metrics. Are there any minimum requirements for a valid distance metric? I guess one thing is what e.g. a Matlab function allow one to call but another thing is what will make sure the clustering algorithm succeds. Do you think a custom distance metric also could make sure that no two observations within a cluster have the same x position? $\endgroup$
    – Andy
    Jun 13, 2012 at 15:33
  • 1
    $\begingroup$ Some algorithms have extra constraints on the distances. For k-means, the mean must be consistent to the distance function. This mostly holds for L_p norms. I read somewhere it only works for Euclidean distance, but I'd assume it works for any distance where the mean is optimal. I don't know about Matlab, as I don't use it. There are lots of crappy implementations around that are tied to particular distance functions. ELKI is pretty interesting here, see the tutorial on their web page on adding custom distance functions. They can be used almost everywhere then. $\endgroup$ Jun 13, 2012 at 15:44
  • $\begingroup$ Hi again Anony-Mousse. I think memory consumption and running time would be a big problem for many hierarchical agglomerate clustering algorithms. I might have 1e6 feature vectors. However I only want to specify the distance matrix for those close in x position. For all others I guess the distance could be considered = 0. Can DBSCAN or OPTICS do this? $\endgroup$
    – Andy
    Jun 16, 2012 at 17:17
  • $\begingroup$ DBSCAN and OPTICS do range queries. So if you have precomputed these neighbors, that should work fine for them and make them essentially O(n). The most expensive part in both is the neighborhood query, which adds n or log n (with a spatial index structure) to this. But precomputed neighbor lists with only the relevant neighbors is ideal for DBSCAN and OPTICS. $\endgroup$ Jun 17, 2012 at 9:10

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.