How to detect in real time a "knee/elbow" (maximal curvature) in a curve

Question

In the following curve (blue line) I'm trying to detect the "knee/elbow" which should be located around x = 2.5

This is the set of values I'm using:

x = {-10, -9, -8, -7, -6, -5, -4, -3, -2 , -1, 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10}

y = {0, 10, 20, 30, 40, 50, 60, 70, 80, 90, 100, 107, 122, 145, 176, 215, 262, 317, 380, 451, 530, 617}

I have tried the Kneedle algorithm and the formal definition of the curvature of a graph (signed curvature). The problem I have with the Kneedle algorithm is that in the real time application (embedded system) I don't know which will be the maximal value of the y axis, so I can't properly normalize the points nor find a slope value which works for all the cases. When using the formal definition of the curvature of a graph I try to fit the curve with an order 5 polynomial (green line) and then get the values of the derivatives to calculate the curvature. Nevertheless, this method finds the curvature around x = -2 as there is a curvature around that point because of the polynomial.

Could somebody suggest me a way to detect the knee/elbow?

Answer 1

This is not a continuous curve, nor an approximation of a continuous curve, so we don't have to waste time here: you have a simple polygon. The polygon equivalent to "radius of curvature" is the difference in angle of incidence and angle of refraction: the smaller that difference is, the larger the radius of "curvature" is.

Provided you sample your data properly, we can compute the angular difference at each data point:

for (i=1, i<p.length -1):
  vector1 = p[i] - p[i-1] // assuming your language of choice has points as primitives
  vector2 = p[1+1] - p[i] // if not, you'll have to extract x/y separately.
  p[i].angle = findAngle(vector1,vector2)

The findAngle function should be easy enough to implement, there are a million tutorials on how to implement this in your favourite language (it's basic arithemtic in many, and even built-in for some, languages). And that's it, we've turned our 2D data into 3D data with (x,y,z) coordinates in which z represents the local angle in travel across the point.

To then find any "knee", we can look at all sets of three data points (x-1), (x), and (x+1) and consider those x where the local angular difference is larger than the angular differences of its neighbours. The x with the biggest difference is the "winner": you've found your knee. (or really, a knee, since point data makes zero promises not to go up and down many many times, resulting in a polygon with lots of inflections)

knee = undefined
max_diff = 0;
for (i=1, i<p.length -1):
  a = p[i].angle

  a1 = p[i-1].angle
  d1 = a1-a

  a3 = p[i+1].angle
  d2 = a3-a

  diff = ... // there's a few ways to compute this
  p[i].diff = diff // always useful if we need to do more work later

  if (diff > max_diff):
    max_diff = diff
    knee = p[i]

I've left the diff computation up to you, because maybe you just want d1+d2, or (d1+d2)/2, or maybe you want to switch based on whether d1 or d2 (but not both) are 0, etc.

The important thing to note here, of course, is that we can do all these things "in one go", as we gather data points, since new points don't affect the position of old points, so at some point n, we already know the angles up to n-1 and we already know the knees up to n-2, so we can compute all these values in a single linear pass. Nice and efficient.

This approach is analogous to finding the root(s) of the derivative of a data-fitting curve, but when all we have is numerical data, sometimes the most correct approach is to work with that data, not with "a presumed reconstruction of the original thing you derived the data from". In this case, you have no idea what your data source is doing between sample points. Maybe it's well behaved, maybe it's not, but you don't know, so it's well worth not wasting cycles on complicating the problem and instead work directly with the properties you know are true (and of course, what they are true for: these are sensor readings and your sensor can absolutely be noisy or even faulty).

Answer 2

I have been working on Knee/Elbow point detection for some time. By no means, I am an expert. Nevertheless, I tested your example on the code I have implemented so far.

Kneedle -> 3.0
L-method -> 4.0
S-method -> 2.0
Menger Curvature -> 1.0
DFDT -> 3.0
DSDT -> 3.0
R-Method -> NA (R2 = NAN)
DK-method -> 1.0

K-needle appears to give a good approximation, as DFDT/DSDT and S-method. As you mentioned kneedle is not relevant to your use case.

DFDT stands for Dynamic First Derivate Threshold. It computes the first derivative and uses a Thresholding algorithm to detect the knee/elbow point. DSDT is similar but uses the second derivative, my evaluation shows that they have similar performances.

S-method is an extension of the L-method. The L-method fits two straight lines to your curve, the interception between the two lines is the knee/elbow point. The best fit is found by looping overall points, fitting the lines and evaluate the MSE (Mean Square Error). The S-method fits 3 straight lines, this improves the accuracy but also requires some more computation.

All of my code is publicly available on github. Furthermore, this article can help you found more information about the topic. It is only four pages long so it should be easy to read. You can use the code, and if you want to discuss any of the methods feel free to do so.