通过 3D x,y,z 散点图数据拟合一条线答案

【问题标题】：Fitting a line through 3D x,y,z scatter plot data通过 3D x,y,z 散点图数据拟合一条线
【发布时间】：2021-04-29 11:33:48
【问题描述】：

我有一些数据点在 3d 空间中沿着一条线聚集。我在要导入的 csv 文件中有 x、y、z 数据。我想找到一个方程来代表那条线，或者垂直于那条线的平面，或者任何数学上正确的东西。这些数据相互独立。也许有比我尝试做的更好的方法来做到这一点，但是......

我试图在这里复制一篇似乎正在做我想做的事情的旧帖子 Fitting a line in 3D

但似乎过去十年的更新导致代码的第二部分无法正常工作？或者，也许我只是做错了什么。我已经在底部包含了我从这个 frankenstein 组合在一起的全部内容。有两行似乎给我带来了问题。

我已经把它们偷偷跑出来了……

import numpy as np

pts = np.add.accumulate(np.random.random((10,3)))
x,y,z = pts.T

# this will find the slope and x-intercept of a plane
# parallel to the y-axis that best fits the data
A_xz = np.vstack((x, np.ones(len(x)))).T
m_xz, c_xz = np.linalg.lstsq(A_xz, z)[0]

# again for a plane parallel to the x-axis
A_yz = np.vstack((y, np.ones(len(y)))).T
m_yz, c_yz = np.linalg.lstsq(A_yz, z)[0]

# the intersection of those two planes and
# the function for the line would be:
# z = m_yz * y + c_yz
# z = m_xz * x + c_xz
# or:
def lin(z):
    x = (z - c_xz)/m_xz
    y = (z - c_yz)/m_yz
    return x,y

#verifying:
from mpl_toolkits.mplot3d import Axes3D
import matplotlib.pyplot as plt

fig = plt.figure()
ax = Axes3D(fig)
zz = np.linspace(0,5)
xx,yy = lin(zz)
ax.scatter(x, y, z)
ax.plot(xx,yy,zz)
plt.savefig('test.png')
plt.show()

他们返回这个，但没有值......

FutureWarning：rcond 参数将更改为机器精度时间的默认值max(M, N)，其中 M 和 N 是输入矩阵维度。要使用未来的默认值并消除此警告，我们建议传递rcond=None，继续使用旧的，明确传递rcond=-1。 m_xz, c_xz = np.linalg.lstsq(A_xz, z)[0] FutureWarning：rcond 参数将更改为默认的机器精度时间 max(M, N) 其中 M 和 N 是输入矩阵维度。要使用未来的默认值并消除此警告，我们建议传递rcond=None，继续使用旧的，明确传递rcond=-1。 m_yz, c_yz = np.linalg.lstsq(A_yz, z)[0]

我不知道从这里去哪里。我什至实际上不需要情节，我只需要一个方程式并且没有能力继续前进。如果有人知道更简单的方法，或者可以指出正确的方向，我愿意学习，但我非常非常迷茫。提前谢谢您！

这是我的整个 frankensteined 代码，以防万一导致问题。

import pandas as pd
import numpy as np
mydataset = pd.read_csv('line1.csv')

x = mydataset.iloc[:,0]
y = mydataset.iloc[:,1]
z = mydataset.iloc[:,2]


data = np.concatenate((x[:, np.newaxis], 
                       y[:, np.newaxis], 
                       z[:, np.newaxis]), 
                      axis=1)


# Calculate the mean of the points, i.e. the 'center' of the cloud
datamean = data.mean(axis=0)

# Do an SVD on the mean-centered data.
uu, dd, vv = np.linalg.svd(data - datamean)

# Now vv[0] contains the first principal component, i.e. the direction
# vector of the 'best fit' line in the least squares sense.

# Now generate some points along this best fit line, for plotting.

# we want it to have mean 0 (like the points we did
# the svd on). Also, it's a straight line, so we only need 2 points.
linepts = vv[0] * np.mgrid[-100:100:2j][:, np.newaxis]

# shift by the mean to get the line in the right place
linepts += datamean

# Verify that everything looks right.

import matplotlib.pyplot as plt
import mpl_toolkits.mplot3d as m3d

