曲面拟合图

曲面拟合用于二维变量与目标函数之间关系的分析；

曲面拟合：

根据实际实验测试数据，求取函数f(x,y)与变量x及y之间的解析式，使其通过或近似通过所有的实验测试点；即使所有实验数据点能近似地分布在函数f(x,y)所表示的空间曲面上；

曲面拟合通常采用方式：

• 插值方式

• 逼近方式

共同点：

均利用曲面上或接近曲面的一组离散点，寻求良好的曲面方程；

区别:

• 插值方式得到的方程，所表示的曲面全部通过这组数据点，例如：lowess曲面拟合；

• 逼近方式只要求在某种准则下，其方程表示的曲面与这组数据点接近即可，例如：多项式曲面拟合；

• 逼近方式一般使用最小二乘法实现；

• 最小二乘法是一种逼近理论，也是采样数据进行拟合时最常用的一种方法；曲面一般不通过已知数据点，而是根据拟合的曲面在取样处的数值与实际值之差的平均和达到最小时求得；使拟合数值与实际数值之间的偏平方差的和达到最小；

绘制曲面拟合图

1. Statsmodels包的ols函数可以实现多元多次曲面拟合；

2. 先使用现有的数据集拟合得到多项式方程z=f(x,y)；

3. 再使用np.meshgird()函数生成x和y的网格数据；

4. 然后使用拟合的多项式预测z数值；

5. 最后使用ax.plot_surface()函数绘制拟合曲面；

一次曲面拟合

from statsmodels.formula.api import ols

import pandas as pd

from mpl_toolkits import mplot3d

import numpy as np

import matplotlib.pyplot as plt

import matplotlib as mpl

df=pd.read_csv('d:\python\out\SurfaceD.csv')

##多项式拟合z=f(x, y)=a+bx+cy+dx2+ey2

formula = 'z~x+y'

est = ols(formula,data=df).fit()

print(est.summary())

N=30

xmar= np.linspace(min(df.x),max(df.x),N)

ymar= np.linspace(min(df.y),max(df.y),N)

X,Y=np.meshgrid(xmar,ymar)

df_grid =pd.DataFrame({'x':X.flatten(),'y':Y.flatten()})

Z=est.predict(df_grid)

fig = plt.figure(figsize=(10,8),dpi =90)

ax = fig.gca(projection='3d')

#ax.set_aspect('equal','box')

ax.view_init(azim=60, elev=20)

##改变绘制图像的视角,即相机的位置,azim沿着z轴旋转，elev沿着y轴

#----------------------------------------------------------------

ax.grid(False)

ax.xaxis._axinfo['tick']['outward_factor'] = 0

ax.xaxis._axinfo['tick']['inward_factor'] = 0.4

ax.yaxis._axinfo['tick']['outward_factor'] = 0

ax.yaxis._axinfo['tick']['inward_factor'] = 0.4

ax.zaxis._axinfo['tick']['outward_factor'] = 0

ax.zaxis._axinfo['tick']['inward_factor'] = 0.4

# =================================================================

ax.xaxis.pane.fill = False

ax.yaxis.pane.fill = False

ax.zaxis.pane.fill = False

# =================================================================

ax.xaxis.pane.set_edgecolor('k')

ax.yaxis.pane.set_edgecolor('k')

ax.zaxis.pane.set_edgecolor('k')

p=ax.plot_surface(X,Y, Z.values.reshape(N,N), rstride=1, cstride=1, cmap='Spectral_r',

alpha=1,edgecolor='k',linewidth=0.25)

ax.set_xlabel( "Gax Mileage (mpg)")

ax.set_ylabel("0-60 mph (sec)")

ax.set_zlabel("Power (KW)")

ax.set_zlim(20,140)

cbar=fig.colorbar(p, shrink=0.5,aspect=10)

cbar.set_label('Power (kW)')

fig.savefig('d:\python\out\三维曲面图3.pdf')

OLS Regression Results

=================================================================

Dep. Variable: z R-squared: 0.614

Model: OLS Adj. R-squared: 0.611

Method: Least Squares F-statistic: 267.7

Date: Wed, 19 Aug 2020 Prob (F-statistic): 2.44e-70

Time: 21:03:14 Log-Likelihood: -1454.1

No. Observations: 340 AIC: 2914.

Df Residuals: 337 BIC: 2926.

Df Model: 2

Covariance Type: nonrobust

=================================================================

coef std err t P>|t| [0.025 0.975]

