适合3D中的一组点的平面:scipy.optimize.minimize vs scipy.linalg.lstsq

Question

给定3D中的一组点,一般问题是以下列形式找到a, b, c平面方程的系数:

z = a*x + b*y + c

这样得到的平面最适合那组点.

1. 在这个SO答案中,函数scipy.optimize.minimize用于解决此问题.

   它依赖于系数的初始猜测,并最小化误差函数,该函数将每个点与平面表面的距离相加.

2. 在此代码中(基于此其他SO答案),scipy.linalg.lstsq函数用于解决相同的问题(当限制为一阶多项式时).

   它C在等式中求解z = A*C,其中点集A的x,y坐标的连接是集合z的z坐标,并且C是a,b,c系数.

   与上述方法中的代码不同,这个代码似乎不需要对平面系数进行初始猜测.

由于该 minimize函数需要初始猜测,这意味着它可能会或可能不会收敛到最优解(取决于猜测的好坏程度).第二种方法是否有类似的警告,或者它会返回一个始终精确的解决方案吗？

Answer 1

scipy.linalg.lstsq保证最小二乘()收敛.事实上,有一种封闭形式的解析解(由(A^T A)^-1 A^Tb(^T矩阵转置和^-1矩阵求逆)

然而,标准优化问题通常无法解决 - 我们无法保证找到最小化值.然而,对于给定的等式,找到一些a, b, c这样的z = a*x + b*y + c,我们有一个线性优化问题(约束和目标在我们试图优化的变量中是线性的).线性优化问题通常是可解决的,因此scipy.optimize.minimize应收敛到最佳值.

注意:即使我们这样做,这在我们的约束中也是线性的z = a*x + b*y + d*x^2 + e*y^2 + f- 我们不必将自己局限于线性空间(x,y),因为我们已经有了这些点(x, y, x^2, y^2).对于我们的算法,这些看起来就像矩阵中的点A.所以我们实际上可以使用最小二乘法获得更高阶的多项式!

简而言之:实际上,所有不使用精确解析解的求解器通常都会在实际答案的某个可接受的范围内停止,因此我们很少得到确切的解决方案,但它往往是如此接近我们在实践中接受它.此外,即使最小二乘解算器很少使用解析解,而是采用像牛顿方法更快的方法.

如果您要更改优化问题,则不会这样.存在某些类型的问题,我们通常可以找到最佳值(这些问题中最大的类称为凸优化问题 - 尽管存在许多非凸问题,我们可以在某些条件下找到最佳值).

如果您有兴趣了解更多信息,请查看Boyd和Vandenberghe的Convex Optimization.第一章不需要太多的数学背景,它概述了一般优化问题以及它与线性和凸规划等可解决优化问题的关系.

Answer 2

我想用另一种方法完成答案，以便找到适合 R^3 中一组点的最佳平面。实际上，该lstsq方法非常有效，除非在（例如）所有点的 x 坐标为 0（或相同）的特定情况下。在这种情况下，使用的 A 矩阵的列lstsq不是线性无关的。例如：

A = [[ 0   y_0    1]
     [ 0   y_1    1]
     ...
     [ 0   y_k    1] 
     ...
     [ 0   y_N    1]]

为了规避这个问题，你可以直接svd在点集的中心坐标上使用。实际上，svd用于lstsq但不在同一个矩阵中。

这是一个python示例，给出了coords数组中点的坐标：

# barycenter of the points
# compute centered coordinates
G = coords.sum(axis=0) / coords.shape[0]

# run SVD
u, s, vh = np.linalg.svd(coords - G)

# unitary normal vector
u_norm = vh[2, :]

使用这种方法，vh矩阵是3x3在其行中包含正交向量的矩阵。前两个向量在平面上形成正交基，第三个向量是垂直于平面的单位向量。

如果你真的需要a, b, c参数，你可以从法向量中得到它们，因为法向量的坐标是(a, b, c)，假设平面的方程是ax + by + cz + d = 0。