机器学习 | 多项式回归处理非线性问题

之前我们学习了一般线性回归，以及加入正则化的岭回归与Lasso，其中岭回归可以处理数据中的多重共线性，从而保证线性回归模型不受多重共线性数据影响。Lasso主要用于高维数据的特征选择，即降维处理。

线性回归中的多重共线性与岭回归

深度理解Lasso回归分析

在使用线性回归时，除了遇到以上问题（数据中存在多重共线性、数据维度过高），还会遇到数据并不总是线性的，若此时仍坚持用线性模型去拟合非线性数据，其结果往往是很糟糕的。庆幸地是除了更换更加复杂的非线性模型外，仍可以使用线性回归方法来拟合非线性数据，只不过在此之前，需要对输入数据进行处理。

多项式回归是一种通过增加自变量上的次数，而将数据映射到高维空间的方法，从而提高模型拟合复杂数据的效果。

---

本文主要运用对比分析方法解释线性数据、非线性数据、线性模型及非线性模型出发，并深入分析多项式回归处理非线性数据的原理、效果及代码实现。

在探究多项式回归之前，先对线性数据、非线性数据、线性模型及非线性模型做一个详细的介绍，以便更加深入地理解多项式回归在非线性数据集上使用线性模型的奥秘。

数据集的线性与非线性

变量之间的线性关系（linear relationship），表示两个变量之间的关系可以展示为一条直线，即可以使用方程 来拟合。绘制散点图通常是表示两个变量之间的线性关系的最简单的方式。如果散点图能够相对均匀地分布在一条直线的两端，则说明这两个变量之间的关系是线性的。而如三角函数（），高次函数（），指数函数（）等等图像不为直线的函数所对应的自变量和因变量之间是非线性关系（non-linear relationship）。

一般情况下，一组数据由多个变量和标签组成。变量分别与标签存在线性关系，则称他们是线性数据。而任意一个变量与标签之间的需要用三角函数、指数函数等来定义，则称其为"非线性数据"。

在回归中，绘制图像的是变量与标签的关系图，横坐标是特征，纵坐标是标签，标签是连续型的，则可以通过是否能够使用一条直线来拟合图像判断数据究竟属于线性还是非线性。

左图可以用 线性方程来进行拟合，称为线性数据；而右图拟合方程为 ， 为非线性方程，因此称之为非线性数据。

在分类中，绘制的是数据分布图，横纵坐标均是数据中的变量，颜色表示标签数据点。则使用 "是否线性可分" (linearly separable)来划分分类数据集。当分类数据的分布上可以使用一条直线来将两类数据分开时，则数据是线性可分的。反之，数据不是线性可分的。

这些数据与上面数据不同，都不能由一条直线来进行拟合，也没有均匀分布在某一条线的周围，但右图可以用一条直线将其分开，而左图却不能。

线性模型与非线性模型

回顾下线性回归模型，可以拟合出一组参数向量 ，从而建立一个模型，即线性回归模型。其建模过程即寻找参数向量的过程，用于拟合线性数据的线性模型。线性回归模型拟合的方程为线性方程，如下

boldsymbol{y=w_0 w_1x_1 w_2x_2 w_3x_3 dots w_nx_n}

而像决策树、支持向量机、各类树的集成模型，以及一切通过三角函数，指数函数等非线性方程来建立的模型。诸如此类的模型并不能像线性回归模型一样，使用形似 的线性方程来拟合数据。但他们被用于拟合非线性数据时，效果出奇的好（甚至过拟合）。

从图像上可以看出，线性回归模型无法拟合出这条带噪音的正弦曲线的真实面貌，只能够模拟出大概的趋势，使用线性回归模型来拟合非线性数据的效果并不好。这是因为线性模型假定自变量和因变量之间总是存在线性关系。这个假设是很弱的，它仅仅是近似。而非线性的决策树模型却通过建立复杂的模型将几乎每个点都拟合出来了，拟合得太细致，但相比之下，决策树在非线性数据上的拟合效果更好一些。

