机器学习之回归算法（多项式扩展）

目标函数（损失函数）分类

多项式介绍

1、传统方法的多项式

多项式回归就是我们要使用的方法，其本质上还是基于线性回归的数学原理，只不过是对损失函数进行巧妙的改变使得原本我们用线性回归求得一条直线，变为求得一条曲线。
所以对与多项式回归而言，就是需要我们为样本数据增加一个特征，使之可以用线性回归的原理更好的拟合样本数据，但是求出的是对于原始样本而言的非线性的曲线。

import numpy as np
import matplotlib.pyplot as plt
# 构建样本数据，从-3到3的100个随机数
# 从一个均匀分布[low,high)中随机采样，注意定义域是左闭右开，即包含low，不包含high.
x = np.random.uniform(-3, 3, size=100)  # 横坐标(100,)
# 求y的方程
# 其中loc表示均值，scale表示方差，size表示输出的size
y = 0.5 * x**2 + x + 2 + np.random.normal(0, 1, size=100)  # 纵坐标100*1
# 将样本数据绘制出来(散点图）
plt.scatter(x, y)  # (100,)

# 转换为100行1列的矩阵，既样本数据只有一个特征
X = x.reshape(-1, 1)  # (100,1)
Y = y.reshape(-1, 1)  # (100,1)
# 导入线性回归类
from sklearn.linear_model import LinearRegression
lr = LinearRegression()
# 训练模型
lr.fit(X, Y)
# 预测出y值
y_predict = lr.predict(X)

# 绘制拟合直线
plt.plot(x, y_predict, color='r')
plt.show()

可以看到这条直线明显无法很好的拟合这些样本数据。下面我们来看看如何使用多项式回归来拟合样本数据：

import numpy as np
import matplotlib.pyplot as plt
# 构建样本数据，从-3到3的100个随机数
# 从一个均匀分布[low,high)中随机采样，注意定义域是左闭右开，即包含low，不包含high.
x = np.random.uniform(-3, 3, size=100)  # 横坐标(100,)
# 求y的方程
# 其中loc表示均值，scale表示方差，size表示输出的size
y = 0.5 * x**2 + x + 2 + np.random.normal(0, 1, size=100)  # 纵坐标100*1
# 将样本数据绘制出来(散点图）
plt.scatter(x, y)  # (100,)

# 转换为100行1列的矩阵，既样本数据只有一个特征
X = x.reshape(-1, 1)  # (100,1)
Y = y.reshape(-1, 1)  # (100,1)
# 导入线性回归类
from sklearn.linear_model import LinearRegression
lr = LinearRegression()

# 根据上文中讲的，我们需要先给样本数据增加一个特征，既X的平方
X2 = X**2    # (100,1)
# 通过np.hstack方法将原始特征和平方后的特征转换为新的特征矩阵，既100行2列的特征矩阵
# np.hstack():在水平方向上平铺
X_new = np.hstack((X, X2))   # (100,2)

# 然后再对新的样本特征用线性回归训练模型
lr.fit(X_new, Y)
y_predict_new = lr.predict(X_new)   # (100,1)
# 绘制出拟合曲线
# np.argsort(x)返回x中数据从小到大的索引x = ([3, 1, 2])返回[1, 2, 0]
plt.plot(np.sort(x), y_predict_new[np.argsort(x)], color='r')
plt.show()

lr2.coef_
# 结果
array([ 1.01051072, 0.45036129])
lr2.intercept_
# 结果
2.2249327738612221

可以看到，第一个特征的系数是1，第二个特征的系数，既对第一个特征做了平方的系数是0.45，截距是2.2。这和我们在之前构建的求y的方程是非常近似的，因为求y的方程里还加了正态分布的噪音，所以系数和截距是有点误差的，但是误差不大。那这就是多项式回归的原理和实现思路。我们在讲PCA的时候知道它的作用主要是给样本数据降维，而多项式回归是给样本数据升维，所以不同的算法，不同的样本数据是有降维，也有升维的，这点需要大家注意。

2、Scikit-Learn中的多项式