ax = m3d.Axes3D(plt.figure())
ax.scatter3D(*data.T)
ax.plot3D(*linepts.T)
plt.show()

# this will find the slope and x-intercept of a plane
# parallel to the y-axis that best fits the data
A_xz = np.vstack((x, np.ones(len(x)))).T
m_xz, c_xz = np.linalg.lstsq(A_xz, z)[0]

# again for a plane parallel to the x-axis
A_yz = np.vstack((y, np.ones(len(y)))).T
m_yz, c_yz = np.linalg.lstsq(A_yz, z)[0]

# the intersection of those two planes and
# the function for the line would be:
# z = m_yz * y + c_yz
# z = m_xz * x + c_xz
# or:
def lin(z):
    x = (z - c_xz)/m_xz
    y = (z - c_yz)/m_yz
    return x,y

print(x,y)

#verifying:
from mpl_toolkits.mplot3d import Axes3D
import matplotlib.pyplot as plt

fig = plt.figure()
ax = Axes3D(fig)
zz = np.linspace(0,5)
xx,yy = lin(zz)
ax.scatter(x, y, z)
ax.plot(xx,yy,zz)
plt.savefig('test.png')
plt.show()

【问题讨论】：

你能给我们一个你的数据集的sn-p吗？尝试使用 mydataset.head(20) 将前 20 行添加到您的问题中
您只是想查看数据吗？我这里有 csv 文件。现在只有 17 个点，但我在不同的线上有几百个点会来来去去，这就是我希望它与 excel 兼容的原因。谢谢你的耐心，哈哈。 drive.google.com/file/d/1I0SkBnpzFTHS7jMWSotegVsjl-pxH2bO/…
@ShelbyJohnston，您应该将数据作为列表或 numpy 数组添加到代码 sn-p 中，这样我们就可以在不下载任何数据的情况下运行您的代码。
我将其保留为导入 csv，因为我认为这可能是问题所在。现在我看到它不是。我复制了原始代码，因为我知道它没有错（而且我不知道自己这样做会搞砸什么）并将其粘贴到 sn-p.xml 中。他们只是使用随机值，它给出了完全相同的错误......

标签： python pandas numpy matplotlib linear-algebra

【解决方案1】：

您可以通过添加rcond=None 来消除来自leastsquares 的投诉，如下所示：

m_xz, c_xz = np.linalg.lstsq(A_xz, z, rcond=None)[0]

这是适合您情况的正确决定吗？我不知道。但在docs 中还有更多相关信息。

当我使用您的输入运行您的代码时，它似乎运行得很好，并且我得到分配给 m_xz、c_xz 等的值。如果您不使用 print('m_xz')（或其他）显式调用它们，那么你不会看到他们。

m_xz
Out[42]: 5.186132604596112

c_xz
Out[43]: 62.5764694106141

此外，您以两种不同的方式引用您的数据。您从 csv 中获取 x、y 和 z，但也将其放入 numpy 数组中。您只需使用 numpy 就可以摆脱重复和熊猫：

data = np.genfromtxt('line1.csv', delimiter=',', skip_header=1)

x = data[:,0]
y = data[:,1]
z = data[:,2]

【讨论】：

【解决方案2】：

正如old post you refer to 中所建议的，您也可以使用主成分分析来代替最小二乘法。为此，我建议 sklearn.decomposition.PCA 来自 sklearn package。

使用您提供的 csv 文件可以在下面找到一个示例。

import pandas as pd
import numpy as np
from sklearn.decomposition import PCA
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D

mydataset = pd.read_csv('line1.csv')

x = mydataset.iloc[:,0]
y = mydataset.iloc[:,1]
z = mydataset.iloc[:,2]

coords = np.array((x, y, z)).T

pca = PCA(n_components=1)
pca.fit(coords)
direction_vector = pca.components_
print(direction_vector)


# Create plot
origin = np.mean(coords, axis=0)
euclidian_distance = np.linalg.norm(coords - origin, axis=1)
extent = np.max(euclidian_distance)

line = np.vstack((origin - direction_vector * extent,
                  origin + direction_vector * extent))

fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
ax.scatter(coords[:, 0], coords[:, 1], coords[:,2])
ax.plot(line[:, 0], line[:, 1], line[:, 2], 'r')

【讨论】：