实际中机器学习中的模型很灵活。线性模型可以用来拟合非线性数据，而非线性模型也可以用来拟合线性数据，有的算法没有模型也可以处理各类数据，而有的模型可以既可以是线性，也可以是非线性模型。

---

非线性模型拟合线性数据

非线性模型拟合或处理线性数据的例子非常多，如随机森林，决策树等分类算法在处理线性可分的数据时，其效果并不逊于线性模型等表现。但与此同时，过拟合现象也更加容易出现在非线性模型拟合线性数据上，如利用随机森林拟合一条直线，因为简单的一条直线对于非线性模型来说过于简单，很容易就把训练集上的训练得很高，训练的很低。

非线性模型拟合非线性数据

诸如随机森林的非线性模型在线性模型上拟合效果不逊于线性模型，但更易过拟合。非线性模型真正的用武之地是非线性数据或非线性可分的数据集，其可以达到线性模型不能逾越的鸿沟。

---

线性模型拟合非线性数据

线性模型若用来拟合非线性数据或者对非线性可分的数据进行分类，那通常都会表现糟糕。

从上面的图中，线性模型们的决策边界都是一条条平行的直线，而非线性模型们的决策边界是交互的直线、曲线或环形等等。对于分类模型，线性模型的决策边界是平行的直线，非线性模型的决策边界是曲线或者交叉的直线。 回归模型中，若自变量上的最高次方为1，则模型是线性的 分类模型中，如果一个分类模型的决策边界上自变量的最高次方为1，则称这个模型是线性模型。

---

既是线性也是非线性的模型

有些模型既可以处理线性模型又可以处理非线性模型，如支持向量机。支持向量机的前身是感知机模型，朴实的感知机模型是线性模型，在线性可分数据上表现优秀，但在非线性可分的数据上基本属于无法使用状态。

而支持向量机通过选用不同的核函数可以在线性和非线性之间自由切换。如选用线性核函数"linear"，数据没有进行变换，就是线性模型。当选用非线性核函数时，数据进行升维变化，就是非线性模型。因此支持向量机在对不同的数据集选用合适的核函数，可以较灵活地高效地处理各种类型的数据。

非模型算法

没有模型的算法，如KNN，其原理不是建模拟合数据，也没有线性和非线性模型之分，但能够直接预测出标签或做出判断。

模型在线性和非线性数据集上的表现为我们选择模型提供了一个思路----当我们获取数据时，我们往往希望使用线性模型来对数据进行最初的拟合（线性回归用于回归，逻辑回归用于分类），如果线性模型表现良好，则说明数据本身很可能是线性的或者线性可分的，如果线性模型表现糟糕，那毫无疑问我们会选择决策树，随机森林等这些更加复杂的非线性模型的，不必浪费时间在线性模型上了。

不过这并不代表着我们就完全不能使用线性模型来处理非线性数据了。在现实中，线性模型有着不可替代的优势----计算速度异常快速，并存在着不得不使用线性回归的情况。

多项式回归

PolynomialFeatures

多项式回归通过增加额外的预测项对简单线性模型进行了拓展，即一个简单的线性回归可以通过从系数构造多项式特征来扩展。在标准线性回归的情况下，对于二维数据，你可能有一个这样的模型：

hat{y}(w, x) = w_0 w_1 x_1 w_2 x_2

如果我们想让数据拟合一个抛物面而不是一个平面，我们可以把这些特征合并成二阶多项式，使模型看起来像这样:

hat{y}(w, x) = w_0 w_1 x_1 w_2 x_2 w_3 x_1 x_2 w_4 x_1^2 w_5 x_2^2

更加一般地，多项式函数拟合数据时，多项式定义为

y(W, x) = w_0 w_1x w_2x^2 cdots w_Nx^N=sum_{j=0}^{N}w_jx^j

其中 为多项式的阶数，即最高次。 是多项式的系数，记做 ，

是关于 非线性函数，但是却是关于多项式系数 的线性函数。因此可以运用均方误差对多项式进行评估