-----------------------------------------------------------------

Intercept 164.2697 5.318 30.888 0.000 153.809 174.731

x -0.2658 0.314 -0.846 0.398 -0.884 0.352

y -3.6913 0.168 -21.974 0.000 -4.022 -3.361

=================================================================

Omnibus: 34.307 Durbin-Watson: 1.457

Prob(Omnibus): 0.000 Jarque-Bera (JB): 55.148

Skew: 0.637 Prob(JB): 1.06e-12

Kurtosis: 4.507 Cond. No. 153.

=================================================================

Warnings:

[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.

二次曲面拟合

from statsmodels.formula.api import ols

import pandas as pd

import numpy as np

import matplotlib.pyplot as plt

df=pd.read_csv('d:\python\out\SurfaceD.csv')

##多项式拟合z=f(x, y)=a+bx+cy+dx2+ey2

formula = 'z~x+np.square(x)+y+np.square(y)'

est = ols(formula,data=df).fit()

print(est.summary())

N=30

xmar= np.linspace(min(df.x),max(df.x),N)

ymar= np.linspace(min(df.y),max(df.y),N)

X,Y=np.meshgrid(xmar,ymar)

df_grid =pd.DataFrame({'x':X.flatten(),'y':Y.flatten()})

Z=est.predict(df_grid)

fig = plt.figure(figsize=(10,8),dpi =90)

ax = fig.gca(projection='3d')

#ax = fig.add_subplot(projection='3d')

#ax.set_aspect('equal','box')

ax.view_init(azim=60, elev=20)

##改变绘制图像的视角,即相机的位置,azim沿着z轴旋转，elev沿着y轴

ax.grid(False)

ax.xaxis._axinfo['tick']['outward_factor'] = 0

ax.xaxis._axinfo['tick']['inward_factor'] = 0.4

ax.yaxis._axinfo['tick']['outward_factor'] = 0

ax.yaxis._axinfo['tick']['inward_factor'] = 0.4

ax.zaxis._axinfo['tick']['outward_factor'] = 0

ax.zaxis._axinfo['tick']['inward_factor'] = 0.4

ax.xaxis.pane.fill = False

ax.yaxis.pane.fill = False

ax.zaxis.pane.fill = False

ax.xaxis.pane.set_edgecolor('k')

ax.yaxis.pane.set_edgecolor('k')

ax.zaxis.pane.set_edgecolor('k')

p=ax.plot_surface(X,Y, Z.values.reshape(N,N), rstride=1, cstride=1, cmap='Spectral_r',

alpha=1,edgecolor='k',linewidth=0.25)

ax.set_xlabel( "Gax Mileage (mpg)")

ax.set_ylabel("0-60 mph (sec)")

ax.set_zlabel("Power (KW)")

ax.set_zlim(50,170)

cbar=fig.colorbar(p, shrink=0.5,aspect=10)

cbar.set_label('Power (kW)')

fig.savefig('d:\python\out\三维曲面图2.pdf')

OLS Regression Results

=================================================================

Dep. Variable: z R-squared: 0.701

Model: OLS Adj. R-squared: 0.697

Method: Least Squares F-statistic: 195.9

Date: Wed, 19 Aug 2020 Prob (F-statistic): 2.28e-86

Time: 21:01:42 Log-Likelihood: -1410.8

No. Observations: 340 AIC: 2832.

Df Residuals: 335 BIC: 2851.

Df Model: 4

Covariance Type: nonrobust

=================================================================

coef std err t P>|t| [0.025 0.975]

-----------------------------------------------------------------

Intercept 290.7007 18.804 15.460 0.000 253.713 327.689

x -4.2189 2.123 -1.987 0.048 -8.396 -0.042

np.square(x) 0.1314 0.064 2.052 0.041 0.005 0.257

y -12.9581 1.002 -12.936 0.000 -14.929 -10.988

np.square(y) 0.2039 0.022 9.395 0.000 0.161 0.247

================================================================

Omnibus: 45.663 Durbin-Watson: 1.814

Prob(Omnibus): 0.000 Jarque-Bera (JB): 123.583

Skew: 0.622 Prob(JB): 1.46e-27

Kurtosis: 5.679 Cond. No. 1.41e+04

=================================================================

Warnings:

[1] Standard Errors assume that the covariance matrix of the errors is correctly specified.

[2] The condition number is large, 1.41e+04. This might indicate that there are

strong multicollinearity or other numerical problems.