如何实时检测曲线中的“膝盖/肘部”（最大曲率）答案

【问题标题】：How to detect in real time a "knee/elbow" (maximal curvature) in a curve如何实时检测曲线中的“膝盖/肘部”（最大曲率）
【发布时间】：2018-05-17 08:39:34
【问题描述】：

在下面的曲线（蓝线）中，我试图检测应该位于 x = 2.5 附近的“膝盖/肘部”

这是我正在使用的一组值：

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}

我已经尝试过Kneedle algorithm 和formal definition of the curvature of a graph（有符号曲率）。我在使用 Kneedle 算法时遇到的问题是，在实时应用程序（嵌入式系统）中，我不知道 y 轴的最大值是哪个，所以我无法正确归一化点，也无法找到斜率值适用于所有情况。当使用图形曲率的正式定义时，我尝试用 5 阶多项式（绿线）拟合曲线，然后获取导数的值来计算曲率。然而，此方法找到 x = -2 附近的曲率，因为该点周围存在多项式曲率。

有人可以建议我一种检测膝盖/肘部的方法吗？

【问题讨论】：

你知道如何数值计算一阶导数吧？（斜率），二阶导数当然只是 THAT 的一阶导数。（斜率变化率）。如果您尝试对完整数据集进行曲线拟合，这实际上会从您需要的信息中“删除”信息，而不是添加信息。另一种方法是“区间二等分”，您将其分成两半，并寻找具有最大“斜率不一致”的“段”，您可以从每个范围内的 10 个样本中采样并获取每个样本的斜率.
顺便说一句，当我说“采取一阶导数”时，我并不是说你最终得到一个数字，而是一个新的 CURVE（点集）。如果不清楚。
Median, best fit of line from mesh skeleton的可能重复
问题：您的数据实际建模的是什么？因为，照原样，没有理由假设函数是严格递增的，所以期望函数回退更有意义，因此（可能）具有任意数量的相同值的最大曲率“点” .同样重要的是：您的信号是平滑的，还是我们只是猜测这是一个平滑的函数？
@Mike'Pomax'Kamermans 我的数据实际上是来自电流传感器的输入。事实上，数据有太多的噪音。尽管如此，使用 Savitzky-Golay 过滤器，我可以使曲线变得相当平滑。为什么你说期望功能回退更有意义？我不清楚。

标签： algorithm curve polynomials

【解决方案1】：

我从事膝/肘点检测已经有一段时间了。无论如何，我是专家。尽管如此，我还是在我目前实现的代码上测试了你的示例。

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 似乎提供了一个很好的近似值，如 DFDT/DSDT 和 S-方法。正如您提到的，knedle 与您的用例无关。

DFDT 代表动态一阶导数阈值。它计算一阶导数并使用阈值算法来检测膝/肘点。 DSDT 类似但使用二阶导数，我的评估表明它们具有相似的性能。

S 方法是 L 方法的扩展。 L 方法将两条直线拟合到曲线上，两条线之间的交点是膝/肘点。通过循环整体点、拟合线并评估 MSE（均方误差）来找到最佳拟合。 S 方法拟合 3 条直线，这提高了准确性，但也需要更多的计算。

我的所有代码都可以在 github 上公开获得。此外，此article 可以帮助您找到有关该主题的更多信息。它只有四页长，所以应该很容易阅读。您可以使用代码，如果您想讨论任何方法，请随时这样做。

【讨论】：

【解决方案2】：

这不是连续曲线，也不是连续曲线的近似值，因此我们不必在这里浪费时间：您有一个简单的多边形。相当于“曲率半径”的多边形就是入射角和折射角的差值：差值越小，“曲率”半径就越大。

如果您正确采样数据，我们可以计算每个数据点的角度差：

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)

findAngle 函数应该很容易实现，有一百万个教程介绍了如何用您最喜欢的语言实现它（它是许多语言的基本算术，甚至是某些语言的内置）。就是这样，我们已将 2D 数据转换为具有 (x,y,z) 坐标的 3D 数据，其中z 表示穿过该点的局部角度。