E(W) = frac12sum_{n=1}^{N}(y(x_n,W)-epsilon_n)^2

狭义线性模型 自变量上不能有高此项，自变量与标签之间不能存在非线性关系。 广义线性模型 只要标签与模型拟合出的参数之间的关系是线性的，模型就是线性的。

线性模型中的升维工具----多项式变化。是一种通过增加自变量上的次数，而将数据映射到高维空间的方法，在sklearn中的类 PolynomialFeatures 设定一个自变量上的次数(大于1)，相应地获得数据投影在高次方的空间中的结果。

语法

代码语言：javascript复制

sklearn.preprocessing.PolynomialFeatures (degree=2, interaction_only=False, include_bias=True)

重要参数

degree : integer 多项式中的次数，默认为2 interaction_only : boolean, default = False 布尔值是否只产生交互项，默认为False。就只能将原有的特征进行组合出新的特征，而不能直接对原特征进行升次。 include_bias : boolean 布尔值，是否产出与截距项相乘的 ，默认True

例

代码语言：javascript复制

>>> from sklearn.preprocessing import PolynomialFeatures
>>> import numpy as np
>>> X = np.arange(6).reshape(3, 2)
>>> X
array([[0, 1],
       [2, 3],
       [4, 5]])
>>> poly = PolynomialFeatures(degree=2)
>>> poly.get_feature_names()
['1', 'x0', 'x1', 'x0^2', 'x0 x1', 'x1^2']
>>> poly.fit_transform(X)
array([[ 1.,  0.,  1.,  0.,  0.,  1.],
       [ 1.,  2.,  3.,  4.,  6.,  9.],
       [ 1.,  4.,  5., 16., 20., 25.]])

>>> # 转换成DataFrame
>>> poly2 = PolynomialFeatures(degree=3)
>>> X_2 = poly2.fit_transform(X)
>>> pd.DataFrame(data = X_2 ,columns=poly2.get_feature_names())

不难发现，数据的维度以一定的规律增加了，而这个规律可以通过poly.get_feature_names() 查看得到。

当我们进行多项式转换的时候，多项式会产出到最高次数为止的所有低高次项。这个最高次数可以通过 degree=n 来控制。

代码语言：javascript复制

>>> poly3 = PolynomialFeatures(degree=3, interaction_only=True)
>>> X_3= poly3.fit_transform(X)
>>> pd.DataFrame(X_3, columns=poly3.get_feature_names()) 

sklearn中存在着控制是否要生成平方和立方项的参数interaction_only ，因为存在只需求产生高次项的情况。其实同样产生一个高次项， 要比更好，因为 与之间的共线性会比与之间的共线性好一些。因为多项式回归模型，在经过多项式转化后仍需要使用线性模型进行拟合数据，若此时因转换数据带来额外的共线性，甚至更加严重的共线性将会严重影响模型拟合的结果。

---

多项式回归处理非线性问题

同样的一个问题，用线性回归模型无法拟合出这条带噪音的正弦曲线的真实面貌，只能够模拟出大概的趋势，而用复杂的决策树模型又拟合地太过细致，即过拟合。此时利用多项式将数据升维，并拟合数据，看看结果如何。

代码

代码语言：javascript复制

from sklearn.preprocessing import PolynomialFeatures as PF
from sklearn.linear_model import LinearRegression
import numpy as np
import matplotlib.pyplot as plt
rnd = np.random.RandomState(42) #设置随机数种子 
X = rnd.uniform(-3, 3, size=100)
y = np.sin(X)   rnd.normal(size=len(X)) / 3
#将X升维，准备好放入sklearn中 
X = X.reshape(-1,1)
#多项式拟合，设定高次项
d=5
#原始特征矩阵的拟合结果
LinearR = LinearRegression().fit(X, y)
#进行高此项转换
X_ = PF(degree=d).fit_transform(X)
LinearR_ = LinearRegression().fit(X_, y)
line = np.linspace(-3, 3, 1000, endpoint=False).reshape(-1, 1) 
line_ = PF(degree=d).fit_transform(line)
#放置画布
fig, ax1 = plt.subplots(1)
#将测试数据带入predict接口，获得模型的拟合效果并进行绘制
ax1.plot(line, LinearR.predict(line), linewidth=2, color='green'
         ,label="linear regression")
