如果数据比简单的直线更为复杂，我们也可以用线性模型来拟合非线性数据。一个简单的方法就是将每一个特征的幂次方添加为一个新的特征，然后在这个拓展的特征集上进行线性拟合，这种方法成为多项式回归。

回归分析的目标是根据自变量（或自变量向量）x 的值来模拟因变量 y 的期望值。在简单的线性回归中，使用模型

\({y=β_{0}+β_{1}x+\varepsilon }\)

其中ε是未观察到的随机误差，其以标量 x 为条件，均值为零。在该模型中，对于 x 值的每个单位增加，y 的条件期望增加 \({β_{1}} \)个单位。

在许多情况下，这种线性关系可能不成立。例如，如果我们根据合成发生的温度对化学合成的产率进行建模，我们可以发现通过增加每单位温度增加的量来提高产率。在这种情况下，我们可能会提出如下所示的二次模型：

\({y=β _{0}+β_{1}x+β_{2}x^{2}+\varepsilon }\)

通常，我们可以将 y 的期望值建模为 n 次多项式，得到一般多项式回归模型：

\({y=β _{0}+β_{1}x+β_{2}x^{2}+β_{3}x^{3}+\cdots +β_{n}x^{n}+\varepsilon }\)

为了方便，这些模型从估计的角度来看都是线性的，因为回归函数就未知参数\(β_0\)、\(β_1\)等而言是线性的。因此，对于最小二乘分析，多项式回归的计算和推理问题可以使用多元回归技术完全解决，这是通过将 \(x\)、\(x^2\) 等视为多元回归模型中的独特自变量来完成的。

一、单个特征值

1.1 生成演示数据

我们举一个简单的例子，我们生成一些符合\(y = 0.5x^2+2x+2+\varepsilon \)数据（\(ε\)为高斯噪声）。

import numpy as np
import matplotlib.pyplot as plt

X = 10 * np.random.rand(100,1) - 5
X = np.sort(X,axis=0) #为了方便作图，X和y的每一个组合的顺序没有关系的
y = 0.5* X ** 2 + 2 * X + 2 + np.random.randn(100,1)
plt.scatter(X,y,c='green', alpha=0.6)

1.2 生成新的特征

可以看出来，直线是不能拟合这个数据的，我们用Scikit-Learn的PolynomialFeatures对特征进行转换，将每一个特征的平方作为新的特征加入训练集。

from sklearn.preprocessing import PolynomialFeatures
pf = PolynomialFeatures(degree = 2,include_bias = False)
#include_bias：默认为True。如果为True的话，那么就会有常量为1的那一项。
#interaction_only： 默认为False，如果指定为True，那么就不会有特征自己和自己结合的项。
X_new = pf.fit_transform(X)

我们看看新的特征X_new

第一列为原来的X，第二列是X的平方。

1.3 用新的特征进行回归

现在新的特征集X_new包含来原来的特征X和该特征的平方，系那种对这个拓展的训练集匹配一个线性回归模型。

from sklearn.linear_model import LinearRegression
reg = LinearRegression()
reg.fit(X_new,y)
reg.intercept_
reg.coef_

输出：

也就是我们拟合得到的模型是\(\hat{y} = 0.47x^2+1.95x+2.37\)，而我们原来的函数是\(y = 0.5x^2+2x+2+ε\)已经很不错的拟合了。

1.4 拟合效果

我们看看拟合的曲线。

plt.scatter(X,y,c='green', alpha=0.6)
plt.plot(X, reg.predict(X_new), color='r')
plt.show()

二、多个特征值

当存在多个特征值的时候，多项式回归能发现特征和特征之间的关系（线性回归做不到这个）。PolynomialFeatures会在给定的多项式的阶数下，添加所有的特征组合。如有两个特征a和b，阶数为3的情况下，PolynomialFeatures会添加\(a^2,b^2,a^3,b^3,ab,a^2b,ab^2\)

警告！！经过 PolynomialFeatures 之后，样本特征呈指数增长，新增的特征包含了所有可能的样式。

接下来我们演示有两个特征的多项式回归模型的拟合。

2.1 生成演示数据

