利用Python中的sklearn函数库的LinearRegression和PolynomialFeatures进行函数拟合

具体程序如下：

import matplotlib.pyplot as plt
import pandas as pd 
import numpy as np 
from sklearn.linear_model import LinearRegression #导入线性回归模块
from sklearn.preprocessing import PolynomialFeatures


x=[[-12,	2.54],[-12,	2.19],[-20,	1.75],[-20,	1.85],[-5,	3.58],[-7,	2.44],[7,		3.20],[25,	2.57],[12,	1.77],[20,	2.77],[15,	3.25],[-20,	2.34],[-30,	1.57],[-20,	1.52],[15,	2.57],[10,	3.75],[15,	1.26],[-12,	2.54],[-5.87,	2.57],[0.08,	2.88],[-11.9,	2.32],[6.1,	1.55],[-0.15,	3.16],[-5.26,	3.82],[-10.29,1.07], [-39.92,1.90],]

y = [41,50,50,50,50,48,44,60,46,60,55,35,35,45,50,60,40,41,31.02,20.98,30.97,11.05,31.16,20.92,21.18,31.09,]
print("==================================================================")
for index in range(1,100):
	data=pd.DataFrame({'IN':x, 'OUT':y})
	data_train=np.array(data['IN']).reshape(data['IN'].shape[0],1)				
	data_test=data['OUT']
	
	poly_reg =PolynomialFeatures(degree = index) 
	X_ploy =poly_reg.fit_transform(x)
	regr=LinearRegression() 													
	regr.fit(X_ploy,data_test)													
	print("vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv")
	print("degree = ", index)
	print("coefficients = ", regr.coef_)
	print("intercept = ", regr.intercept_)
	print("R^2 = ",regr.score(X_ploy,data_test))								
	print("^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^")
	if(regr.score(X_ploy,data_test) >= 0.99):
		break
print("==================================================================")

计算结果如下：

==================================================================
vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv
degree =  1
coefficients =  [0.         0.23417502 2.49000913]
intercept =  35.655517534736624
R^2 =  0.13062845852051208
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv
degree =  2
coefficients =  [ 0.00000000e+00 -8.67998076e-01  2.85867915e+01  5.05996188e-03
  5.41937770e-01 -5.16031113e+00]
intercept =  2.8917936495706513
R^2 =  0.36680282974834966
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv
degree =  3
coefficients =  [ 0.00000000e+00  2.32337035e+00  8.62837626e+01  7.32845678e-02
 -2.76400095e+00 -3.16485345e+01  1.25383404e-03 -1.86885458e-02
  7.02626264e-01  3.89746542e+00]
intercept =  -43.565823068235204
R^2 =  0.595036736861785
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv
degree =  4
coefficients =  [ 0.00000000e+00 -9.27359226e+00  7.50519878e+02  7.95413552e-01
  8.23353562e+00 -4.11409690e+02  1.71490391e-02 -6.59155370e-01
 -2.36676840e+00  9.91069353e+01  5.28441826e-05 -7.07256550e-03
  1.39851401e-01  2.41607515e-01 -8.88958626e+00]
intercept =  -469.6580685131609
R^2 =  0.6830703676187968
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv
degree =  5
coefficients =  [ 0.00000000e+00  9.28501268e+01 -1.34960520e+04 -4.52985493e+00
 -9.95780624e+01  1.14314413e+04 -2.77506890e-01  8.06281464e+00
  3.22856832e+01 -4.68204475e+03 -4.20565631e-03  2.50969496e-01
 -4.00613710e+00 -2.50785085e+00  9.33885800e+02 -2.17313708e-05
  1.67482250e-03 -5.36753698e-02  6.08933509e-01 -2.07588174e-01
 -7.28846679e+01]