ax1.plot(line, LinearR_.predict(line_), linewidth=2, color='orange'
,label="Polynomial regression") #将原数据上的拟合绘制在图像上
ax1.plot(X[:, 0], y, 'o', c='k')
#其他图形选项
ax1.legend(loc="best")
ax1.set_ylabel("Regression output") 
ax1.set_xlabel("Input feature") 
ax1.set_title("Linear Regression ordinary vs poly") 
ax1.text(0.8,-1,f"LinearRegression score:n{LinearR.score(X,y)}",fontsize=15)
ax1.text(0.8,-1.3,f"LinearRegression score after poly :n{LinearR_.score(X_,y)}",fontsize=15)

plt.tight_layout()
plt.show()

输出结果

这里可以看出，多项式回归能够较好地拟合非线性数据，还不容易发生过拟合，可以说是保留了线性回归作为线性模型所带的"不容易过拟合"和"计算快速"的性质，同时又实现了优秀地拟合非线性数据。

这里我们使用了 degree=5，在实际应用中，我们并不能一次搞定degree的值。其实，不同的最高次取值，对模型拟合效果有重要的影响。

通过选择不同的 degree 对数据进行拟合：

代码

代码语言：javascript复制

# 生产数据函数
def uniform(size):
    x = np.linspace(0,1,size)
    return x.reshape(size,1)

def create_data(size):
    x = uniform(size)
    np.random.seed(42) #设置随机数种子 
    y = sin_fun(x) np.random.normal(scale=0.25, size=x.shape)
    return x,y

def sin_fun(x):
    return np.sin(2*np.pi*x)

代码语言：javascript复制

X_train,y_train = create_data(20)
X_test = uniform(200)
y_test = sin_fun(X_test)
fig = plt.figure(figsize=(12,8))
for i,degree in enumerate([0,1,3,6,9,12]):  
    plt.subplot(2,2,i 1)

    poly = PolynomialFeatures(degree)
    X_train_ploy = poly.fit_transform(X_train)
    X_test_ploy = poly.fit_transform(X_test)
    
    lr = LinearRegression()
    lr.fit(X_train_ploy,y_train)
    y_pred = lr.predict(X_test_ploy)

    plt.scatter(X_train,y_train,facecolor="none", edgecolor="g", s=25, label="training data")
    plt.plot(X_test,y_pred,c="orange",label="fitting")
    plt.plot(X_test,y_test,c="k",label="$sin(2pi x)$")
    plt.title("N={}".format(order))
    plt.legend(loc="best")
    plt.ylabel("Regression output") 
#     plt.xlabel("Input feature") 
    plt.legend()
plt.show()

输出结果

如图可见，黄色曲线是预测函数，黑色曲线是无噪声的原始曲线。 当 degree=0 与 degree=1 时多项式拟合数据效果极差，即欠拟合。其实就是线性模型拟合非线性数据。 当 degree=3、6、9 时，其实已经拟合很好了，当 degree=12 时，多项式函数对噪声数据有一定的敏感性，即过拟合。

前面有提到使用均方误差对拟合出的多项式进行评估，拟合数据的目的是最小化误差函数，因为误差函数是多项式系数 的二次函数，因此它关于系数 的导数是线性函数，所以误差函数的最小值有一个唯一解，我们记作 。

代码

代码语言：javascript复制