然后要找到任何“拐点”，我们可以查看所有三个数据点 (x-1)、(x) 和 (x+1) 的集合，并考虑那些局部角度差较大的x比其邻居的角度差异。 x 最大的不同是“赢家”：你找到了你的膝盖。（或者真的，a 膝盖，因为点数据使零保证不会上下多次，从而导致多边形有很多拐点）

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]

我已将差异计算留给您，因为您可能只想要 d1+d2 或 (d1+d2)/2，或者您可能想根据 d1 或 d2（但不是两者）进行切换0等

当然，这里要注意的重要一点是，当我们收集数据点时，我们可以“一次性”完成所有这些事情，因为新点不会影响旧点的位置，所以在某些时候n，我们已经知道到 n-1 的角度，并且我们已经知道到 n-2 的膝盖，所以我们可以在一次线性传递中计算所有这些值。好用又高效。

这种方法类似于寻找数据拟合曲线的导数的根，但是当我们只有数值数据时，有时最正确的方法是使用该数据，而不是“对您从中获取数据的原始事物的假定重建”。在这种情况下，您不知道数据源在样本点之间做了什么。也许它表现得很好，也许不是，但是你不知道，所以不要浪费周期来使问题复杂化，而是直接使用你知道是真的属性（当然， 它们的真实用途：这些是传感器读数，您的传感器绝对可能有噪音甚至有故障）。

【讨论】：

【解决方案3】：

答案取决于您如何获得这些值。下面的答案假设您在输入时获得一组这些 XY 数字（如您的示例中所示）。

还有你是如何定义曲率的。

对于曲率的几何意义（与轴方向无关），您可以将 B 点的曲率估计为abs( crossproduct(B-A, C-B) ) / length(C-A)，其中 A 和 C 是两个相邻点。因为平方根（长度需要）可能很昂贵，尤其是。在嵌入时，您可以改用该值的平方。

对于曲率的代数定义（二阶导数，很大程度上取决于轴的方向），请参见这个公式：https://mathformeremortals.wordpress.com/2013/01/12/a-numerical-second-derivative-from-three-points/

更新：对于曲率的几何意义，上面的公式只有在点分布均匀的情况下才可以。对于不同的点密度，最好使用a precise formula 作为曲率半径。 4*K = 2*abs( crossproduct(B-A, C-B) )。不幸的是，您需要一个平方根（您可以先将所有 3 个边的平方相乘，然后从乘积中取一个根）。

【讨论】：

通过曲率的代数定义，对于我提供的点集，我得到了很好的结果，但对于平滑的真实数据却没有。即使在过滤后使用真实数据时，曲率的几何含义似乎也能正常工作。只是为了确定我在做什么是正确的，abs(crossproduct(B-A, C-B)) 是否意味着 (B-A) x (C-B) 的叉积的范数？
@Marco 二维向量的叉积本质上是一个标量：allenchou.net/2013/07/cross-product-of-2d-vectors 所以我的意思是这个绝对值：cplusplus.com/reference/cstdlib/abs

【解决方案4】：

这是一个算法：

1) 将范围分成三等份。 A,B,C

2) 对于每个范围 (A,B,C)，计算任意点的最大斜率减去任意点的最小斜率，并将其称为 SlopeDiversity。

3) 无论哪三分之一（A、B 或 C）具有最大的坡度多样性，取其并从步骤 #1 开始重复整个过程。这就像二分法，但又是三分法。

也可能有一些二分算法有效，但划分为 3rds 是有道理的，因为您想丢弃不好的部分并保留最好的部分，这需要每次分成 3 个范围。

另一种“陈述”或执行的方法是，当您将其中的所有点与简单的线性回归模型进行比较时，只选择具有“最大总误差”的第 3 个（仅适用于第 3 个的线性公式） .

【讨论】：