intercept =  6113.3609609005725
R^2 =  0.9155071848237998
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
vvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvvv
degree =  6
coefficients =  [ 1.49368226e-01  1.32073179e+03 -1.45048362e+02 -2.26176086e+02
 -1.25241633e+03 -7.13337162e+02 -1.02801495e+01  5.40473325e+02
 -6.11134917e+02 -6.92340948e+02 -2.88372238e-01  1.72578578e+01
 -4.38748168e+02  9.77187097e+02  1.17156755e+03 -3.37666499e-03
  2.83940755e-01 -8.79360290e+00  1.47098644e+02 -3.39152229e+02
 -4.67229808e+02 -1.24293556e-05  1.53024369e-03 -6.66028716e-02
  1.40134218e+00 -1.74726208e+01  3.71705693e+01  5.87943145e+01]
intercept =  1316.821170263324
R^2 =  0.9999999999995062
^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
==================================================================

在多项式最高阶数达到6阶后，其R^2超过99%。

生成多项式如下：

二元二阶：

X1：

X2：

 生成矩阵为：

则多元多项式为：

从多项式可以看出，整个矩阵没有全部使用，具体规则如下：

1         2         4

3         5         0

6         0         0

其中1标识的是，其他的以此类推。0代表没有使用。

二元三阶的也按同样的方式即可，如下：

1        2        4        7

3        5        8        0

6        9        0        0

10        0        0        0

（一）使用matlab画多项式函数图像

在用matlab画五次多项式的时候，发现使用y=(0.2771*t^3 - 0.081*t^4 + 0.0063*t^5);总是会报错。搜了搜资料，发现可以这么画：

例如：五次多项式（quintic polynomial） y=0.2771t³-0.081t⁴+0.0063t⁵

.m脚本文件

t=0:0.01:5;
p=[0.0063,-0.081,0.2771,0,0,0];
y=polyval(p,t);

figure(1)
plot(t,y,'r','linewidth',2)
grid on
hold on;

（二）使用python画多项式函数图像

例如：五次多项式（quintic polynomial） y=0.2771t³-0.081t⁴+0.0063t⁵
.py脚本文件

# coding=utf-8
import numpy as np
import matplotlib.pyplot as plt

plt.rcParams['font.sans-serif'] = ['Arial']   # 如果要显示中文字体,则在此处设为：SimHei
plt.rcParams['axes.unicode_minus'] = False    # 显示负号

x = np.linspace(0,5)
y = 0.2771*x*x*x-0.081*x*x*x*x+0.0063*x*x*x*x*x

plt.figure(figsize=(10, 5))
plt.grid(linestyle="--")                      # 设置背景网格线为虚线
ax = plt.gca()
ax.spines['top'].set_visible(False)           # 去掉上边框
ax.spines['right'].set_visible(False)         # 去掉右边框

plt.plot(x, y, marker='o', color="blue",  label="V0=10m/s", linewidth=1.5)
plt.yticks(fontsize=12, fontweight='bold')
plt.title("Quintic Polynomial Graph", fontsize=12, fontweight='bold')
plt.xlabel("Time/s", fontsize=13, fontweight='bold')
plt.ylabel("Vehicle lateral distance/m", fontsize=13, fontweight='bold')

plt.legend(loc=0, numpoints=1)
leg = plt.gca().get_legend()
ltext = leg.get_texts()
plt.setp(ltext, fontsize=12, fontweight='bold')  # 设置图例字体的大小和粗细
plt.show()

更新：2020年4月1日
偶然间看到matlab还有一种画多项式函数的方法，而且就是之前报错的那种方法。
.m脚本

t = 0:0.01:5;
y = 0.2771 * t.^3 - 0.081 * t.^4 + 0.0063 * t.^5;  %注意.和^之间没有空格
plot(t, y, 'r')

数据分析中经常会使用到数据拟合，本文中将阐述如何实现一元以及多元的线性拟合以及多项式拟合，本文中只涉及实现方式，不涉及理论知识。

模型拟合中涉及的误差评估方法如下所示：

import numpy as np
def stdError_func(y_test, y):
  return np.sqrt(np.mean((y_test - y) ** 2))


def R2_1_func(y_test, y):
  return 1 - ((y_test - y) ** 2).sum() / ((y.mean() - y) ** 2).sum()