from sklearn.pipeline import Pipeline # 一种具有最终估计器的转换管道。文末有介绍
from sklearn.metrics import mean_squared_error
import math
import seaborn as sns
import matplotlib.patches as mpatches
# 数据还是上面例子的数据
mse1 = []
mse2 = []
plt.figure(figsize=(12,6))
degrees = [0,1,3,6,9,12]
for degree in degrees:
    pipe = Pipeline([("poly", PolynomialFeatures(degree=degree)),
                     ("lin_reg", LinearRegression())])
    pipe.fit(X_train,y_train)
    pred1 = pipe.predict(X_train)
    train_mse = mean_squared_error(y_train,pred1)
    
    pred2 = pipe.predict(X_test)
    test_mse = mean_squared_error(y_test,pred2)   
    
    mse1.append(train_mse)
    mse2.append(test_mse)
# 可视化
sns.pointplot(degrees,mse1,color='orange',label='Train')
sns.pointplot(degrees,mse2,color='k',label='Test')
orange_patch = mpatches.Patch(color='orange', label='Train MSE')
black_patch = mpatches.Patch(color='k', label='Test MSE')
plt.legend(handles = [red_patch,blue_patch])
plt.show()

输出结果

结合均方误差评价模型拟合效果的得分，可以更加精确地看出，在degree>3时MSE已经接近0了。而在degree=12 时，训练得分更加接近于0，若此时加上数据量过少，在模型训练时，模型更加专注于训练数据，导致模型过拟合。

由于数据量较少导致模型过拟合，可通过增加数据量，可同时增加模型复杂度（提高幂次degree的值）。但当我们增加幂次的值时，曲线开始高频震荡。这导致曲线的形状过于复杂，最终引起过拟合现象。此时可通过增加正则项来约束模型复杂度，即在最小化误差上增加惩罚项。此处可参见前面章节，岭回归于Lasso回归。