首先我们模拟两个正态分布的随机特征，通过函数\(y = 2x_1+3x_2+4x_1^2+5x_1x_2+6x_2^2+7+ ε\)生成因变量y的数据，我们把X1和X2合并成X，并把X的里两个特征命名为x1和x2

X1 = 10 * np.random.rand(100,1) - 5
X2 = 2 * np.random.rand(100,1) - 9
y = 2 * X1 + 3 * X2 + 4 * X1**2 + 5 * X1*X2 + 6 * X2**2 + 7 + np.random.randn(100,1)
X = pd.DataFrame(np.c_[X1,X2])
X.columns=['x1','x2']

2.2 生成新的特征

针对x1和x2，生成包含\(x_1x_2,x_1^2,x_2^2\)特征的df数据框

pf = PolynomialFeatures(degree=2).fit(X)
X_new = pf.transform(X)
clumns_list = pf.get_feature_names(X.columns)
features = pd.DataFrame(X_new, columns=clumns_list)

X_new是这样的

获取X的列名

通过PolynomialFeatures的get_feature_names得到遍历的符合PolynomialFeatures顺序的列名。

我们通过pd的DataFrame转换为df数据框，列名为clumns_list，最终的数据如下图：

2.3 用新的特征进行回归

reg = LinearRegression()
reg.fit(features,y)
reg.intercept_
reg.coef_

回归得到的公式是\(\hat{y} = 1.78x_1+9.96x_2+4.01x_1^2+4.97x_1x_2+6.44x_2^2+34.41\)，
而我们原来的函数是：\(y = 2x_1+3x_2+4x_1^2+5x_1x_2+6x_2^2+7+ ε\)还是有点差距的，但这里仅仅是说明这样的方法，不代表方法一定非常准确。

!!多个特征我们无法在平面画出拟合图。

2.4 数据标准化

因为新的拓展的特征的指数导致新的特征的绝对值非常大，数据的差距也会变动非常的大，所以我们要用标准化后的数据进行回归。

from sklearn.preprocessing import StandardScaler

scaler = StandardScaler().fit(X_new)
X_new_std = scaler.transform(X_new)
features_std = pd.DataFrame(X_new_std, columns=clumns_list)
reg = LinearRegression()
reg.fit(features_std,y)
reg.intercept_

reg.coef_

如果用回归模型进行预测，输入的特征是包含拓展了的特征，并且也是标准化后的特征，如果有新的数据进入需要用拟合的方程进行预测，要把新的数据按同样的规则生成新的特征，并把对应的特征做标准化处理。对新的特征的标准化处理方法：(数值-scaler.mean_)/np.sqrt(scaler.var_)

具体结果不演示了。

三、Pipeline

为了简化这个过程，Scikit-Learn中的 Pipeline 将这三个模型封装起来串联操作，让模型接口更加简洁，使用起来方便的减少代码量同时让机器学习的流程变得直观和更加的优雅。函数 PolynomialRegression() 中传入的超参数 degree 是用来指定所得的多项式回归中所用多项式的阶次。

pipeline最后一步如果有predict()方法我们才可以对pipeline使用fit_predict()，同理，最后一步如果有transform()方法我们才可以对pipeline使用fit_transform()方法。

3.1 生成模拟数据

我们同样用会第一部分的模拟数据

X = 10 * np.random.rand(100,1) - 5
X = np.sort(X,axis=0) ##为了方便作图，X和y的每一个组合的顺序没有关系的
y = 0.5* X ** 2 + 2 * X + 2 + np.random.randn(100,1)

3.2 建立模型

用6阶转换以及标准化后用线性模型

poly_reg = Pipeline([
("poly", PolynomialFeatures(degree=6)),
("std_scaler", StandardScaler()),
("lin_reg", LinearRegression())
])

3.3 拟合曲线

poly_reg.fit(X, y)
y_predict = poly_reg.predict(X)
plt.scatter(X,y, alpha=0.6)
plt.plot(X, y_predict, color='r')
plt.show()

我们也可以得到模型的参数，可以通过参数把拟合方程描述处理，但如果用这个拟合方程预测，输入的新数据都要做特征转换和标准化，比较麻烦，直接用上面的poly_reg.predict(X)更方便。

多项式回归的内容很多，过拟合、模型交叉检验、学习曲线等内容下篇再介绍。