def R2_2_func(y_test, y):
  y_mean = np.array(y)
  y_mean[:] = y.mean()
  return 1 - stdError_func(y_test, y) / stdError_func(y_mean, y)

一元线性拟合

import pandas as pd
import numpy as np
from sklearn.preprocessing import PolynomialFeatures
from sklearn import linear_model


filename = "E:/data.csv"
df= pd.read_csv(filename)
x = np.array(df.iloc[:,0].values)

y = np.array(df.iloc[:,5].values)

cft = linear_model.LinearRegression()
cft.fit(x[:,np.newaxis], y) #模型将x变成二维的形式, 输入的x的维度为[None, 1]

print("model coefficients", cft.coef_)
print("model intercept", cft.intercept_)


predict_y =  cft.predict(x[:,np.newaxis])
strError = stdError_func(predict_y, y)
R2_1 = R2_1_func(predict_y, y)
R2_2 = R2_2_func(predict_y, y)
score = cft.score(x[:,np.newaxis], y) ##sklearn中自带的模型评估，与R2_1逻辑相同

print(' strError={:.2f}, R2_1={:.2f},  R2_2={:.2f}, clf.score={:.2f}'.format(
    strError,R2_1,R2_2,score))

结果输出为：
model coefficients [-31.2375]
model intercept 7.415750000000001
strError=1.11, R2_1=0.28, R2_2=0.15, clf.score=0.28

模型拟合的表达式为：
y = 7.415750000000001 +(-31.2375) * x
从拟合的均方误差和得分来看效果不佳

多元线性拟合

 import pandas as pd
import numpy as np
from sklearn.preprocessing import PolynomialFeatures
from sklearn import linear_model




filename = "E:/data.csv"
df= pd.read_csv(filename)
x = np.array(df.iloc[:,0:4].values)

y = np.array(df.iloc[:,5].values)

cft = linear_model.LinearRegression()
print(x.shape)
cft.fit(x, y) #

print("model coefficients", cft.coef_)
print("model intercept", cft.intercept_)


predict_y =  cft.predict(x)
strError = stdError_func(predict_y, y)
R2_1 = R2_1_func(predict_y, y)
R2_2 = R2_2_func(predict_y, y)
score = cft.score(x, y) ##sklearn中自带的模型评估，与R2_1逻辑相同

print('strError={:.2f}, R2_1={:.2f},  R2_2={:.2f}, clf.score={:.2f}'.format(
    strError,R2_1,R2_2,score))

结果输出为：
model coefficients [-31.2375 17.74375 44.325 5.7375 ]
model intercept 0.5051249999999978
strError=0.58, R2_1=0.80, R2_2=0.56, clf.score=0.80

模型拟合的表达式为：
y = 0.5051249999999978 +(-31.2375) * x11 + 17.74375 *x2 + 44.325 * x3 + 5.7375 * x4
从拟合的均方误差和得分来看在之前的基础上有所提升

一元多项式拟合

以三次多项式为例

import pandas as pd
import numpy as np
from sklearn.preprocessing import PolynomialFeatures
from sklearn import linear_model

filename = "E:/data.csv"
df= pd.read_csv(filename)
x = np.array(df.iloc[:,0].values)

y = np.array(df.iloc[:,5].values)

poly_reg =PolynomialFeatures(degree=3) #三次多项式
X_ploy =poly_reg.fit_transform(x[:, np.newaxis])
print(X_ploy.shape)
lin_reg_2=linear_model.LinearRegression()
lin_reg_2.fit(X_ploy,y)
predict_y =  lin_reg_2.predict(X_ploy)
strError = stdError_func(predict_y, y)
R2_1 = R2_1_func(predict_y, y)
R2_2 = R2_2_func(predict_y, y)
score = lin_reg_2.score(X_ploy, y) ##sklearn中自带的模型评估，与R2_1逻辑相同