我们知道多项式回归其实就是对样本数据进行预处理，然后用线性回归的算法进行模型训练，所以多项式回归的核心在于多样本数据做预处理。那么在Scikit Learn中针对多项式回归的类PolynomialFeatures是在preprocessing包中：

import numpy as np
import matplotlib.pyplot as plt

# 构建样本数据
x = np.random.uniform(-3, 3, size=100)  # (100,)
X = x.reshape(-1, 1)  #(100, 1)
y = 0.5 * x ** 2 + x + 2 + np.random.normal(0, 1, size=100)  ## (100,)
Y = y.reshape(-1, 1)  #(100, 1)

# 导入多项式处理的类
from sklearn.preprocessing import PolynomialFeatures
# degree参数表示将样本数据处理为至少包含有几次幂的特征
'''
在实例化PolynomialFeatures时构造函数需要的degree参数表示将样本数据处理为至少包含有几次幂特征的数据。
比如degree=2，那么转换后的样本数据中的特征数量至少是2，其中一个是原始的特征，
另一个是将原始特征平方后的特征。

我们看到转换后的样本数据X2是一个100行3列的矩阵，说明通过PolynomialFeatures(degree=2)转换后，
我们期望的2个特征变成了3个，我们来看看前5行的样本数据：X2[:5, :],
array([[ 1.00000000e+00, -2.67091119e+00, 7.13376657e+00],
[ 1.00000000e+00, 9.26429770e-01, 8.58272118e-01],
[ 1.00000000e+00, 2.19969933e+00, 4.83867715e+00],
[ 1.00000000e+00, 3.13320069e-01, 9.81694657e-02],
[ 1.00000000e+00, -8.02234471e-02, 6.43580146e-03]])
可以看到转换后的第一列相当于是对特征做了0次方，值都是1，第二列是原始特征，
第三列是原始特征平方后的特征。
'''
poly = PolynomialFeatures(degree=2)
poly.fit(X)
X2 = poly.transform(X)   #(100, 3)

# 线性回归对转换后的样本数据进行处理
from sklearn.linear_model import LinearRegression
lr = LinearRegression()
lr.fit(X2, Y)
y_predict = lr.predict(X2)  #(100, 1)
plt.scatter(x, y)
plt.plot(np.sort(x), y_predict[np.argsort(x)], color='r')
plt.show()

以上就是在Scikit Learn中使用多项式回归的方式。下面我们再来仔细看一下PolynomialFeatures这个类。

# 构造简单的样本数据
X_simple = np.arange(1, 11).reshape(-1, 2)
X_simple
# 结果
array([[ 1, 2],
[ 3, 4],
[ 5, 6],
[ 7, 8],
[ 9, 10]])

构造一个简单的样本数据，元素从1到10，一共5行2列，既是一个有两个特征的样本数据，然后使用PolynomialFeatures对这个样本数据进行处理：

poly = PolynomialFeatures(degree=2)
poly.fit(X_simple)
X2_simple = poly.transform(X_simple)
X2_simple
# 结果
array([[ 1., 1., 2., 1., 2., 4.],
[ 1., 3., 4., 9., 12., 16.],
[ 1., 5., 6., 25., 30., 36.],
[ 1., 7., 8., 49., 56., 64.],
[ 1., 9., 10., 81., 90., 100.]])

我们看到当原始样本数据有2个特征，然后希望转换后的样本数据里增加对每个原始特征进行平方的新特征，通常情况下我们的预期可能是转换后有样本数据有4个特征，2个原始特征，2个对应做了平方处理的新特征。但是实际的情况是转换后的样本数据是一个5行6列的矩阵，既有6个特征的样本数据。

• 第一列全部为1。
• 第二列和第三列是原始特征。
• 第四列和第六列分别是第二列原始特征和第三列原始特征的平方。
• 第五列是第二列原始特征和第三列原始特征的乘积。

我们再来看一下如果将degree参数设置为3，转换后的样本数据会有几个特征：

poly = PolynomialFeatures(degree=3)
poly.fit(X_simple)
X3_simple = poly.transform(X_simple)
X3_simple.shape
# 结果
(5, 10)

3、pipeline

