用于拟合具有两条线的曲线的亚二次算法

Question

问题是找到具有两条线组成的折线的实值二维曲线（由点集给出）的最佳拟合。

一种蛮力方法是为曲线的每个点找到“左”和“右”线性拟合，并选择误差最小的对。我可以在迭代曲线的点时增量计算两个线性拟合，但我找不到增量计算误差的方法。因此，这种方法产生二次复杂度。

问题是是否有一种算法可以提供次二次复杂度？

第二个问题是是否有用于此类算法的方便的 C++ 库？

---

编辑 对于单线拟合，有公式：

m = (?xiyi - ?xi?yi/N) / (?xi2 - (?xi)2/N)
b = ?yi/N - m * ?xi/N

其中m是斜率，b是线的偏移量。拥有这样一个拟合误差公式可以以最好的方式解决问题。

Answer 1

免责声明：我不想弄清楚如何在 C++ 中执行此操作，因此我将使用 Python (numpy) 表示法。这些概念是完全可以转移的，因此您应该可以毫无困难地翻译回您选择的语言。

\n

假设您有一对数组 和 ，x其中y包含数据点，并且x是单调递增的。我们还假设您始终选择一个划分点，该划分点在每个划分中至少留下两个元素，因此方程是可解的。

\n

现在您可以计算一些相关量：

\n

N = len(x)\n\nsum_x_left = x[0]\nsum_x2_left = x[0] * x[0]\nsum_y_left = y[0]\nsum_y2_left = y[0] * y[0]\nsum_xy_left = x[0] * y[0]\n\nsum_x_right = x[1:].sum()\nsum_x2_right = (x[1:] * x[1:]).sum()\nsum_y_right = y[1:].sum()\nsum_y2_right = (y[1:] * y[1:]).sum()\nsum_xy_right = (x[1:] * y[1:]).sum()\n

\n

我们需要这些量（要初始化）的原因O(N)是您可以直接使用它们来计算线性回归参数的一些众所周知的公式。例如，最优m和b的y = m * x + b计算公式为

\n

\n\xce\xbc x = \xce\xa3x i /N\n\xce\xbc y = \xce\xa3y i /N\nm = \xce\xa3(x i - \xce\xbc x )(y i - \xce\xbc y ) / \xce\xa3(x i - \xce\xbc x ) 2 \nb = \xce\xbc y - m * \xce\xbc x \n

\n

误差平方和由下式给出

\n

\ne = \xce\xa3(y i - m * x i - b) 2 \n

\n

这些可以使用简单代数扩展为以下形式：

\n

\nm = (\xce\xa3x i y i - \xce\xa3x i \xce\xa3y i /N) / (\xce\xa3x i 2 - (\xce\xa3x i ) 2 /N)\nb = \xce \xa3y i /N - m * \xce\xa3x i /N\ne = \xce\xa3y i 2 + m 2 * \xce\xa3x i 2 + N * b 2 - 2 * m * \xce\xa3x i y i - 2 * b * \xce\xa3y i + 2 * m * b * \xce\xa3x i \n

\n

因此，您可以遍历所有可能性并记录最小的e：

\n

for p in range(1, N - 3):\n    # shift sums: O(1)\n    sum_x_left += x[p]\n    sum_x2_left += x[p] * x[p]\n    sum_y_left += y[p]\n    sum_y2_left += y[p] * y[p]\n    sum_xy_left += x[p] * y[p]\n\n    sum_x_right -= x[p]\n    sum_x2_right -= x[p] * x[p]\n    sum_y_right -= y[p]\n    sum_y2_right -= y[p] * y[p]\n    sum_xy_right -= x[p] * y[p]\n\n    # compute err: O(1)\n    n_left = p + 1\n    slope_left = (sum_xy_left - sum_x_left * sum_y_left * n_left) / (sum_x2_left - sum_x_left * sum_x_left / n_left)\n    intercept_left = sum_y_left / n_left - slope_left * sum_x_left / n_left\n    err_left = sum_y2_left + slope_left * slope_left * sum_x2_left + n_left * intercept_left * intercept_left - 2 * (slope_left * sum_xy_left + intercept_left * sum_y_left - slope_left * intercept_left * sum_x_left)\n\n    n_right = N - n_left\n    slope_right = (sum_xy_right - sum_x_right * sum_y_right * n_right) / (sum_x2_right - sum_x_right * sum_x_right / n_right)\n    intercept_right = sum_y_right / n_right - slope_right * sum_x_right / n_right\n    err_right = sum_y2_right + slope_right * slope_right * sum_x2_right + n_right * intercept_right * intercept_right - 2 * (slope_right * sum_xy_right + intercept_right * sum_y_right - slope_right * intercept_right * sum_x_right)\n\n    err = err_left + err_right\n    if p == 1 || err < err_min\n        err_min = err\n        n_min_left = n_left\n        n_min_right = n_right\n        slope_min_left = slope_left\n        slope_min_right = slope_right\n        intercept_min_left = intercept_left\n        intercept_min_right = intercept_right\n

\n

您可能还可以进行其他简化，但这足以拥有一个O(n)算法。

\n