拟合平面数据的三条相交线



假设平面上有三组点。分别称它们为红色组、蓝色组和黄色组。

假设每个集合至少包含两个不同的点。

我想找到三条线——分别叫它们红线、蓝线和黄线——这样:

i)红、蓝、黄线有一个共同点

ii)红线"合身";红色设置尽可能接近

iii)蓝线"fits"蓝色的集合尽可能接近

iv)黄线"fits"黄色设置尽可能接近

注意:

  1. 没有第一个条件,这是线性代数中的标准练习,通常被称为"最小二乘最小化"。

  2. 我想到了某种最小二乘拟合,但如果它们导致更简单的解决方案(或任何解决方案),我愿意接受其他最佳拟合标准。

这实际上是我遇到的实际问题的稍微简化,这是一个更高维度的版本。在这个版本中,三组数据位于空间中,而不仅仅是一个平面,我试图将数据与共享同一条线的三个平面相匹配。

所以你必须选择一个公共向量c,和G向量d[g[i]],其中g[i]表示点v[i]的组,你想最小化|d[g[i]] * t[i] - v[i] - c|,通过选择d,tc。注意两个变量d[g[i]] * t[i]的乘积,一般来说这个问题不再是凸的了。

你可以简单地把它放到非线性优化中

也许如果你假设t[i] > 0,你可以找到一个凸公式(在这里猜测)。

你可以尝试一些启发式的方法,例如,你可以开始解决个别线性回归问题,给定的解决方案找到点c,使所有的和的平方残差最小,然后你把这个点添加到所有的点的集合,并再次拟合,重复添加点的权重开始增加,希望你会找到一个好的解决方案。

虽然是线性拟合,但这听起来更像是正交距离问题。我同意2D解决方案应该提供如何在3D中做到这一点的提示。对于正交距离,我用的是非线性拟合。如下所示:

import matplotlib.pyplot as plt
import numpy as np
from scipy.optimize import least_squares
"""
just some noise
"""
def norm2( s ):
return np.random.normal( size=2, scale=s )
"""
definition of a line that allows arbitrary slope
actually, for a single line this has an unnecessary parameter,
as 2 should be enough, but with the upcoming restriction in mind, 
that's the one
"""
def lin( t, x0, y0, alpha, disorder=0 ):
d = np.array( [ x0, y0 ] )
m = np.array( [ np.cos( alpha ), np.sin( alpha ) ] )
if isinstance( t, ( float, int ) ):
out =  d + t * m + norm2( disorder )
elif isinstance( t, (tuple, list, np.ndarray ) ):
out = [ d + item * m + norm2( disorder ) for item in t ]
else:
raise TypeError
return np.asarray( out )
"""
test data
"""
tl = np.linspace( -7, 7, 25 )
pl1 = lin( tl, 3, 4, +0.133, disorder=0.1 )
pl2 = lin( tl, 3, 4, -1.150, disorder=0.1 )
pl3 = lin( tl, 3, 4, +1.513, disorder=0.1 )
"""
distance from line to point whre line is given by r + m * t 
and ||m|| = 1
(actually the least sqaure function is sqaring that again later, but ok )
"""
def dist_lin_pnt( r, m, p ):
q = p - r
q2 = np.dot( q, q )
mq = np.dot( m, q )
d2 = q2 - mq**2
return np.sqrt( d2 )
"""
residual function
no weighting introduced here, but if different lines or points
have different weights, that should be easy to change
"""
def res( params, line1, line2, line3 ):
x0, y0, a, b, g = params
r = np.array( [ x0, y0 ] )
m1 = np.array( [ np.cos( a ), np.sin( a ) ] )
m2 = np.array( [ np.cos( b ), np.sin( b ) ] )
m3 = np.array( [ np.cos( g ), np.sin( g ) ] )
dist1 = [ dist_lin_pnt( r, m1, p ) for p in line1 ]
dist2 = [ dist_lin_pnt( r, m2, p ) for p in line2 ]
dist3 = [ dist_lin_pnt( r, m3, p ) for p in line3 ]
return np.array( dist1 + dist2 + dist3 )
"""
Solution is quite stable w.r.t starting conditions
"""
solution = least_squares(
res,
x0=( 0, 0, 0, 0, 0 ),
args=( pl1, pl2, pl3 )
)
print( solution.x )
xf, yf, af, bf, gf = solution.x
pl1f = lin( tl, xf, yf, af )
pl2f = lin( tl, xf, yf, bf )
pl3f = lin( tl, xf, yf, gf )
"""
plotting the data
"""
fig = plt.figure()
ax = fig.add_subplot( 1, 1, 1 )
ax.plot( pl1[ ::, 0 ], pl1[ ::, 1 ], ls="", marker="o", color="r" )
ax.plot( pl2[ ::, 0 ], pl2[ ::, 1 ], ls="", marker="o", color="b" )
ax.plot( pl3[ ::, 0 ], pl3[ ::, 1 ], ls="", marker="o", color="y" )
ax.plot( pl1f[ ::, 0 ], pl1f[ ::, 1 ], color="r" )
ax.plot( pl2f[ ::, 0 ], pl2f[ ::, 1 ], color="b" )
ax.plot( pl3f[ ::, 0 ], pl3f[ ::, 1 ], color="y" )
plt.show()

最新更新