【问题标题】:How to fit a plane to a 3D dataset in Python如何在 Python 中将平面拟合到 3D 数据集
【发布时间】:2019-08-06 00:02:45
【问题描述】:

我有一组 x、y 和 z 点,并且正在尝试将平面拟合到这个 3D 数据,以便可以为任何 x 和 y 计算 z=f(x,y)。

我希望得到一个平面方程并在 Jupyter 笔记本中绘制图形以进行可视化。

这是我用来绘制数据的(工作)代码

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

x = np.arange(-12, 1)
y = np.arange(-40,-25) 
Z = array([[402., 398., 395., 391., 387., 383., 379., 375., 371., 367., 363.,358., 354.],
       [421., 417., 413., 409., 406., 402., 398., 393., 389., 385., 381.,
        376., 372.],
       [440., 436., 432., 429., 425., 421., 416., 412., 408., 404., 399.,
        395., 391.],
       [460., 456., 452., 448., 444., 440., 436., 432., 427., 423., 419.,
        414., 410.],
       [480., 477., 473., 469., 465., 460., 456., 452., 447., 443., 438.,
        434., 429.],
       [501., 498., 494., 490., 485., 481., 477., 472., 468., 463., 459.,
        454., 449.],
       [523., 519., 515., 511., 507., 502., 498., 494., 489., 484., 480.,
        475., 470.],
       [545., 541., 537., 533., 529., 524., 520., 515., 511., 506., 501.,
        496., 492.],
       [568., 564., 560., 556., 551., 547., 542., 538., 533., 528., 523.,
        518., 513.],
       [592., 588., 583., 579., 575., 570., 565., 561., 556., 551., 546.,
        541., 536.],
       [616., 612., 607., 603., 598., 594., 589., 584., 579., 575., 569.,
        564., 559.],
       [640., 636., 632., 627., 623., 618., 613., 609., 604., 599., 593.,
        588., 583.],
       [666., 662., 657., 653., 648., 643., 638., 633., 628., 623., 618.,
        613., 607.],
       [692., 688., 683., 679., 674., 669., 664., 659., 654., 649., 643.,
        638., 632.],
       [ nan, 714., 710., 705., 700., 695., 690., 685., 680., 675., 669.,
        664., 658.]])

X, Y = np.meshgrid(x, y)

fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')

print (X.shape, Y.shape, Z.shape)

ax.plot_surface(X, Y, Z)
ax.set_xlabel('X')
ax.set_ylabel('Y')
ax.set_zlabel('Z')

plt.show()

我已尝试实施这些解决方案:

https://gist.github.com/amroamroamro/1db8d69b4b65e8bc66a6

http://inversionlabs.com/2016/03/21/best-fit-surfaces-for-3-dimensional-data.html

但是,由于我的 x 和 y 数组的长度不同,我收到以下错误消息:

ValueError: all the input array dimensions except for the concatenation axis must match exactly

非常欢迎任何帮助。谢谢!

【问题讨论】:

    标签: python numpy matplotlib jupyter-notebook


    【解决方案1】:

    您分享的数据在绘图过程中似乎对我有用。您的 X, Y, Z 都具有相同的大小。您的 Z 数组中有一个 nan 值。您可以在估计平面方程时删除该点。

    您希望将数据拟合到 3D 计划中。因此,这是一个线性回归问题。您可以使用 scikit-learn 包中的多元回归来估计平面方程的系数。

    平面方程如下:

    Z = a1 * X + a2 * Y + c

    您可以按如下方式展平数据并使用scikit-learnlinear_model 将平面拟合到数据。请参考以下:

    # your data is stored as X, Y, Z
    print(X.shape, Y.shape, Z.shape)
    
    x1, y1, z1 = X.flatten(), Y.flatten(), Z.flatten()
    
    X_data = np.array([x1, y1]).reshape((-1, 2))
    Y_data = z1
    
    from sklearn import linear_model
    
    reg = linear_model.LinearRegression().fit(X_data, Y_data)
    
    print("coefficients of equation of plane, (a1, a2): ", reg.coef_)
    
    print("value of intercept, c:", reg.intercept_)
    

    上面的代码将把一个平面拟合到给定的线性数据。

    要适合二度曲面,请进一步阅读。

    您将获得以下形式的二度曲面方程:

    Z = a1*X + a2*Y + a3*X*Y + a4*X*X + a5*Y*Y + c

    要使用线性回归拟合这条曲线,您必须按以下方式修改上述代码:

    # your data is stored as X, Y, Z
    print(X.shape, Y.shape, Z.shape)
    
    x1, y1, z1 = X.flatten(), Y.flatten(), Z.flatten()
    x1y1, x1x1, y1y1 = x1*y1, x1*x1, y1*y1
    
    X_data = np.array([x1, y1, x1y1, x1x1, y1y1]).T  # X_data shape: n, 5
    Y_data = z1
    
    from sklearn import linear_model
    
    reg = linear_model.LinearRegression().fit(X_data, Y_data)
    
    print("coefficients of equation of plane, (a1, a2, a3, a4, a5): ", reg.coef_)
    
    print("value of intercept, c:", reg.intercept_)
    

    【讨论】:

    • 非常感谢您的帮助 Aditya。您的解决方案有效,但 Z 存在较大误差(至少 50%)。如何更改代码以适应多项式曲面而不是平面?
    • @rad189 我已更新代码以包含二阶曲线拟合。我希望它能解决你的问题。
    • 非常感谢阿迪亚!很好的解决方案
    猜你喜欢
    • 2019-01-17
    • 1970-01-01
    • 1970-01-01
    • 2019-10-19
    • 2013-03-07
    • 2017-01-02
    • 2020-01-24
    • 1970-01-01
    • 1970-01-01
    相关资源
    最近更新 更多