print("model coefficients", lin_reg_2.coef_)
print("model intercept", lin_reg_2.intercept_)
print('degree={}: strError={:.2f}, R2_1={:.2f},  R2_2={:.2f}, clf.score={:.2f}'.format(
    3, strError,R2_1,R2_2,score))

输出结果
model coefficients [ 0. 990.64583333 -11906.25 44635.41666667]
model intercept -20.724999999999117
degree=3: strError=1.08, R2_1=0.32, R2_2=0.17, clf.score=0.32

对应的函数表达式为 -20.724999999999117 + [0, 990.64583333, -11906.25, 44635.41666667] *[1， x, x^2, x.^3].T = -20.724999999999117 + 990.64583333 * x + ( -11906.25) * x^2 + 44635.41666667 * x^3

多元多项式拟合

import pandas as pd
import numpy as np
from sklearn.preprocessing import PolynomialFeatures
from sklearn import linear_model

filename = "E:/data.csv"
df= pd.read_csv(filename)
x = np.array(df.iloc[:,0:4].values)

y = np.array(df.iloc[:,5].values)


poly_reg =PolynomialFeatures(degree=2) #三次多项式
X_ploy =poly_reg.fit_transform(x)
lin_reg_2=linear_model.LinearRegression()
lin_reg_2.fit(X_ploy,y)
predict_y =  lin_reg_2.predict(X_ploy)
strError = stdError_func(predict_y, y)
R2_1 = R2_1_func(predict_y, y)
R2_2 = R2_2_func(predict_y, y)
score = lin_reg_2.score(X_ploy, y) ##sklearn中自带的模型评估，与R2_1逻辑相同

print("coefficients", lin_reg_2.coef_)
print("intercept", lin_reg_2.intercept_)
print('degree={}: strError={:.2f}, R2_1={:.2f},  R2_2={:.2f}, clf.score={:.2f}'.format(
    3, strError,R2_1,R2_2,score))

函数输出结果为：

coefficients [ 0. 332.28129937 -19.9240981 -9.10607925
-191.05593023 -287.93919929 -912.11402936 -1230.21922184
-207.90033986 99.03441748 190.26204994 433.25169929
273.13674555 257.66550523 344.92652936]
intercept 4.35175537840722
degree=3: strError=0.23, R2_1=0.97, R2_2=0.82, clf.score=0.97

代码中输入的自变量是一个包含四个变量的输入， 对应coefficients输出的是长度为15的向量， 其中对应到的变量分别为 variable_X = [1, x1, x2, x3, x4, $x1\ast x1$, $x1\ast x2$, $x1\ast x3$, $x1\ast x4$, $x2\ast x2$, $x2\ast x3$, $x2\ast x4$, $x3\ast x3$, $x3\ast x4$, $x4\ast x4$]

对应的方程式为： $intercept+coefficients\ast variabl{e}_X.T$

代码中涉及到的数据集如下：

a,b,c,d,e
0.06,0.2,0.02,0.1,0.340 
0.1,0.28,0.02,0.14,0.370 
0.12,0.32,0.02,0.16,0.377 
0.08,0.24,0.02,0.12,0.383 
0.08,0.32,0.04,0.1,0.383 
0.12,0.28,0.03,0.1,0.393 
0.1,0.24,0.05,0.1,0.385 
0.06,0.32,0.05,0.14,0.362 
0.12,0.2,0.05,0.12,0.320 
0.06,0.28,0.04,0.12,0.393 
0.08,0.28,0.05,0.16,0.402 
0.08,0.2,0.03,0.14,0.349 
0.1,0.2,0.04,0.16,0.335 
0.1,0.32,0.03,0.12,0.387 
0.12,0.24,0.04,0.14,0.390 
0.06,0.24,0.03,0.16,0.315 

指数函数和幂函数拟合参照网址：
https://blog.csdn.net/kl28978113/article/details/88818885
参考链接：
https://blog.csdn.
net/weixin_44794704/article/details/89246032
https://blog.csdn.net/bxg1065283526/article/details/80043049
https://www.cnblogs.com/Lin-Yi/p/8975638.html