现在我们知道，使用多项回归，首先就是对样本数据进行预处理，但是当degree值比较大的时候，我们的样本数据的特征分布会非常分散，比如1的1次方和100的100次方，这使得样本数据的特征值分布非常不均衡，这会使训练出的模型误差比较大。在第三篇笔记中，我们讲过数据归一化，就是解决这个问题的工具，所以在多项式回归中，数据归一化也是我们经常需要处理的步骤，再然后才是线性回归的处理。

那么整个多项式回归至少就有三个步骤，分别为对样本数据进行多项式的预处理、对样本数据进行归一化、对样本数据进行线性回归处理。对不同的样本数据，每次都要做这三次步骤还是比较繁琐的，所以Scikit Learn给我们提供了一个Pipeline的工具，既可以将若干步骤打包成一个对象，相当于制作一个流程模板，这样面对不同的样本数据，我们只需要处理一次既可，Pipeline内部帮我们处理了打包在里面的其他步骤，极大提高了效率。

# 导入Pipeline和其他需要打包进Pipeline的类
from sklearn.pipeline import Pipeline
from sklearn.preprocessing import StandardScaler
from sklearn.preprocessing import PolynomialFeatures
# 封装Pipeline
poly_pipeline = Pipeline([
("poly", PolynomialFeatures(degree=2)),
("std_scalar", StandardScaler()),
("lr", LinearRegression())
])

构建Pipeline的过程就是封装步骤的过程，Pipeline构造函数的参数需要传入一个**数组，数组元素的类型为元组，元组的第一个元素为字符串，既标识这一步是做什么事的，第二个元素是处理对应步骤的类。**从上面构造的Pipeline可以看出，是将三个步骤进行了打包，第一个步骤的对样本数据进行多项式的预处理，第二个步骤是对样本数据进行归一化，第三个步骤是对样本数据进行线性回归处理。封装好Pipeline后，使用的方式和机器学习的算法保持一致：

poly_pipeline.fit(X, y)
y_predict_pipeline = poly_pipeline.predict(X)
plt.scatter(x, y)
plt.plot(np.sort(x), y_predict_pipeline[np.argsort(x)], color='r')
plt.show()

时间与电压的多项式扩展

原始的6个特征结果拟合不够好