此外为了克服多项式回归的缺点（为了更好地拟合数据，增加多项式函数的复杂性，特征数量的增长很难控制），有学者提出使用样条回归法来克服多项式的众多缺点。这种方法没有将模型应用到整个数据集中，而是将数据集划分到多个区间，为每个区间中的数据单独拟合一个模型。具体可参见[https://zhuanlan.zhihu.com/p/36535032]

---

多项式回归的原理

由上面可知，用多项式回归模型拟合出这条带噪音的正弦曲线，其原理是泰勒级数展开。

泰勒公式泰勒公式是一个用函数在某点的信息描述其附近取值的公式。数学家们在柯西中值定理的基础上，推导出了泰勒中值定理（泰勒公式）。 若函数 在包含 的某个开区间 上具有 阶的导数，那么对于任一 , 有

f(x)=frac{f(x_0)}{0!} frac{f^prime(x_0)}{1!}(x-x_0) frac{f^{primeprime}(x_0)}{2!}(x-x_0)^2 cdots frac{f^{(n)}(x_0)}{n!}(x-x_0)^n R_n(x)

f(x)=frac{f(x_0)}{0!} frac{f^prime(x_0)}{1!}(x-x_0) frac{f^{primeprime}(x_0)}{2!}(x-x_0)^2 cdots frac{f^{(n)}(x_0)}{n!}(x-x_0)^n R_n(x)

其中，R_n(x)=frac{f^{n 1}(epsilon)}{(n 1)!}(x-x_0)^{n 1},epsilon为x_0与x之间的某个值。f(x)称为n阶泰勒公式

其中，R_n(x)=frac{f^{n 1}(epsilon)}{(n 1)!}(x-x_0)^{n 1},epsilon为x_0与x之间的某个值。f(x)称为n阶泰勒公式

一般情况下，泰勒公式在处展开。这样可以简化了泰勒公式得到在处的 阶泰勒公式，也称麦克劳林（Maclaurin）公式

f(x)=f(0) f^prime(0)(x) frac{f^{primeprime}(0)}{2!}(x)^2 cdots frac{f^{(n)}(0)}{n!}(x)^n R_n(x)

由于为与之间的某个值，可令 则其余项写为

o(x^n)=frac{f^{n 1}(theta x)}{(n 1)!}x^{n 1}(0<theta<1)

泰勒公式的几何意义是利用多项式函数来逼近原函数。

根据泰勒级数展开定理可知，任一函数可以用多项式表示，如

sin(x) = x-frac{x^3}{3!} frac{x^5}{5!}-,cdots,(-1)^{(n-1)}frac{x^{2n-1}}{(2n-1)!} R_{2n}

同理， 的泰勒展开为

e^x = 1 x frac{x^2}{2!} frac{x^3}{3!} ,cdots, frac{x^{n}}{n!} frac{e^{theta x}}{(n 1)!}x^{n 1}

当线性回归的模型太简单导致欠拟合时，我们就可以通过增加特征多项式来让线性回归模型更好的拟合数据。

---

多项式回归的可解释性

线性回归是一个具有高解释性的模型，它能够对每个特征拟合出参数以帮助我们理解每个特征对于标签的作用。当进行了多项式转换后，随着数据维度和多项式次数的上升，方程也变得异常复杂，但多项式回归的可解释性依然是存在的，我们可以使用接口get_feature_names来调用生成的新特征矩阵的各个特征上的名称。

案例

加利福利亚房屋数据集

代码语言：javascript复制

from sklearn.datasets import fetch_california_housing as fch
import pandas as pd
housevalue = fch()
X = pd.DataFrame(housevalue.data)
y = housevalue.target
housevalue.feature_names
X.columns = ["住户收入中位数","房屋使用年代中位数","平均房间数目" ,"平均卧室数目","街区人口","平均入住率","街区的纬度","街区的经度"]
poly = PolynomialFeatures(degree=2).fit(X,y)
poly.get_feature_names(X.columns)
X_poly = poly.transform(X) 
reg = LinearRegression().fit(X_poly,y)
coef = reg.coef_
[*zip(poly.get_feature_names(X.columns),reg.coef_)]
#放到dataframe中进行排序
coeff = pd.DataFrame([poly.get_feature_names(X.columns),reg.coef_.tolist()]).T 
coeff.columns = ["feature","coef"]
coeff.sort_values(by="coef").head(10)

输出结果

多项式回归虽然拟合了多项式曲线，但其本质仍然是线性回归，只不过我们将输入的特征做了些调整，增加了它们的多次项数据作为新特征。其实除了多项式回归，我们还可以使用这种方法拟合更多的曲线，我们只需要对原始特征作出不同的处理即可。

Pipeline

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

Pipeline原理示意图

pipeline最后一步如果有predict()方法我们才可以对 pipeline 使用 fit_predict()，同理，最后一步如果有transform()方法我们才可以对pipeline使用fit_transform()方法。 更多内容请参见https://scikit-learn.org/stable/modules/generated/sklearn.pipeline.Pipeline.html

代码语言：javascript复制

X,y = create_data(200) # 利用上面的生产数据函数
degree = 6
# 利用Pipeline将三个模型封装起来串联操作n
poly_reg = Pipeline([
                    ("poly", PolynomialFeatures(degree=degree)),
                    ("std_scaler", StandardScaler()),
                    ("lin_reg", LinearRegression())
                    ])
fig = plt.figure(figsize=(10,6)) 
poly = PolynomialFeatures(order)
poly_reg.fit(X, y)
y_pred = poly_reg.predict(X)
# 可视化结果
plt.scatter(X,y,facecolor="none", edgecolor="g", s=25, label="training data")
plt.plot(X,y_pred,c="orange",label="fitting")
# plt.plot(X,y,c="k",label="$sin(2pi x)$")
plt.title("degree={}".format(degree))
plt.legend(loc="best")
plt.ylabel("Regression output") 
plt.xlabel("Input feature") 
plt.legend()
plt.show()

输出结果

可以运用如下方法获取对应模型的接口。

代码语言：javascript复制

>>> poly_reg.named_steps['lin_reg'].coef_
array([[  0.        ,  -0.14987173,  16.4142791 , -65.75198762,
         87.4087241 , -45.6444121 ,   7.36587511]])

>>> poly_reg.named_steps['poly'].get_feature_names()
['1', 'x0', 'x0^2', 'x0^3', 'x0^4', 'x0^5', 'x0^6']

线性回归 决策树 编程算法 机器学习 神经网络

0 人点赞