多项式的扩展(阶数过大，有可能导致过拟合）

#导包
import numpy as np
import matplotlib as mpl
import matplotlib.pyplot as plt
import pandas as pd
import time

import sklearn
from sklearn.model_selection import train_test_split
from sklearn.linear_model import LinearRegression   
from sklearn.preprocessing import StandardScaler   # 标准化
from sklearn.metrics import mean_absolute_error,r2_score  #评估标准

## 设置字符集，防止中文乱码
mpl.rcParams['font.sans-serif']=[u'simHei']
mpl.rcParams['axes.unicode_minus']=False

# 加载数据
path='C:/Users/hubert/Desktop/Python零基础人工智能就业25G视频教程/5、机器学习/data/household_power_consumption_1000.txt' ## 1000行数据
df = pd.read_csv(path, sep=';', low_memory=False)   #当文件中的数据类型为单一类型时，low_memory能提高效率
# 异常数据处理(异常数据过滤)
new_df = df.replace('?', np.nan)    #用nan替换 有问题的数据,方便进行dropna操作
datas = new_df.dropna(axis=0,how = 'any') # any只要该行有数据为nan，就进行删除操作;all是全为nan时才删除行

"""
Pipeline:管道的意思，将多个操作合并在一起
Pipeline总可以给定多个操作，给定操作的名字即可，执行的时候从前到后

eg.如果执行fit_transform():
PolynomialFeatures中fit_transform对数据进行转换并构建模型，
然后对转换的数据调用LinearRegression的fit方法构建模型
"""
from sklearn.pipeline import Pipeline
from sklearn.preprocessing import PolynomialFeatures
# 给定进行多项式扩展操作：
# 1：多项式扩展
# 2：线性回归
models = [
    Pipeline([
        ('Poly',PolynomialFeatures()), 
        ('Linear',LinearRegression(fit_intercept=False))  # 是否计算偏置
    ])
]
model = models[0]

# 提取特征数据
## 创建一个时间字符串格式化字符串
def date_format(dt):
    #dt在这里是一行,是Series类型，这一行有两个数据
    #接下来要把这两个数据拼成字符串，再提取其中的‘年月日时分秒’
    import time
    t = time.strptime(' '.join(dt), '%d/%m/%Y %H:%M:%S')  #函数根据指定的格式把一个时间字符串解析为时间元组。
    return pd.Series([t.tm_year, t.tm_mon, t.tm_mday, t.tm_hour, t.tm_min, t.tm_sec])

## 时间和电压之间的关系(Linear)
# 获取x和y变量, 并将时间转换为数值型连续变量
X = datas.iloc[:,0:2]  # 选取datas中前三列
X = X.apply(date_format, axis=1)   #1代表按行，默认的是0代表列
Y = datas['Voltage'] 

# 划分数据集，标准化，模型训练，并使用模型预测
# 对数据集进行测试集和训练集划分
X_train,X_test,Y_train,Y_test = train_test_split(X, Y, test_size=0.2, random_state=0)  # 训练集80%

# 数据标准化，必须先fit_transform(X_train)，然后transform(X_test)，保证归一化标准一样
# fit,训练模型，没有返回值
#transform，对数据进行转化
#fit_transform，先训练，再转换
ss = StandardScaler()
X_train = ss.fit_transform(X_train) # 对训练集进行训练并标准化
X_test = ss.transform(X_test) #在训练集的fit标准下，用X_test进行的标准化模型训练，用训练好的模型，再对X_test进行转换

# 模型训练
t=np.arange(len(X_test))  # 0-199
N=5
d_pool=np.arange(1,N,1)  # 阶数1-4   [1, 2, 3, 4]
m=d_pool.size   #4

clrs=[]
# 颜色构建
for c in np.linspace(1671160,255,m):  # [1.67116000e+06, 1.11419167e+06, 5.57223333e+05, 2.55000000e+02]
    clrs.append('#%06x'%int(c))
line_width=3

# 窗口大小，颜色
plt.figure(figsize=(12,6),facecolor='w')

for i,d in enumerate(d_pool):  #[(0, 1), (1, 2), (2, 3), (3, 4)]
    plt.subplot(N-1,1,i+1)  # 411,412,413,414
    plt.plot(t,Y_test,'r-',label='真实值',ms=10,zorder=N)
    
    # 设置管道对象中的参数值，Poly是在管道对象中定义的操作名称，后面跟参数名称，设置多项式的阶乘；中间是两个下划线
    model.set_params(Poly__degree=d)  # 设置多项式的阶乘
    model.fit(X_train,Y_train)  # 模型训练
    # linear是管道中定义的操作名称
    # 获取线性回归算法模型对象
    lin = model.get_params('Linear')['Linear']
    output =u'%d阶，系数为：'%d
    # 判断lin对象中是否有对应的属性
    if hasattr(lin,'alpha_'):
        idx = output.find(u'系数')
        output = output[:idx]+(u'alpha=%.6f,'%lin.alpha_)+output[idx:]
    if hasattr(lin,'l1_ratio_'):
        idx = output.find(u'系数')
        output = output[:idx]+(u'l1_ratio=%.6f,'%lin.l1_ratio_)+output[idx:]
    print(output,lin.coef_.ravel())
    
    # 模型结果预测
    y_hat = model.predict(X_test)
    # 计算评估值
    s = model.score(X_test, Y_test)
    
    # 画图
    z = N - 1 if (d == 2) else 0
    label = u'%d阶, 准确率=%.3f' % (d,s)
    plt.plot(t, y_hat, color=clrs[i], lw=line_width, alpha=0.75, label=label, zorder=z)
    plt.legend(loc = 'upper left')
    plt.grid(True)
    plt.ylabel(u'%d阶结果' % d, fontsize=12)
# 预测值和实际值画图比较
plt.suptitle(u"线性回归预测时间和电压之间的多项式关系", fontsize=20)
plt.grid(b=True)
plt.show()