Python机器学习基础教程学习笔记

笔记学习来源：Python机器学习基础教程

学习时间：2021年4月15日开始学习

第一章 引言

1.1 概念

监督学习

从输入输出对中进行学习的算法。

示例：识别信封上手写的邮政编码，基于医学影像判断肿瘤是否为良性，检测信用卡交易中的诈骗行为。

无监督学习

只有输入数据是已知的，没有为算法提供输出数据。

示例：确定一系列博客文章的主题，将客户分成具有相似偏好的群组，检测网站的异常访问模式。

1.2 环境搭建和工具介绍

scikit-learn

Scikit-learn（以前称为scikits.learn，也称为sklearn）是针对Python 编程语言的免费软件机器学习库。它具有各种分类，回归和聚类算法，包括支持向量机，随机森林，梯度提升，k均值和DBSCAN，并且旨在与Python数值科学图书馆NumPy和SciPy。

1
pip install numpy scipy matplotlib ipython scikit-learn pandas

必要的库和工具

  • NumPy

NumPy(Numerical Python) 是 Python 语言的一个扩展程序库，支持大量的维度数组与矩阵运算，此外也针对数组运算提供大量的数学函数库。NumPy 是一个运行速度非常快的数学库，主要用于数组计算，核心功能是ndarray类，即n维数组，数组中的所有元素必须是同一类型。

1
2
3
4
import numpy as np
# 创建一个 ndarray 只需调用 NumPy 的 array 函数即可
x = np.array([[1,2,3],[4,5,6]])
print("x:\n{}".format(x))

image-20210415145430297

  • SciPy

Scipy是一个用于数学、科学、工程领域的常用软件包，可以处理插值、积分、优化、图像处理、常微分方程数值解的求解、信号处理等问题。它用于有效计算Numpy矩阵，使Numpy和Scipy协同工作，高效解决问题。最重要的是scipy.sparse

1
2
3
4
5
from scipy import sparse
# 创建一个二维数组NumPy，对角线为1，其余都为0
# eye是scipy包中的一个创建特殊矩阵(单位矩阵E)的方法
eye = np.eye(4)
print("NumPy array:\n{}".format(eye))

image-20210415145815904

  • matplotlib

Matplotlib 是 Python 的绘图库。 它可与 NumPy 一起使用，提供了一种有效的 MatLab 开源替代方案。 它也可以和图形工具包一起使用，如 PyQt 和 wxPython。

1
2
3
4
5
6
7
8
import matplotlib.pyplot as plt
# 在-10和10之间生成一个数列，共100个数
x = np.linspace(-10,10,100)
# 用正弦函数创建另一个数组
y = np.sin(x)
# plot函数绘制一个数组关于另一个数组的折线图
plt.plot(x,y,marker="x")
plt.show()

image-20210415150214522

  • pandas

pandas，python+data+analysis的组合缩写，是python中基于numpy和matplotlib的第三方数据分析库，与后两者共同构成了python数据分析的基础工具包，享有数分三剑客之名。正因为pandas是在numpy基础上实现，其核心数据结构与numpy的ndarray十分相似，但pandas与numpy的关系不是替代，而是互为补充。

默认导入环境

1
2
3
4
5
import numpy as np
from scipy import sparse
import matplotlib.pyplot as plt
import pandas as pd
import mglearn

1.3 第一个应用：鸢尾花分类

  • 目标：构建一个机器学习模型，可以从这些已知品种的鸢尾花测量数据中进行学习，从而能够预测新鸢尾花的品种，这一模型为监督学习。
  • 核心术语：单个数据点（一朵鸢尾花）的预期输出是这朵花的品种，对于一个数据点来说，他的品种叫做标签。

1.3.1 初识数据

iris数据集

1
2
3
# 经典数据集，可以调用load_iris()函数来加载数据
from sklearn.datasets import load_iris
iris_dataset = load_iris()

iris_dataset对象是一个Bunch对象，里面包括键和值。

1
2
3
print("Keys of iris_dataset:\n{}".format(iris_dataset.keys()))
# 注：format()函数
# 一种格式化字符串的函数 str.format()，它增强了字符串格式化的功能。基本语法是通过{}和:来代替以前的%

image-20210415154827716

DESCR键

DESCR键对应的值是数据集的简要说明。

1
print(iris_dataset['DESCR'][:193] + "\n...")

image-20210415155400917

target_names键

target_names键对应的值是一个字符串数组，里面包含要预测的花的品种。

1
print("Target names: {}".format(iris_dataset['target_names']))

image-20210415155544844

feature_names键

feature_names键对应的值是一个字符串列表，对每一个特征进行了说明。

1
print("Feature names: \n{}".format(iris_dataset['feature_names']))

image-20210415155723911

data键

data键对应的值是二维NumPy数组，每一行对应一朵花，列代表每朵花的四个测量数据。

1
print("Shape of data: {}".format(iris_dataset['data'].shape))

image-20210415155959901

可以看出，数组中包含150朵不同的花的测量数据。data数组的形状（shape）是样本数乘以特征数，这是scikit-learn的约定。

1
print("First five rows of data:\n{}".format(iris_dataset['data'][:5]))

image-20210415160207996

以上为data数组的前五行数据。

target键

target键对应的值是一个一维NumPy数组，包含的是测量过的每朵花的品种（标签）。品种被转换为0到2的整数。

1
print("Target:\n{}".format(iris_dataset['target']))

image-20210415160427681

1.3.2 训练数据和测试数据

我们将数据集分成训练集（training set）和测试集（test set），前者用于构建模型，后者用于评估模型对于前所未见的新数据的泛化能力。

  • 训练数据：用于构建机器学习模型，又叫做训练集。
  • 测试数据：用于评估机器学习模型，又叫做测试集，留出集。

train_test_split函数

train_test_split函数可以将数据集打乱并进行拆分，这个函数将75%的行数据及对应标签作为训练集，剩下的25%作为测试集。通常，数据用大写的X表示，而标签用小写的y表示，X为输入的二维数组，y则为输出的一维数组。

1
2
3
4
5
6
7
8
9
from sklearn.datasets import load_iris
from sklearn.model_selection import train_test_split
iris_dataset = load_iris()
# X_train包含了75%的行数据，X_test包含了25%的行数据
X_train,X_test,y_train,y_test = train_test_split(iris_dataset['data'],iris_dataset['target'],random_state=0)
print("X_train shape: {}".format(X_train.shape))
print("y_train shape: {}".format(y_train.shape))
print("X_test shape: {}".format(X_test.shape))
print("y_test shape: {}".format(y_test.shape))

image-20210415161708657

1.3.3 观察数据

暂略

1.3.4 构建第一个模型：k近邻算法

k近邻算法含义

scikit-learn中有许多可用的分类算法。要对一个新的数据点做出预测，k近邻算法会在训练集中寻找与这个新数据点距离最近的数据点，然后将找到的数据点的标签赋给这个新数据点。

k近邻算法中的k的含义就是考虑训练集中与新数据点最近的任意k个邻居。

k近邻算法实在neighbors模块的KNeighborsClassifier类中实现的，最重要的参数就是这个k，即邻居的数量。

代码实现

1
2
3
4
5
6
from sklearn.neighbors import KNeighborsClassifier
# knn对象对算法进行了封装，既包括了用训练数据构建模型的算法，也包括对新数据点进行预测的算法
knn = KNeighborsClassifier(n_neighbors=1)
# 想要基于训练集来构建模型，需要调用fit方法，输入参数为X_train和y_train，二者都是NumPy数组
# 前者为训练数据，后者为相应的训练标签
knn.fit(X_train,y_train)

1.3.5 做出预测

新的数据点

1
2
# 注意X_new为二维数组，作为其输入
X_new = np.array([[5,2.9,1,0.2]])

做出预测——predict

1
2
3
4
# 调用knn对象的predict方法进行预测
prediction = knn.predict(X_new)
print("Prediction: {}".format(prediction))
print("Predicted target name :{}".format(iris_dataset['target_names'][prediction]))

image-20210415170753481

1.3.6 评估模型

可以对测试数据中的每朵鸢尾花进行预测，并将预测结果与标签进行对比。

代码实现

1
2
3
4
5
6
7
iris_dataset = load_iris()
X_train,X_test,y_train,y_test = train_test_split(iris_dataset['data'],iris_dataset['target'],random_state=0)
knn = KNeighborsClassifier(n_neighbors=1)
# 训练
knn.fit(X_train,y_train)
# 调用knn对象的score方法计算测试集的精度
print("Test set score: {:.2f}".format(knn.score(X_test,y_test)))

image-20210415172433706

可知该模型的准确率为97%。

1.3.7 小结

下面汇总了整个训练和评估过程所必需的代码：

1
2
3
4
5
6
7
8
# 划分原始数据集为训练集和测试集
X_train,X_test,y_train,y_test = train_test_split(iris_dataset['data'],iris_dataset['target'],random_state=0)
# 构造k近邻算法模型
knn = KNeighborsClassifier(n_neighbors=1)
# 训练该模型
knn.fit(X_train,y_train)
# 调用knn对象的score方法计算测试集的精度
print("Test set score: {:.2f}".format(knn.score(X_test,y_test)))

这个代码片段包含了应用scikit-learn中任何机器学习算法的核心代码。fit，predict和score方法是scikit-learn监督学习模型中最常用的接口。

第二章 监督学习

2.1 分类与回归

监督机器学习问题主要有两种：分类（classification）和回归（regression）

分类

分类问题的目标是预测类别标签，可分为二分类和多分类问题。

回归

回归问题的目标是预测一个连续值。

区分分类和回归，就是问一个问题，输出是否具有某种连续性。

2.2 泛化，过拟合和欠拟合

泛化 generalization

如果一个模型能够对没见过的数据做出准确预测，我们就说它能够从训练集泛化到测试集。

过拟合 overfitting

如果在拟合模型时过分关注训练集的细节，得到一个在训练集上表现很好，但不能泛化到新数据上的模型，那么就存在过拟合。这种情况通常是构建了一个对现有信息（训练集）过于复杂的模型。越靠近训练集，越过拟合。或者，训练集和测试集的效率差异太大，也是过拟合。

欠拟合 underfitting

和过拟合相反，构建了过于简单的模型，没有充分考虑到训练集中的各种有效信息。或者，训练集和测试集的效率相似且很差，就是欠拟合。

总结

我们的模型越复杂，在训练数据上的预测结果就越好，但是，如果我们的模型过于复杂，我们开始过多关注训练集中每个单独的数据点，模型就不能很好的泛化到新数据上去。

经验：收集更多的数据，适当构建更复杂的模型，对监督学习任务往往特别有用，在现实世界中，你往往能够决定收集多少数据，这可能比模型调参更为有效，永远不要低估更多数据的力量。

2.3 监督学习算法

许多算法都有分类和回归两种模式。

2.3.1 一些样本数据集

二分类数据集——forge

1
2
3
4
5
6
7
8
9
10
11
import mglearn
import matplotlib.pyplot as plt
# 生成数据集
X, y = mglearn.datasets.make_forge()
# 数据集绘图
mglearn.discrete_scatter(X[:,0], X[:,1], y)
plt.legend(["Class 0", "Class 1"], loc=4)
plt.xlabel("First feature")
plt.ylabel("Second feature")
plt.show()
print("X.shape: {}".format(X.shape))

image-20210416153203952

图像以第一个特征为x轴，第二个特征为y轴，每个数据点对应图像中的一点。

image-20210416153300067

还可知，X数据集总共26个数据，2个特征。

回归用数据集——wave

wave数据集只有一个输入特征和一个连续的目标变量（或==响应==），后者是模型想要预测的对象（一个连续值）。

1
2
3
4
5
6
X, y = mglearn.datasets.make_wave(n_samples=40)
plt.plot(X, y, 'o')
plt.ylim(-3,3)
plt.xlabel("Feature")
plt.ylabel("Target")
plt.show()

image-20210416153737341

图中x轴表示特征，y轴表示回归目标（对应的响应）。

现实中的二分类数据集——cancer

1
2
3
from sklearn.datasets import load_breast_cancer
cancer = load_breast_cancer()
print("keys of cancer: \n{}".format(cancer.keys()))

image-20210416154018036

该数据集包含569个数据，每个数据有30个特征。212个被标记为恶性（malignant），357个被标记为良性（benign）。

现实中的回归数据集——boston

与这个数据集相关的任务是，利用犯罪率，是否临近查尔斯河，公路可达性等信息，预测20世纪70年代波士顿地区房屋价格的中位数，这个数据集包含506个数据点和13个特征。

1
2
3
from sklearn.datasets import load_boston
boston = load_boston()
print("Data shape: {}".format(boston.data.shape))

image-20210416154442447

2.3.2 k近邻算法

k近邻分类算法

分为单近邻和多近邻。可用第一章的知识分析两者的效率区别。

分析KNeighborsClassifier

查看决策边界，即算法对类别0和类别1的分界线。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
import numpy as np
import matplotlib.pyplot as plt
import pandas as pd
import mglearn
from sklearn.neighbors import KNeighborsClassifier
X, y = mglearn.datasets.make_forge()
fig, axes = plt.subplots(1,3, figsize=(10,3))
for n_neighbors, ax in zip([1,3,9], axes):
    # fit 方法返回对象本身，所以我们可以将实例化和拟合放在一行代码中
    clf = KNeighborsClassifier(n_neighbors=n_neighbors).fit(X,y)
    mglearn.plots.plot_2d_separator(clf, X, fill=True, eps=0.5, ax=ax, alpha=.4)
    mglearn.discrete_scatter(X[:,0], X[:,1], y, ax=ax)
    ax.set_title("{} neighbor(s)".format(n_neighbors))
    ax.set_xlabel("feature 0")
    ax.set_ylabel("feature 1")
axes[0].legend(loc=3)
plt.show()

image-20210416160626726

从图中可以看出，使用单一邻居绘制的决策边界紧跟着训练数据，模型越复杂，随着邻居的个数越来越多，决策边界也越来越平滑，更平滑的边界对应着更简单的模型。

  • 应用到cancer数据集的结果：
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
from sklearn.datasets import load_breast_cancer
from sklearn.model_selection import train_test_split

cancer = load_breast_cancer()
X_train, X_test, y_train, y_test = train_test_split(cancer.data, cancer.target, stratify=cancer.target, random_state=66)

training_accuracy=[]
test_accuracy=[]

# neighbors取值从1到10
neighbors_settings=range(1,11)

for n_neighbors in neighbors_settings:
    # 构建模型
    clf = KNeighborsClassifier(n_neighbors=n_neighbors)
    clf.fit(X_train, y_train)
    # 记录训练集精度
    training_accuracy.append(clf.score(X_train, y_train))
    # 记录泛化精度
    test_accuracy.append(clf.score(X_test, y_test))

plt.plot(neighbors_settings, training_accuracy, label="training accuracy")
plt.plot(neighbors_settings, test_accuracy, label="test accuracy")
plt.ylabel("Accuracy")
plt.xlabel("n_neighbors")
plt.legend()
plt.show()

image-20210416161223780

从图中可知，最佳点在邻居为6的附近，6之前的越靠近训练集，则过拟合，导致泛化能力下降，6以后的越远离训练集，则欠拟合，泛化能力也不好。

k近邻回归算法

单一近邻时，利用单一邻居的预测结果就是最近邻居的目标值。

image-20210416161536580

使用多个近邻时，预测结果为这些邻居的平均值。

image-20210416161607389

利用k近邻回归算法在scikit-learn的KNeighborsRegressor类中实现。

1
2
3
4
5
6
7
8
9
10
from sklearn.neighbors import KNeighborsRegressor
X, y=mglearn.datasets.make_wave(n_samples=40)
# 划分训练测试集
X_train, X_test, y_train, y_test=train_test_split(X,y, random_state=0)
# 实例化模型，邻居设定为3
reg=KNeighborsRegressor(n_neighbors=3)
# 拟合模型
reg.fit(X_train, y_train)
print("Test set predictions:\n{}".format(reg.predict(X_test)))
print("Test set R^2: {:.2f}".format(reg.score(X_test, y_test)))

image-20210416161917275

可知准确率为0.83。

分析KNeighborsRegressor

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
from sklearn.neighbors import KNeighborsRegressor
X, y=mglearn.datasets.make_wave(n_samples=40)
# 划分训练测试集
X_train, X_test, y_train, y_test=train_test_split(X,y, random_state=0)
fig, axes = plt.subplots(1,3, figsize=(15,4))
#创建1000个数据点，在-3和3之间均匀分布
line=np.linspace(-3,3,1000).reshape(-1,1)
for n_neighbors,ax in zip([1,3,9], axes):
    # 利用1，3，9个邻居分别进行预测
    reg=KNeighborsRegressor(n_neighbors=n_neighbors)
    reg.fit(X_train, y_train)
    ax.plot(line, reg.predict(line))
    ax.plot(X_train, y_train, '^', c=mglearn.cm2(0), markersize=8)
    ax.plot(X_test, y_test, 'v', c=mglearn.cm2(1), markersize=8)
    ax.set_title(
        "{} neighbor(s)\n train score: {:.2f} test score: {:.2f}".format(
            n_neighbors, 
            reg.score(X_train, y_train),
            reg.score(X_test, y_test)))
    ax.set_xlabel("Feature")
    ax.set_ylabel("Target")
axes[0].legend(["Model predictions", "Training data/target",
                "Test data/target"], loc="best")
plt.show()

image-20210416162210898

从图中可知，仅使用单一邻居，训练集中的每个点对预测结果都有显著影响，预测结果的图像经过所有数据点，这导致测试结果非常不稳定，考虑更多邻居后，预测结果变得更加平滑，但对训练数据的拟合也不太好。

k近邻算法优缺点分析

  • 优点：容易理解，模型参数调节简单
  • 缺点：训练集很大时，预测速度会变慢，这一算法对有很多特征的数据集往往效果不好。

综上，k近邻算法由于预测速度慢且不能处理多特征的数据集，所以实践中往往不会用到，但可作为基准测试方法使用。

2.3.3 线性模型

线性模型利用输入特征的线性函数进行预测。

用于回归的线性模型

对于回归问题，线性模型预测的一般公示如下：

1
y = w[0]*x[0] + w[1]*x[1] + ... + w[p]*x[p] +b

x表示单个数据点的特征，w和b是学习模型的参数，y是模型的预测结果。

下面线性模型在wave数据集上学习的参数w和b：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
import numpy as np
import matplotlib.pyplot as plt

from sklearn.linear_model import LinearRegression
from sklearn.model_selection import train_test_split
from mglearn.datasets import make_wave
from mglearn.plot_helpers import cm2

def plot_linear_regression_wave():
  	# 获得wave数据集
    X, y = make_wave(n_samples=60)
    # 拆分数据集为训练集和测试集
    X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=42)

    line = np.linspace(-3, 3, 100).reshape(-1, 1)
		
    # 利用线性回归模型进行拟合
    lr = LinearRegression().fit(X_train, y_train)
    print("w[0]: %f  b: %f" % (lr.coef_[0], lr.intercept_))

    plt.figure(figsize=(8, 8))
    plt.plot(line, lr.predict(line))
    plt.plot(X, y, 'o', c=cm2(0))
    ax = plt.gca()
    ax.spines['left'].set_position('center')
    ax.spines['right'].set_color('none')
    ax.spines['bottom'].set_position('center')
    ax.spines['top'].set_color('none')
    ax.set_ylim(-3, 3)
    ax.legend(["model", "training data"], loc="best")
    ax.grid(True)
    ax.set_aspect('equal')

if __name__ == "__main__":
    plot_linear_regression_wave()
    plt.show()

image-20210417143438717

image-20210417143514365

有许多线性回归模型，区别在于如何从训练集中学习到参数w和b，以及如何控制模型的复杂度。

标准线性回归（普通最小二乘法）

标准线性回归寻找参数w和b，使得对训练集的预测值与真实的回归目标值y之间的均方误差最小，均方误差是指预测值和真实值之差的平方和除以样本数。线性回归没有额外的控制参数，因此无法控制模型的复杂度。

  • w：斜率参数，或者权重，系数，保存在coef_属性中，为NumPy数组。sklearn将从训练数据中得出的值保存在以下划线结尾的属性中，与用户设置的参数区分开。

  • b：偏移，或者截距，保存在intercept_属性中，为浮点数。

1
2
3
4
5
6
7
8
from sklearn.linear_model import LinearRegression
import mglearn
from sklearn.model_selection import train_test_split
X, y = mglearn.datasets.make_wave(n_samples=60)
X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=42)
lr = LinearRegression().fit(X_train, y_train)
print("lr.coef_: {}".format(lr.coef_))
print("lr.intercept_: {}".format(lr.intercept_))

image-20210417144639249

因为wave数据集的数据只有一个特征，则coef_的元素个数为1。

在boston数据集上的准确率测试：

1
2
3
4
5
X, y = mglearn.datasets.load_extended_boston()
X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=0)
lr = LinearRegression().fit(X_train, y_train)
print("Traning set score: {:.2f}".format(lr.score(X_train, y_train)))
print("Test set score: {:.2f}".format(lr.score(X_test, y_test)))

image-20210417145641931

可以看出，训练集和测试集的效率差异过大，是过拟合的标志。

岭回归——ridge regression

岭回归特点：

  1. 对系数w的拟合附加约束
  2. L2正则化：希望w尽可能小，都接近于0，这样意味着每个特征对输出的影响应尽可能小，模型越简单，避免过拟合。用数学术语来讲，Ridge惩罚了系数w的L2范数或者w的欧式长度。惩罚系数（penalty）为alpha。
  3. 岭回归在linear_model.Ridge中实现。

在boston数据集上的准确率测试：

1
2
3
4
5
6
from sklearn.linear_model import Ridge
X, y = mglearn.datasets.load_extended_boston()
X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=0)
ridge = Ridge().fit(X_train, y_train)
print("Traning set score: {:.2f}".format(ridge.score(X_train, y_train)))
print("Test set score: {:.2f}".format(ridge.score(X_test, y_test)))

image-20210417150004961

可以看出，训练集上的准确率不及线性回归，但是测试集准确率高，表示泛化能力强。

岭回归模型通过用户设置alpha惩罚参数来在模型的简单性和训练集性能之间做出权衡。

==增大惩罚系数alpha，会使得系数coef_(或w)越趋近于0，正则化约束越强，模型越简单，从而降低了训练集性能，但是泛化能力可能会增强==。

1
ridge = Ridge(alpha = 10).fit(X_train, y_train)

lasso回归

lasso回归模型特点：

  1. 对系数w的拟合附加约束
  2. L1正则化：使用lasso时，某些系数刚好为0，这说明某些特征被模型完全忽略。lasso惩罚了系数向量的L1范数，即系数的绝对值之和。惩罚系数记也为alpha。
  3. 岭回归在linear_model.Lasso中实现。

lasso在boston数据集上的表现：

1
2
3
4
5
6
7
8
from sklearn.linear_model import Lasso
import numpy as np
X, y = mglearn.datasets.load_extended_boston()
X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=0)
lasso = Lasso().fit(X_train, y_train)
print("Traning set score: {:.2f}".format(lasso.score(X_train, y_train)))
print("Test set score: {:.2f}".format(lasso.score(X_test, y_test)))
print("Number of features userd: {}".format(np.sum(lasso.coef_ != 0)))

image-20210417151338608

可以看出，lasso在训练集和测试集上的表现都很差，说明存在欠拟合，并且只使用了4个特征，而单个数据点具有105个特征。

lasso默认alpha的值为1.0，为了提高两者表现，需要减小alpha的值，使其拟合一个更复杂的模型。

1
2
# 减小alpha的值，同时增加max_iter的值，增加运行迭代的次数
lasso = Lasso(alpha=0.01,max_iter=100000).fit(X_train, y_train)

image-20210417151653278

如果把alpha设置的太小，那么就会消除正则化的约束，使其变得和标准线性回归的结果一样。

在实践中，一般选择岭回归模型，如果考虑到模型的可解释性，即存在某几个特征的影响相比于其它的更显著，那么选择lasso回归模型。

用于二分类的线性模型

线性模型也可用于分类问题。有二分类公式：

1
y = w[0]*x[0] + w[1]*x[1] + ... + w[p]*x[p] + b > 0

上式为预测设置了阈值（0），当函数值小于0时，则贴上标签-1，否则为+1。

最常见的两种线性分类算法为Logistic回归和线性支持向量机（linear support vector machine，简称线性SVM），注意前者虽然叫做“回归”，其实是用于分类算法。

两者特点：

  1. 使用L2正则化。权衡系数记为C，与惩罚系数alpha的效果==相反==。
  2. 前者在linear_model.LogisticRegression中实现，后者在svm.LinearSVC中实现。

两者在forge数据集上的决策边界可视化结果：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
from sklearn.linear_model import LogisticRegression
from sklearn.svm import LinearSVC
import matplotlib.pyplot as plt
X, y = mglearn.datasets.make_forge()
fig, axes = plt.subplots(1, 2, figsize=(10, 3))
for model, ax in zip([LinearSVC(),LogisticRegression()], axes):
  	# 拟合模型
    clf = model.fit(X,y)
    mglearn.plots.plot_2d_separator(clf, X, fill=False, eps=0.5, ax=ax, alpha=.7)
    mglearn.discrete_scatter(X[:,0], X[:,1], y, ax=ax)
    ax.set_title("{}".format(clf.__class__.__name__))
    ax.set_xlabel("Feature 0")
    ax.set_ylabel("Feature 1")
axes[0].legend()
plt.show()

image-20210417153247677

注意：两个模型中都有两个点的分类是错误的。

参数C对线性SVM的决策边界影响结果：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
def plot_linear_svc_regularization():
    X, y = make_blobs(centers=2, random_state=4, n_samples=30)
    fig, axes = plt.subplots(1, 3, figsize=(12, 4))

    # a carefully hand-designed dataset lol
    y[7] = 0
    y[27] = 0
    x_min, x_max = X[:, 0].min() - .5, X[:, 0].max() + .5
    y_min, y_max = X[:, 1].min() - .5, X[:, 1].max() + .5

    for ax, C in zip(axes, [1e-2, 10, 1e3]):
        discrete_scatter(X[:, 0], X[:, 1], y, ax=ax)
				# 对SVM设置权衡参数C
        svm = LinearSVC(C=C, tol=0.00001, dual=False).fit(X, y)
        w = svm.coef_[0]
        a = -w[0] / w[1]
        xx = np.linspace(6, 13)
        yy = a * xx - (svm.intercept_[0]) / w[1]
        ax.plot(xx, yy, c='k')
        ax.set_xlim(x_min, x_max)
        ax.set_ylim(y_min, y_max)
        ax.set_xticks(())
        ax.set_yticks(())
        ax.set_title("C = %f" % C)
    axes[0].legend(loc="best")

if __name__ == "__main__":
    plot_linear_svc_regularization()
    plt.show()

image-20210417153551456

Logistic回归在cancer数据集上的测试：

1
2
3
4
5
6
7
8
from sklearn.datasets import load_breast_cancer
from sklearn.linear_model import LogisticRegression
cancer = load_breast_cancer()
X_train, X_test, y_train, y_test = train_test_split(cancer.data, cancer.target, stratify=cancer.target, random_state=42)
# 设置C为100，默认为1，则正则化约束变弱，模型更复杂
logreg = LogisticRegression(C=100).fit(X_train, y_train)
print("Traning set score: {:.3f}".format(logreg.score(X_train, y_train)))
print("Test set score: {:.3f}".format(logreg.score(X_test, y_test)))

image-20210417154427555

用于多分类的线性模型

许多线性分类模型只适用于二分类问题，将其推广到多分类算法的一种常见方法是“一对其余”。对每个类别都学习一个二分类模型，将这个类别与其所有其他类别尽量区分开，这样就生成了与类别个数相同的二分类模型。在测试点上运行所有二分类器来进行预测，在对应类别上分数最高的分类器胜出，将这个类别标签返回作为预测结果。

以下是一个三分类的数据集：

1
2
3
4
5
6
7
8
from sklearn.datasets import make_blobs
import matplotlib.pyplot as plt
X, y = make_blobs(random_state=42)
mglearn.discrete_scatter(X[:,0],X[:,1],y)
plt.xlabel("Feature 0")
plt.ylabel("Feature 1")
plt.legend(["Class 0", "Class 1", "Class 2"])
plt.show()

image-20210417163703497

在这个数据集上训练一个线性SVC分类器

1
2
3
4
5
6
from sklearn.datasets import make_blobs
from sklearn.svm import LinearSVC
X, y = make_blobs(random_state=42)
linear_svm = LinearSVC().fit(X, y)
print("Coeffient shape: ",linear_svm.coef_.shape)
print("Intercept shape: ",linear_svm.intercept_.shape)

image-20210417164023296

coef_的形状为（3,2），说明每行包含三个之一的系数向量，每列包含2个特征的系数值。可视化结果为：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
from sklearn.datasets import make_blobs
from sklearn.svm import LinearSVC
import matplotlib.pyplot as plt
import numpy as np
X, y = make_blobs(random_state=42)
linear_svm = LinearSVC().fit(X, y)
mglearn.discrete_scatter(X[:,0],X[:,1],y)
line = np.linspace(-15,15)
for coef, intercept, color in zip(linear_svm.coef_, linear_svm.intercept_, ['b','r','g']):
    plt.plot(line, -(line * coef[0] + intercept) / coef[1], c=color)
plt.ylim(-10,15)
plt.xlim(-10,8)
plt.xlabel("Feature 0")
plt.ylabel("Feature 1")
plt.legend(["Class 0", "Class 1", "Class 2","line class 0", "line class 1","line class 2"],loc=(1.01,0.3))
plt.show()

image-20210417164716398

image-20210417165104817

总结

线性模型的主要参数时正则化参数，或称为惩罚参数，权衡参数，在回归模型中为alpha，分类模型中为C。

  • 优点：训练速度快，预测速度也快，可以推广到非常大的数据集，或者特征数量很多的数据集。

  • 缺点：在低维数据集中（特征数较少），表现不是很好。

补充：方法链

sklearn中的所有模型的fit方法返回的都是self，因此可以编写下面的代码：

1
2
# 用一行代码初始化模型并拟合
logreg = LogisticRegression().fit(X_train, y_train)

这种方式被称为方法链。

1
2
3
4
# 初始化模型
logreg = LogisticRegression()
# 另一个常见的用法是在一行代码中同时实现拟合和预测
y_pred = logreg.fit(X_train, y_train).predict(X_test)

甚至可以在一行代码中实现模型初始化，拟合和预测，但不推荐，因为代码可读性差。此外，拟合后的模型也没有保存在任何变量中，形成了“浪费”。

1
y_pred = LogisticRegression().fit(X_train, y_train).predict(X_test)

2.3.4 朴素贝叶斯分类器

特点

与线性模型分类器相似，但训练速度更快，因此泛化能力相比较弱。

三类分类器

  • GaussianNB：可应用于任意连续数据，高维数据，文本数据分类。
  • BernoulliNB：它假定输入数据为二分类数据。
  • MultinomiaNB：它假定输入数据为计数数据，例如统计单词出现次数，可用于文本数据分类。

2.3.5 决策树

决策树是广泛应用于分类和回归任务的模型。本质上，它从一层层的if/else问题中进行学习，并得出结论。

决策树的目标就是通过提出尽可能少的if/else问题来得到正确答案。

image-20210425144035732

上图，用机器语言来讲，就是为了区分四种动物，利用三个特征来构建一个模型。

构造决策树

数据通常并不是上例中的那样具有二元性（是或否）的形式，而是表示为连续特征。用于连续数据的测试形式为：==特征i的值是否大于a？==

对数据进行反复递归划分，直到划分后的每个区域（决策树的每个叶节点）只包含单一目标值（单一类别或单一回归值）。如果树中某个叶节点所包含数据点的目标值都相同，那么这个叶节点就是纯的。

决策树也可用于回归任务，预测的方法是：基于每个节点的测试对树进行遍历，最终找到新数据点所属的叶节点。这一数据点的输出就是此叶节点中所有训练点的平均目标值。

控制决策树的复杂度——防止过拟合

通常来说，构造决策树直到所有叶节点都是纯的叶节点，会导致模型非常复杂（对训练集数据的预测率将会是100%），泛化能力会很低。

防止过拟合的策略有两种：

  • 预剪枝（pre-pruning）：及时停止树的生长，可以限制树的最大深度，限制叶节点的最大数目，或者规定一个节点中数据点的最小数目来防止继续划分。
  • 后剪枝（post-pruning）：先构造树，然后删除或折叠信息量很少的节点。

sklearn的决策树在DecisionTreeRegressor和DecisionTreeClassifier类中实现，但只实现了预剪枝。

决策树在cancer数据集上的表现：

1
2
3
4
5
6
7
8
9
10
from sklearn.tree import DecisionTreeClassifier
from sklearn.datasets import load_breast_cancer
from sklearn.model_selection import train_test_split

cancer = load_breast_cancer()
X_train, X_test, y_train, y_test = train_test_split(cancer.data, cancer.target, stratify=cancer.target, random_state=42)
tree = DecisionTreeClassifier(random_state=0)
tree.fit(X_train, y_train)
print("Accuracy on training set: {:.3f}".format(tree.score(X_train, y_train)))
print("Accuracy on test set: {:.3f}".format(tree.score(X_test, y_test)))

image-20210425145940981

可以看出，训练集准确率为100%，出现了过拟合。如果不限制决策树的深度，它的深度和复杂度会非常大，因此未剪枝的树容易出现过拟合，对新数据的泛化能力不佳。

可以设置参数max_depth的值，来限制树的最大深度：

1
tree = DecisionTreeClassifier(max_depth=4,random_state=0)

image-20210425150206666

可以看出，训练集下降了，但是测试集的准确率上升，泛化能力得到了增强。

分析决策树——树的可视化

1
2
3
4
5
6
rom sklearn.tree import export_graphviz
export_graphviz(tree,out_file="tree.dot",class_names=["magligant","benign"],feature_names=cancer.feature_names,impurity=False,filled=True)
import graphviz
with open("tree.dot") as f:
    dot_graph=f.read()
graphviz.Source(dot_graph)

树的特征重要性——fearture importance

特征重要性：决策树为每个特征对树的决策的重要性进行排序，是一个0到1之间的数字，其中0表示根本没有用到，1表示完美预测目标值，特征重要性之和始终为1。

1
print("Feature importance: \n {}".format(tree.feature_importances_))

image-20210425152037005

柱状图可视化：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
from sklearn.tree import DecisionTreeClassifier
from sklearn.datasets import load_breast_cancer
from sklearn.model_selection import train_test_split
import matplotlib.pyplot as plt
import numpy as np
cancer = load_breast_cancer()
X_train, X_test, y_train, y_test = train_test_split(cancer.data, cancer.target, stratify=cancer.target, random_state=42)
tree = DecisionTreeClassifier(max_depth=4,random_state=0)
tree.fit(X_train, y_train)
def plot_feature_importances_cancer(model):
    n_features = cancer.data.shape[1]
    plt.barh(range(n_features), model.feature_importances_, align='center')
    plt.yticks(np.arange(n_features), cancer.feature_names)
    plt.xlabel("Feature importance")
    plt.ylabel("Feature")
plot_feature_importances_cancer(tree)
plt.show()

image-20210425152649499

可以看出，顶部划分用到的特征（worst radius）是最重要的特征，第一层划分已经将两个类别区分得很好。

特征重要性告诉我们worst radius特征最重要，但是并没有告诉我们半径大表示样本是良性还是恶性，在特征和类别之间没有这样简单的正比反比关系。

用于回归任务的决策树

注意：==基于树的回归模型不能外推，也不能在训练数据范围之外进行预测。==

示例：以下为ram历史价格的数据集：

1
2
3
4
5
6
import pandas as pd
ram_prices = pd.read_csv("data/ram_price.csv")
plt.semilogy(ram_prices.date, ram_prices.price)
plt.xlabel("Year")
plt.ylabel("Price in $/Mbyte")
plt.show()

image-20210425153508602

我们将利用线性回归模型和决策树模型和2000年以前的历史价格来预测2000年以后的价格。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
import pandas as pd
from sklearn.tree import DecisionTreeRegressor
from sklearn.linear_model import LinearRegression
ram_prices = pd.read_csv("data/ram_price.csv")
# 利用历史数据策瑜2000年以后的ram价格
data_train = ram_prices[ram_prices.date < 2000 ]
data_test = ram_prices[ram_prices.date >= 2000]

# 基于日期来预测价格
X_train = data_train.date[:,np.newaxis]

# 利用对数变换来得到数据和目标之间更简单的关系
y_train = np.log(data_train.price)

tree = DecisionTreeRegressor().fit(X_train, y_train)
linear_reg = LinearRegression().fit(X_train, y_train)

# 对所有数据进行预测
X_all = ram_prices.date[:,np.newaxis]

pred_tree = tree.predict(X_all)
pred_lr = linear_reg.predict(X_all)

# 对数变换逆运算
price_tree = np.exp(pred_tree)
price_lr = np.exp(pred_lr)

# 绘制对比图像
plt.semilogy(data_train.date, data_train.price, label="Training data")
plt.semilogy(data_test.date, data_test.price, label="Test data")
plt.semilogy(ram_prices.date, price_tree, label="Tree prediction")
plt.semilogy(ram_prices.date, price_lr, label="linear prediction")
plt.legend()
plt.show()

image-20210425154538978

可以看出，两个模型的差异非常明显。对于线性模型，它丢失了训练集的某些信息，但是泛化能力很强，对2000年以后的价格预测较好。对于决策树模型，在训练集的预测上，与原始数据完全一致（图中蓝线和绿线重合在了一起），但是一旦超出训练数据的范围，即对2000年以后的价格预测上，则只能持续预测最后一个已知数据点。树不能在训练数据的范围之外生成新的响应，是所有基于树的回归模型的缺点。

决策树的优缺点

  • 优点：得到的模型很容易可视化；算法完全不受数据缩放的影响，不需要数据预处理，例如归一化或者标准化。
  • 缺点：即使做了预剪枝，也经常会过拟合，泛化能力很差。实践中常常使用下一节介绍的决策树集成来代替单一的决策树模型。
  • 重要参数：预剪枝参数，例如max_depth,max_leaf_nodes,max_samples_leaf等。

2.3.6 决策树集成

集成（ensemble）是合并多个机器学习模型来构建更强大模型的方法，已证明有两种集成模型对大量分类和回归的数据集都是有效的，二者都以决策树为基础，分别是随机森林和梯度提升决策树。

随机森林——random forest

  • 思想：每棵树的预测都可能相对较好，但可能对部分数据存在过拟合，所以可以构造许多棵树，对这些树的结果取平均值来降低过拟合。

在sklearn中，随机森林在RandomForestRegressor和RandomForestClassifier中实现，首先要确定用于构造的树的个数，取决于参数n_estimators

  • 树的随机化有两种：

  1. 通过选择用于构造树的数据点：首先要对数据进行自助采样（bootstrap sample），从n_samples个数据点中有放回地重复随机抽取样本。共抽取n_samples次，这样会创建一个与原数据集想用大小的数据集，但是有些数据点（大约三分之一）会被丢失，另一些则会重复。例如：创建数据集[a,b,c,d]的自助采样，一种可能的是采样是[a,a,b,b]，或者[d,c,,d,a]等等。

  2. 通过选择每次划分测试的特征：在每个决策节点处，算法随机选择特征的一个子集（特征子集），并对其中一个特征寻找最佳测试。选择特征的个数（特征子集的大小）由参数max_features决定。

由以上两种方法，随机森林中构造每棵树的数据集都是略有不同的，由于每个节点的特征选择，每棵树中的每次划分都是基于特征的不同子集，这两种随机化方法保证随机森林中的每棵树都是不同的。

  • max_features：关键参数。如果设置max_features=n_featrues，那么每次划分都会考虑数据集的所有特征，在特征选择的过程中就没有添加随机性。如果设置为1，那么划分时将无法选择对哪个特征进行测试。max_features越大，森林中的树会越相似，越小则差异越大。

  • two_moons数据集介绍：这个数据集有两个类别，每个类别包含50个数据点。

  • 分析随机森林：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
from sklearn.datasets import make_moons
from sklearn.ensemble import RandomForestClassifier
from sklearn.model_selection import train_test_split
import matplotlib.pyplot as plt
import mglearn
X, y = make_moons(n_samples=100, noise=0.25, random_state=3)
X_train, X_test, y_train, y_test = train_test_split(X, y, stratify=y, random_state=42)

forest = RandomForestClassifier(n_estimators=5, random_state=2)
forest.fit(X_train, y_train)
fig, axes = plt.subplots(2,3,figsize=(20,10))
for i,(ax,tree) in enumerate(zip(axes.ravel(), forest.estimators_)):
    ax.set_title("Tree {}".format(i))
    mglearn.plots.plot_tree_partition(X_train, y_train, tree, ax=ax)
mglearn.plots.plot_2d_separator(forest, X_train, fill=True, ax=axes[-1,-1], alpha=.4)
axes[-1,-1].set_title("Random Forest")
mglearn.discrete_scatter(X_train[:,0], X_train[:,1], y_train)
plt.show()

image-20210425165223411

可以看出，每棵树都犯了一些错误，因为这里画出的一些训练点实际上并没有包含在这些树的训练集中，原因在于自助采样。

在cancer数据集上的表现（n_estimators=100）：

1
2
3
4
5
6
cancer = load_breast_cancer()
X_train, X_test, y_train, y_test = train_test_split(cancer.data, cancer.target, random_state=0)
forest = RandomForestClassifier(n_estimators=100, random_state=2)
forest.fit(X_train, y_train)
print("Accuracy on training set: {:.3f}".format(forest.score(X_train, y_train)))
print("Accuracy on test set: {:.3f}".format(forest.score(X_test, y_test)))

image-20210425165701974

  • 与决策树相似，随机森林也给出了特征重要性，一般来说，后者的可靠性要比单棵树的要好。

image-20210425165926359

此处，随机森林认为worst perimeter是信息量最大的特征。

  • 优缺点：对于维度非常高的稀疏数据（例如文本数据），随机森林表现不是很好，这种情况通常采用线性模型处理。优点是鲁棒性很好，不需要反复调节参数，也不需要数据预处理就可以给出很好的结果。
  • 重要参数：
  1. n_estimators：总是越大越好，鲁棒性越强，不过收益是递减的，对计算机硬件的要求就越高，训练时间也越长。
  2. max_features：决定每棵树的随机性大小，最好使用默认值，为sqrt(n_features)；对于回归，则为n_features；有时调整该参数可能提高性能。
  3. 预剪枝参数：预剪枝参数可能影响随机森林的性能。

梯度提升决策树（梯度提升机）——gradient boosted decision tree

梯度提升采用连续的方式构造树，每棵树都视图纠正前一棵树的错误。默认情况下，梯度提升回归树没有随机化，而是用到了预强剪枝。梯度提升树通常采用深度很小（1到5之间）的树。

  • learn_rate：学习率。用于控制每棵树纠正前一棵树的错误的强度。越大，意味着做出强的修正，模型更为复杂。

  • 在cancer数据集上的表现：（默认100棵树，最大深度3，学习率0.1）

1
2
3
4
5
6
7
from sklearn.ensemble import GradientBoostingClassifier
cancer = load_breast_cancer()
X_train, X_test, y_train, y_test = train_test_split(cancer.data, cancer.target, random_state=0)
gbrt = GradientBoostingClassifier(random_state=0)
gbrt.fit(X_train, y_train)
print("Accuracy on training set: {:.3f}".format(gbrt.score(X_train, y_train)))
print("Accuracy on test set: {:.3f}".format(gbrt.score(X_test, y_test)))

image-20210426161050813

出现了过拟合，可以限制最大深度，或者降低学习率：

1
gbrt = GradientBoostingClassifier(random_state=0, max_depth=1)

image-20210426161154819

1
gbrt = GradientBoostingClassifier(random_state=0, learning_rate=0.01)

image-20210426161218159

  • 特征重要性：

image-20210426161439175

  • 优点：鲁棒性强；监督学习中最强大最常用的模型之一
  • 缺点：需要仔细调参，训练时间也可能比较长

2.3.7 核支持向量机

核支持向量机：kernelized support vector machine

可用于分类和回归，分类在SVC中实现，回归在SVR中实现。

线性模型和非线性特征

线性模型在低维空间中可能非常受限，因为线和平面的灵活性有限。

对二分类blobs数据集：

1
2
3
4
5
6
X, y = make_blobs(centers=4, random_state=8)
y = y % 2
mglearn.discrete_scatter(X[:,0],X[:,1],y)
plt.xlabel("Feature 0")
plt.ylabel("Feature 1")
plt.show()

image-20210426162208420

这样的类别并不是线性可分的，如果采用一条直线的线性模型来划分数据点，将无法给出较好的结果：

1
2
3
4
5
6
7
8
9
10
# 引入线性支持向量机
from sklearn.svm import LinearSVC
X, y = make_blobs(centers=4, random_state=8)
y = y % 2
linear_svm = LinearSVC().fit(X, y)
mglearn.plots.plot_2d_separator(linear_svm, X)
mglearn.discrete_scatter(X[:,0],X[:,1],y)
plt.xlabel("Feature 0")
plt.ylabel("Feature 1")
plt.show()

image-20210426162547187

现在对输入特征进行扩展，例如添加第二个特征的平方作为一个新特征，将每个数据点表示为三维点。

2.3.8 神经网络（深度学习）

第三章 无监督学习

无监督学习包括没有已知输出，没有老师指导学习算法的各种机器学习，在无监督学习中，学习算法只有输入数据，并需要从这些数据中提取知识。

3.1 无监督学习的类型

主要研究两种类型的无监督学习：

  • 数据集变换
  • 聚类

常见应用：降维：接受包括许多特征的数据的高维表示，并找到表示该数据的一种新方法，用较少的特征就可以概括其重要特性；聚类：将数据划分成不同的组，每组包括相似的物项。作为监督算法的预处理步骤

3.2 无监督学习的挑战

无监督学习算法常用于不包含任何标签信息的数据，评估无监督算法结果的唯一方法就是人工检查。

3.3 预处理和缩放

通常来说，这是对数据的一种简单的按特征的缩放和移动。

3.3.1 不同类型的预处理

image-20220414201612681

左图显示的是一个有两个特征的二分类数据集，第一个特征x轴位于10-15之间，第二个特征y轴大约位于1-9之间。

StandardScaler确保每个特征的平均值为0，方差为1，使所有特征都位于同一量级。

RobustScaler与上者类似，但是利用中位数和四分位数。

MinMaxScaler移动数据，使所有特征都刚好位于0-1之间。

Normalizer对每个数据点进行缩放，使得特征向量的欧式长度为1，换言之，它将一个数据点投射到半径为1的圆上（或者球面），这意味着每个数据点的缩放比例都不同。

3.3.2 应用数据变换

将SVM应用在cancer数据集上，并使用MinMaxScaler来预处理数据。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
from sklearn.datasets import load_breast_cancer
from sklearn.model_selection import train_test_split
from sklearn.preprocessing import MinMaxScaler
cancer = load_breast_cancer()
X_train, X_test, y_train, y_test = train_test_split(cancer.data, cancer.target, stratify=cancer.target, random_state=1)
scaler = MinMaxScaler()
# 使用fit方法拟合缩放器scaler，并将其应用于训练数据
# fit方法计算每个特征的最大值和最小值
scaler.fit(X_train)
# 进行实际缩放
# 使用transform函数，得到模型返回数据的新表示
# 变换数据：
X_train_scaled = scaler.transform(X_train)
# 打印出缩放前后的数据集属性
print("transformed shape: {}".format(X_train_scaled.shape))
print("per-feature minimum before scaling:\n {}".format(X_train.min(axis=0)))
print("per-feature maximum before scaling:\n {}".format(X_train.max(axis=0)))
print("per-feature mininum after scaling:\n {}".format(X_train_scaled.min(axis=0)))
print("per-feature mininum after scaling:\n {}".format(X_train_scaled.max(axis=0)))

image-20220414201637140

可以看出没有改变数据集的形状，但是所有特征的值都位于0-1之间。

然后对测试集进行缩放：

1
2
3
4
X_test_scaled = scaler.transform(X_test)
# 打印出缩放后的测试数据集属性
print("per-feature mininum after scaling:\n {}".format(X_test_scaled.min(axis=0)))
print("per-feature mininum after scaling:\n {}".format(X_test_scaled.max(axis=0)))

image-20220414201658596

可以看出对测试集进行缩放后其最大最小值并不是0和1，这是因为，所有的缩放器总是对训练集和测试集进行相同的变换，也就是说，transform方法总是减去训练集的最小值，然后除以训练集的范围，这两个值可能与测试集中的不同。

3.3.3 对训练数据和测试数据进行相同的缩放

略

快捷方式与高效的替代方法

通常来说，你需要在某个数据集上拟合（fit）一个模型，然后再将其transform。所有具有transform方法的模型都具有fit_transform方法：

1
2
3
4
5
6
from sklearn.preprocessing import StandardScaler
scaler = StandardScaler()
# 依次调用fit和transform（使用方法链）
X_scaled = scaler.fit(X).transform(X)
# 方法2
X_scaled = scaler.fit_transform(X)

3.3.4 预处理对监督学习的作用

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
from sklearn.datasets import load_breast_cancer
from sklearn.model_selection import train_test_split
from sklearn.preprocessing import MinMaxScaler
from sklearn.svm import SVC
cancer = load_breast_cancer()
X_train, X_test, y_train, y_test = train_test_split(cancer.data, cancer.target, random_state=0)
svm = SVC(C=100)
svm.fit(X_train, y_train)
print("Test set accuracy: {:.2f}".format(svm.score(X_test, y_test)))
# 对数据进行预处理后的准确率：
scaler = MinMaxScaler()
scaler.fit(X_train)
X_train_scaled = scaler.transform(X_train)
X_test_scaled = scaler.transform(X_test)
svm.fit(X_train_scaled, y_train)
print("Scaled test set accuracy: {:.2f}".format(svm.score(X_test_scaled, y_test)))

image-20220414201711287

可见数据经过预处理后，准确率提高了。

3.4 降维，特征提取与流形学习

利用无监督学习进行数据变换可能有多种目的，最常见的目的就是可视化，压缩数据，以及寻找信息量更大的数据表示以用于进一步的处理。常见算法有主成分分析，非负矩阵分解和t-SNE。

3.4.1主成分分析PCA

简介

PCA是一种基本的数据降维技术。 PCA是一种著名的数据降维算法，它应用的条件是数据/特征之间具有明显的线性相关性，它的两个主要的指导思想是抓主要矛盾和方差即信息，它的基本应用是数据降维，以此为基础还有数据可视化、数据压缩存储、异常检测、特征匹配与距离计算等。从数学上理解，它是一种矩阵分解算法；从物理意义上理解，它是线性空间上的线性变换。

主成分分析：principal component analysis。是一种旋转数据集的方法，旋转后的特征在统计上不相关。在做完这种旋转后，通常是根据新特征对解释数据的重要性来选择它的一个子集。

image-20220414201722970

  • 左上图显示的是原始数据点。算法首先找到方差最大的方向，将其标记为“成分1”（component 1）。这是数据中包含信息最多的方向（或向量），沿着这个方向的特征之间最为相关。然后算法找到与第一个方向正交且包含信息最多的方向，在二维空间中只有一个，但在高维空间中则有许多个。这些方向被称为主成分，一般来说，主成分个数与原始特征相同。

  • 右上图将其左上图旋转得来，使得成分1与x轴平行，成分2与y轴平行。在旋转之前，从数据中减去平均值，使得变换后的数据以0为中心。

  • 可以通过仅保留一部分主成分来使用pca进行降维，对于左下图，保留成分1，使得数据变成1维数据集，但要注意，并不是保留了原始特征中的某一个，而是保留了主成分的某一方向。
  • 右下图由左下图反向旋转并重新加回平均值得来，这种变换有时用于去除数据中的噪声影响，或将主成分中保留的那部分信息可视化。

原理概述

见打印资料

将PCA应用于cancer数据集并可视化

PCA最常见的应用之一是将高维数据集可视化。

在应用PCA之前，利用StandardScaler缩放数据。

1
2
3
4
5
6
from sklearn.datasets import load_breast_cancer
from sklearn.preprocessing import StandardScaler
cancer = load_breast_cancer()

scaler = StandardScaler()
X_scaled = scaler.fit_transform(cancer.data)

首先将PCA对象实例化，然后调用fit找到主成分，然后调用transform来旋转并降维，在这个过程中默认保留所有主成分，可以指定保留的个数。

1
2
3
4
5
6
7
8
9
10
from sklearn.decomposition import PCA
# 保留数据的前两个主成分
pca = PCA(n_components=2)
# 对cancer数据集拟合PCA模型
pca.fit(X_scaled)

# 将数据变换到前两个主成分的方向上
X_pca = pca.transform(X_scaled)
print("Original shape: {}".format(str(X_scaled.shape)))
print("Reduced shape: {}".format(str(X_pca.shape)))

image-20220414201736852

然后利用前两个主成分绘制cancer数据集的二维散点图：

1
2
3
4
5
6
7
plt.figure(figsize=(8,8))
mglearn.discrete_scatter(X_pca[:,0], X_pca[:,1],cancer.target)
plt.legend(cancer.target_names,loc="best")
plt.gca().set_aspect("equal")
plt.xlabel("First principal component")
plt.ylabel("Second principal component")
plt.show()

image-20220414201747284

注意，pca是无监督算法，在寻找旋转方向时没有用到任何类别信息，但是可以看出两个类别被很好地区分开来了。

pca的缺点是解释性不好，主成分对应于原始数据中的方向，是原始特征的组合，这些组合往往非常复杂。

在拟合过程中，主成分被保存在PCA对象的components_属性中：

1
2
print("PCA component shape: {}".format(pca.components_.shape))
print("PCA components:\n{}".format(pca.components_))

image-20220414201759498

components_中的每一行对应于一个主成分，他们按照重要性排序（第一主成分排在第一位，以此类推），列对应于PCA的原始特征属性，共30个，例如mean radius，mean texture等。

用热图进行可视化：

1
2
3
4
5
6
7
plt.matshow(pca.components_, cmap='viridis')
plt.yticks([0,1],["First component","Second component"])
plt.colorbar()
plt.xticks(range(len(cancer.feature_names)), cancer.feature_names, rotation=60, ha='left')
plt.xlabel("Feature")
plt.ylabel("Principal components")
plt.show()

image-20220414201812234

这些所有特征的混合使得解释坐标轴含义变得十分困难。

PCA用于人脸识别

PCA的另一个应用是特征提取。

Wild数据集一共有3023张图像，每张大小为87*65像素，分别属于62个不同的人。但这个数据集有数据偏斜，例如乔治布什的有近530张图像，为降低数据偏斜，我们对每个人最多只取50张图像。

1
2
3
4
5
6
7
8
9
10
from sklearn.datasets import fetch_lfw_people
import matplotlib.pyplot as plt
import numpy as np
people = fetch_lfw_people(min_faces_per_person=20, resize=0.7)
image_shape = people.images[0].shape
fix, axes = plt.subplots(2, 5, figsize=(15,8), subplot_kw={'xticks': (), 'yticks': ()})
for target, image, ax in zip(people.target, people.images, axes.ravel()):
    ax.imshow(image)
    ax.set_title(people.target_names[target])
plt.show()

image-20220414201824757

1
2
3
4
5
6
7
8
9
mask = np.zeros(people.target.shape, dtype=np.bool_)
for target in np.unique(people.target):
    mask[np.where(people.target == target)[0][:50]] = 1
X_people = people.data[mask]
y_people = people.target[mask]

# 将灰度值缩放到0-1之间，而不是在0-255之间
# 以得到更好的数据稳定性
X_people = X_people / 255.

利用单一近邻算法对人脸进行预测：

1
2
3
4
5
6
7
from sklearn.neighbors import KNeighborsClassifier
from sklearn.model_selection import train_test_split
X_train, X_test, y_train, y_test = train_test_split(X_people, y_people, stratify=y_people,random_state=0)
# 使用1个邻居构建分类器
knn = KNeighborsClassifier(n_neighbors=1)
knn.fit(X_train, y_train)
print("Test set score of 1-nn: {:.2f}".format(knn.score(X_test, y_test)))

image-20220414201836271

利用PCA进行数据预处理：这里可以使用PCA的白化选项

image-20220414201848526

上图是经过PCA白化（数据预处理）后的图像。

现在对训练数据拟合PCA对象，并提取前100个主成分，然后对训练数据和测试数据进行变换：

1
2
3
4
5
6
7
from sklearn.decomposition import PCA
# 调用方法链，启用白化选项
pca = PCA(n_components=100, whiten=True ,random_state=0).fit(X_train)
# 进行数据变换
X_train_pca = pca.transform(X_train)
X_test_pca = pca.transform(X_test)
print("X_train_pca.shape: {}".format(X_train_pca.shape))

image-20220414201858134

可以看出新数据集有100个新特征，即前100个主成分。现在利用k近邻算法进行分类：

1
2
3
knn = KNeighborsClassifier(n_neighbors=1)
knn.fit(X_train_pca,y_train)
print("Test set score of 1-nn: {:.2f}".format(knn.score(X_test_pca, y_test)))

image-20220414201910291

可见准确率得到了提升。

现在查看前几个主成分：

1
2
3
4
5
fix, axes = plt.subplots(3, 5, figsize=(15,12), subplot_kw={'xticks': (), 'yticks': ()})
for i, (component, ax) in enumerate(zip(pca.components_, axes.ravel())):
    ax.imshow(component.reshape(image_shape), cmap='viridis')
    ax.set_title("{}. component".format((i+1)))
plt.show()

image-20220414202415022

理解主成分很难，但是可以猜测，例如左上图似乎主要编码的是人脸与背景的对比，第二主成分编码的可能是脸左半部分和右半部分的明暗程度差异。

3.4.2 非负矩阵分解NMF

非负矩阵分解：non-negative matrix factorization。目的在于提取有用特征，这种方法只能应用于每个特征都是非负的数据，于PCA相比，NMF得到的分量更容易解释。

NMF原理简介

对于任意给定的一个非负矩阵V，其能够寻找到一个非负矩阵W和一个非负矩阵H，满足条件V=WH,从而将一个非负的矩阵分解为左右两个非负矩阵的乘积。*其中，V矩阵中每一列代表一个观测(observation)，每一行代表一个特征(feature)；W矩阵称为基矩阵，H矩阵称为系数矩阵或权重矩阵。这时用系数矩阵H代替原始矩阵，就可以实现对原始矩阵进行降维，得到数据特征的降维矩阵，从而减少存储空间。过程如下图所示：

image-20220414202452142

其物理意义为每个基向量W(如眼睛、鼻子等等)乘以相应权重H最终组合而成V(脸)。虽然NMF是一个很厉害的算法，但其实质是加权和，我们可以在原理上等效为基本的线性方程。

sklearn中封装的代码

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
from sklearn.decomposition import NMF
from sklearn.datasets import load_iris


X, _ = load_iris(True) # X为原始矩阵V

# can be used for example for dimensionality reduction, source separation or topic extraction
# 个人认为最重要的参数是n_components、alpha、l1_ratio、solver
nmf = NMF(n_components=2,  # k value,默认会保留全部特征
          init=None,  # W H 的初始化方法，包括'random' | 'nndsvd'(默认) |  'nndsvda' | 'nndsvdar' | 'custom'.
          solver='cd',  # 'cd' | 'mu'
          beta_loss='frobenius',  # {'frobenius', 'kullback-leibler', 'itakura-saito'}，一般默认就好
          tol=1e-4,  # 停止迭代的极限条件
          max_iter=200,  # 最大迭代次数
          random_state=None,
          alpha=0.,  # 正则化参数
          l1_ratio=0.,  # 正则化参数
          verbose=0,  # 冗长模式
          shuffle=False  # 针对"cd solver"
          )

# -----------------函数------------------------
print('params:', nmf.get_params())  # 获取构造函数参数的值，也可以nmf.attr得到，所以下面我会省略这些属性

# 下面四个函数很简单，也最核心，例子中见
nmf.fit(X)
W = nmf.fit_transform(X)
W = nmf.transform(X)
nmf.inverse_transform(W)
# -----------------属性------------------------
H = nmf.components_  # H矩阵 系数矩阵
print('reconstruction_err_', nmf.reconstruction_err_)  # 损失函数值
print('n_iter_', nmf.n_iter_)  # 实际迭代次数

将NMF应用于模拟数据

image-20220414202515413

以上的图给出了NMF在二维模拟数据上的结果。对于两个分量的NMF，显然所有的数据点都可以写成这两个分量的正数组合，如果有足够的分量能够完美地重建数据（分量个数等于特征个数），那么算法会指向数据极值的方向。

NMF的主要参数是我们想要提取的分量个数。

将NMF应用于信号源区分

假设有一个信号由三个不同的信号源合成的。

1
2
3
4
5
6
7
8
9
import mglearn
import matplotlib.pyplot as plt

S = mglearn.datasets.make_signals()
plt.figure(figsize=(6,1))
plt.plot(S, '-')
plt.xlabel("Time")
plt.ylabel("Signal")
plt.show()

image-20220414202530024

我们想要将混合信号分解为原始分量，假设有100种方法来观测混合信号，每种方法都为我们提供了一系列的观测结果。

1
2
3
4
# 将数据混合成100维的状态
A = np.random.RandomState(0).uniform(size=(100,3))
X = np.dot(S, A.T)
print("Shape of measurements: {}".format(X.shape))

image-20220414202543634

我们可以用NMF来还原这三个信号：

1
2
3
nmf = NMF(n_components=3, random_state=42)
S_ = nmf.fit_transform(X)
print("Recovered signal shape: {}".format(S_.shape))

image-20220414202555098

添加了PCA方法用于比较：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
# 将数据混合成100维的状态
A = np.random.RandomState(0).uniform(size=(100,3))
X = np.dot(S, A.T)
nmf = NMF(n_components=3, random_state=42)
S_ = nmf.fit_transform(X)
pca = PCA(n_components=3)
H = pca.fit_transform(X)
models = [X, S, S_, H]
names = ['Observations (first three measurements)', 'True sources', 'NMF recovered signals', 'PCA recovered signals']
fig, axes = plt.subplots(4, figsize=(8,4), gridspec_kw={'hspace': .5}, subplot_kw={'xticks':(), 'yticks':()})
for model, name, ax in zip(models, names, axes):
    ax.set_title(name)
    ax.plot(model[:, :3],'-')
plt.show()

image-20220414202606483

可以看出NMF在发现原始信号源时得到了不错的结果，而PCA则失败了，仅使用第一个成分来解释数据中的大部分变化。

3.4.3 用t-SNE进行流形学习

简介

有一类用于可视化的算法叫做流形学习算法，其中特别有用的一个是t-SNE算法。

t-SNE全称为t-distributed Stochastic Neighbor Embedding，翻译为t-随机邻近嵌入，它是一种嵌入模型，能够将高维空间中的数据映射到低维空间中，并保留数据集的局部特性，该算法在论文中非常常见，主要用于高维数据的降维和可视化。

t-SNE可以算是目前效果最好的数据降维和可视化方法之一，当我们想对高维数据集进行分类，但又不清楚这个数据集有没有很好的可分性（同类之间间隔小、异类之间间隔大）时，可以通过t-SNE将数据投影到2维或3维空间中观察一下：如果在低维空间中具有可分性，则数据是可分的；如果在低维空间中不可分，则可能是因为数据集本身不可分，或者数据集中的数据不适合投影到低维空间。

t-SNE将数据点之间的相似度转化为条件概率，原始空间中数据点的相似度由高斯联合分布表示，嵌入空间中数据点的相似度由学生t分布表示。通过原始空间和嵌入空间的联合概率分布的KL散度（用于评估两个分布的相似度的指标，经常用于评估机器学习模型的好坏）来评估嵌入效果的好坏，即将有关KL散度的函数作为损失函数（loss function），通过梯度下降算法最小化损失函数，最终获得收敛结果。要注意t-SNE的缺点很明显：占用内存较多、运行时间长。

此外：流形学习算法主要用于可视化，对t-SNE算法，不允许变换新数据，这意味着不能用于测试集，只能变换用于训练的数据。

应用1——降维

输入4个5维的数据，通过t-SNE将其降维成2维的数据，代码如下：

1
2
3
4
5
6
7
8
9
10
11
12
13
import numpy as np
from sklearn.manifold import TSNE

# 将3维数据降维2维

# 4个3维的数据
x = np.array([[0, 0, 0, 1, 2], [0, 1, 1, 3, 5], [1, 0, 1, 7, 2], [1, 1, 1, 10, 22]])
# 嵌入空间的维度为2，即将数据降维成2维
ts = TSNE(n_components=2)
# 训练模型
ts.fit_transform(x)
# 打印结果
print(ts.embedding_)

image-20220414202817135

应用2——S型曲线的降维与可视化

S型曲线中的数据是高维的数据，不同的颜色表示不同的数据点。当我们通过t-SNE将数据嵌入到2维空间中后，可以看到数据点之间的类别信息被完整地保留了下来。代码如下：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
import matplotlib.pyplot as plt
from sklearn import manifold, datasets

# 生成1000个S型曲线数据
x, color = datasets._samples_generator.make_s_curve(n_samples=1000, random_state=0)		# x是[1000,2]的2维数据，color是[1000,1]的一维数据

n_neighbors = 10
n_components = 2

# 创建自定义图像
fig = plt.figure(figsize=(8, 8))		# 指定图像的宽和高
plt.suptitle("Dimensionality Reduction and Visualization of S-Curve Data ", fontsize=14)		# 自定义图像名称

# 绘制S型曲线的3D图像
ax = fig.add_subplot(211, projection='3d')		# 创建子图
ax.scatter(x[:, 0], x[:, 1], x[:, 2], c=color, cmap=plt.cm.Spectral)		# 绘制散点图，为不同标签的点赋予不同的颜色
ax.set_title('Original S-Curve', fontsize=14)
ax.view_init(4, -72)		# 初始化视角

# t-SNE的降维与可视化
ts = manifold.TSNE(n_components=n_components, init='pca', random_state=0)
# 训练模型
y = ts.fit_transform(x)
ax1 = fig.add_subplot(2, 1, 2)
plt.scatter(y[:, 0], y[:, 1], c=color, cmap=plt.cm.Spectral)
ax1.set_title('t-SNE Curve', fontsize=14)
# 显示图像
plt.show()

image-20220414202828161

应用3——手写数字数据集的降维与可视化

手写数字数据集是一个经典的图片分类数据集，数据集中包含0-9这10个数字的灰度图片，每张图片以8*8共64个像素点表示。

首先利用PCA降维到二维可视化。对前两个主成分作图，并按类别对数据点着色。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
from logging import disable
import matplotlib.pyplot as plt
from sklearn.decomposition import PCA
from sklearn.datasets import load_digits
digits = load_digits()

# 构建一个PCA模型
pca = PCA(n_components=2)
pca.fit(digits.data)
# 将digits数据变换到前两个主成分的方向上
digits_pca = pca.transform(digits.data)
colors = ["#476A2A","#7851B8","#BD3430","#4A2D4E","#875525","#A83683","#4E655E","#853541","#3A3120","#535D8E"]
plt.figure(figsize=(10,10))
plt.xlim(digits_pca[:,0].min(), digits_pca[:,0].max())
plt.ylim(digits_pca[:,1].min(), digits_pca[:,1].max())
for i in range(len(digits.data)):
    plt.text(digits_pca[i,0], digits_pca[i,1], str(digits.target[i]), color=colors[digits.target[i]], fontdict={'weight': 'bold', 'size':9})
plt.xlabel("First principal component")
plt.ylabel("Second principal component")
plt.show()

image-20220414202841091

可以看出各个类别的界限不太好。现在利用t-sne算法进行降维并可视化，由于其算法不支持变换新数据，即没有transform方法，但可以调用fit_transform方法，其会构建模型并立即返回变换后的数据。

1
2
3
4
5
6
7
8
9
10
11
from sklearn.manifold import TSNE
tsne = TSNE(random_state=42)
digits_tsne = tsne.fit_transform(digits.data)
plt.figure(figsize=(10,10))
plt.xlim(digits_tsne[:,0].min(), digits_tsne[:,0].max() + 1)
plt.ylim(digits_tsne[:,1].min(), digits_tsne[:,1].max() + 1)
for i in range(len(digits.data)):
    plt.text(digits_tsne[i,0], digits_tsne[i,1], str(digits.target[i]), color=colors[digits.target[i]], fontdict={'weight': 'bold', 'size':9})
plt.xlabel("t-SNE feature 0")
plt.ylabel("t-SNE feature 1")
plt.show()

image-20220414202858313

所有类别都被明确分开，要记住，这种方法并不知道类别标签，它完全是无监督的，但它能够找到数据的一种二维表示，仅根据原始空间中数据点之间的靠近程度就能够将各个类别明确分开。

3.5 聚类

聚类（clustering）是将数据集划分成组的任务，这些组称为簇。

分类和聚类

分类：分类其实是从特定的数据中挖掘模式，作出判断的过程。比如Gmail邮箱里有垃圾邮件分类器，一开始的时候可能什么都不过滤，在日常使用过程中，我人工对于每一封邮件点选“垃圾”或“不是垃圾”，过一段时间，Gmail就体现出一定的智能，能够自动过滤掉一些垃圾邮件了。这是因为在点选的过程中，其实是给每一条邮件打了一个“标签”，这个标签只有两个值，要么是“垃圾”，要么“不是垃圾”，Gmail就会不断研究哪些特点的邮件是垃圾，哪些特点的不是垃圾，形成一些判别的模式，这样当一封信的邮件到来，就可以自动把邮件分到“垃圾”和“不是垃圾”这两个我们人工设定的分类的其中一个。

聚类：聚类的目的也是把数据分类，但是事先我是不知道如何去分的，完全是算法自己来判断各条数据之间的相似性，相似的就放在一起。在聚类的结论出来之前，我完全不知道每一类有什么特点，一定要根据聚类的结果通过人的经验来分析，看看聚成的这一类大概有什么特点。

聚类和分类最大的不同在于：分类的目标是事先已知的，而聚类则不一样，聚类事先不知道目标变量是什么，类别没有像分类那样被预先定义出来。

3.5.1 k均值聚类

聚类算法有很多种（几十种），K-Means是聚类算法中的最常用的一种，算法最大的特点是简单，好理解，运算速度快，但是只能应用于连续型的数据，并且一定要在聚类前需要手工指定要分成几类。

算法步骤

下面，我们描述一下K-means算法的过程，为了尽量不用数学符号，所以描述的不是很严谨，大概就是这个意思，“物以类聚、人以群分”：

  1. 首先输入k的值，即我们希望将数据集经过聚类得到k个分组。
  2. 从数据集中随机选择k个数据点作为初始大哥（质心，Centroid）
  3. 对集合中每一个小弟，计算与每一个大哥的距离（距离的含义后面会讲），离哪个大哥距离近，就跟定哪个大哥。
  4. 这时每一个大哥手下都聚集了一票小弟，这时候召开人民代表大会，每一群选出新的大哥（其实是通过算法选出新的质心）。
  5. 如果新大哥和老大哥之间的距离小于某一个设置的阈值（表示重新计算的质心的位置变化不大，趋于稳定，或者说收敛），可以认为我们进行的聚类已经达到期望的结果，算法终止。
  6. 如果新大哥和老大哥距离变化很大，需要迭代3~5步骤。

image-20220414202912671

上图中，三角形表示簇中心，圆形表示数据点。我们指定了要寻找三个簇。下图为算法的决策边界：

image-20220414202922790

算法的简单实现

下面给出它的Python实现，其中中心点的选取是手动选择的。在代码中随机产生了一个样本，用于测试K-means算法：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
# K-means Algorithm is a clustering algorithm
import numpy as np
import matplotlib.pyplot as plt
import random
 
 
def get_distance(p1, p2):
    diff = [x-y for x, y in zip(p1, p2)]
    distance = np.sqrt(sum(map(lambda x: x**2, diff)))
    return distance
 
 
# 计算多个点的中心
# cluster = [[1,2,3], [-2,1,2], [9, 0 ,4], [2,10,4]]
def calc_center_point(cluster):
    N = len(cluster)
    m = np.matrix(cluster).transpose().tolist()
    center_point = [sum(x)/N for x in m]
    return center_point
 
 
# 检查两个点是否有差别
def check_center_diff(center, new_center):
    n = len(center)
    for c, nc in zip(center, new_center):
        if c != nc:
            return False
    return True
 
 
# K-means算法的实现
def K_means(points, center_points):
 
    N = len(points)         # 样本个数
    n = len(points[0])      # 单个样本的维度
    k = len(center_points)  # k值大小
 
    tot = 0
    while True:             # 迭代
        temp_center_points = [] # 记录中心点
 
        clusters = []       # 记录聚类的结果
        for c in range(0, k):
            clusters.append([]) # 初始化
 
        # 针对每个点，寻找距离其最近的中心点（寻找组织）
        for i, data in enumerate(points):
            distances = []
            for center_point in center_points:
                distances.append(get_distance(data, center_point))
            index = distances.index(min(distances)) # 找到最小的距离的那个中心点的索引，
 
            clusters[index].append(data)    # 那么这个中心点代表的簇，里面增加一个样本
 
        tot += 1
        print(tot, '次迭代   ', clusters)
        k = len(clusters)
        colors = ['r.', 'g.', 'b.', 'k.', 'y.']  # 颜色和点的样式
        for i, cluster in enumerate(clusters):
            data = np.array(cluster)
            data_x = [x[0] for x in data]
            data_y = [x[1] for x in data]
            plt.subplot(2, 3, tot)
            plt.plot(data_x, data_y, colors[i])
            plt.axis([0, 1000, 0, 1000])
 
        # 重新计算中心点（该步骤可以与下面判断中心点是否发生变化这个步骤，调换顺序）
        for cluster in clusters:
            temp_center_points.append(calc_center_point(cluster))
 
        # 在计算中心点的时候，需要将原来的中心点算进去
        for j in range(0, k):
            if len(clusters[j]) == 0:
                temp_center_points[j] = center_points[j]
 
        # 判断中心点是否发生变化：即，判断聚类前后样本的类别是否发生变化
        for c, nc in zip(center_points, temp_center_points):
            if not check_center_diff(c, nc):
                center_points = temp_center_points[:]   # 复制一份
                break
        else:   # 如果没有变化，那么退出迭代，聚类结束
            break
 
    plt.show()
    return clusters # 返回聚类的结果
 
# 随机获取一个样本集，用于测试K-means算法
def get_test_data():
 
    N = 1000
 
    # 产生点的区域
    area_1 = [0, N / 4, N / 4, N / 2]
    area_2 = [N / 2, 3 * N / 4, 0, N / 4]
    area_3 = [N / 4, N / 2, N / 2, 3 * N / 4]
    area_4 = [3 * N / 4, N, 3 * N / 4, N]
    area_5 = [3 * N / 4, N, N / 4, N / 2]
 
    areas = [area_1, area_2, area_3, area_4, area_5]
    k = len(areas)
 
    # 在各个区域内，随机产生一些点
    points = []
    for area in areas:
        rnd_num_of_points = random.randint(50, 200)
        for r in range(0, rnd_num_of_points):
            rnd_add = random.randint(0, 100)
            rnd_x = random.randint(area[0] + rnd_add, area[1] - rnd_add)
            rnd_y = random.randint(area[2], area[3] - rnd_add)
            points.append([rnd_x, rnd_y])
 
    # 自定义中心点，目标聚类个数为5，因此选定5个中心点
    center_points = [[0, 250], [500, 500], [500, 250], [500, 250], [500, 750]]
 
    return points, center_points
 
 
if __name__ == '__main__':
 
    points, center_points = get_test_data()
    clusters = K_means(points, center_points)
    print('#######最终结果##########')
    for i, cluster in enumerate(clusters):
        print('cluster ', i, ' ', cluster)

结果：

image-20220414202945532

sklearn中的算法使用

1
2
3
4
5
6
7
8
from sklearn.datasets import make_blobs
from sklearn.cluster import KMeans
X, y = make_blobs()

# 构建聚类模型
# 指定簇的个数为3，默认为8
kmeans = KMeans(n_clusters=3)
kmeans.fit(X)

算法会为X中的每个数据点分配一个簇标签

1
print("Cluster memberships:\n{}".format(kmeans.labels_))

image-20220414203009442

调用predict方法会为新数据点分配簇标签，返回同上的结果。

1
print(kmeans.predict(X))

可以看出，每个数据点都有标签，但这些标签本身没有任何含义，算法唯一给出的信息就是同一个标签下的数据点非常相似而已。

当然我们也可以指定其他的簇的数量。

k均值聚类算法失败的案例

即使知道给定数据集中簇的正确个数，k均值可能也不是总能找到他们。k均值只能找到相对简单的形状，此外他还假设所有的簇在某种程度上具有相同的“直径”，所以他总是将簇之间的边界刚好画在簇的中心的中间位置。这会产生令人惊讶的结果：

1
2
3
4
5
6
7
X_varied, y_varied = make_blobs(n_samples=200, cluster_std=[1.0,2.5,0.5],random_state=170)
y_pred = KMeans(n_clusters=3, random_state=0).fit_predict(X_varied)
mglearn.discrete_scatter(X_varied[:,0], X_varied[:,1], y_pred)
plt.legend(["cluster 0", "cluster 1", "cluster 2"], loc='best')
plt.xlabel("Feature 0")
plt.ylabel("Feature 1")
plt.show()

image-20220414203028375

可以看出，“应该在”中间零散的数据点组成的簇中的一些数据点，分布在了左下角和右上角密集簇中。

又如：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
X, y = make_blobs(n_samples=600, random_state=170)
rng = np.random.RandomState(74)
# 变换数据使其拉长
transformation = rng.normal(size=(2,2))
X = np.dot(X, transformation)

# 将数据聚类成3簇
kmeans = KMeans(n_clusters=3)
kmeans.fit(X)
y_pred = kmeans.predict(X)

# 画出簇分配和簇中心
plt.scatter(X[:,0], X[:,1], c=y_pred, cmap=mglearn.cm3)
plt.scatter(kmeans.cluster_centers_[:,0], kmeans.cluster_centers_[:,1], marker='^', c=[0,1,2], s=100, linewidths=2, cmap=mglearn.cm3)
plt.xlabel("Feature 0")
plt.ylabel("Feature 1")
plt.show()

image-20220414203038877

可以看出，数据中包含明确分开的三部分，但这三部分被沿着对角线方向拉长，由于k均值仅考虑到最近簇中心的距离，所以他无法处理这种类型的数据。

如果簇的形状更复杂，那么k均值的表现会更差。

矢量量化，或者将k均值看作分解

k均值尝试利用簇中心来表示每个数据点，可以将其看成仅用一个分量来表示每个数据点，该分量由簇中心给出，这种观点将k均值看作是一种分解方法，其中每个点用单一分量来表示，这种观点称为矢量量化。

对于k均值，重建就是在训练集中找到的最近的簇中心。

总结

k均值是非常流行的聚类算法，缺点在于它依赖于随机初始化，此外它对簇形状的假设的约束性较强，而且还要求指定所要寻找的簇的个数（现实应用中可能并不知道这个数字）。

3.5.2 凝聚聚类

简介

层次聚类：是一种很直观的算法。顾名思义就是要一层一层地进行聚类，可以从下而上地把小的cluster合并聚集，也可以从上而下地将大的cluster进行分割。自下而上地进行聚类称为凝聚式层次聚类，自上而下地进行聚类称为分裂式层次聚类。似乎一般用得比较多的是从下而上地聚集。

  • 凝聚式：从点作为个体簇开始，每一步合并两个最接近的簇。
  • 分裂式：从包含所有点的某个簇开始，每一步分裂一个簇，直到仅剩下单点簇。

凝聚聚类（agglomerative clustering）指的是许多基于相同原则构建的聚类算法，这一原则是：算法首先声明每个点是自己的簇，然后合并两个最相似的簇，直到满足某种停止准则为止。scikit-learn中实现的停止准则是簇的个数，因此相似的簇被合并，直到仅剩下指定个数的簇。还有一些链接准则，规定如何度量“最相似的簇”。这种度量总是定义在两个现有的簇之间。

scikit-learn中提供三种选项：

  • ward：默认选项，ward挑选两个簇来合并，使得所有簇中的方差增加最小。这通常会得到大小差不多相等的簇。
  • average：将簇中所有点之间平均距离最小的两个簇合并。
  • complete：也称为最大链接，将簇中点之间最大距离最小的两个簇合并。

ward适用于大多数数据集，如果簇中的成员个数非常不同，那么可以使用后两者效果可能更好。

算法步骤

  1. 如果需要，计算邻近度矩阵

  2. repeat

  3. 合并最接近的两个簇

  4. 更新邻近度矩阵，以反映新的簇与原来的簇之间的邻近性

  5. until 仅剩下一个簇/达到某个终止条件

其中的概念有：

  • 邻近度矩阵

邻近度有许多种定义方式，比如欧氏距离，曼哈顿距离，马氏距离，余弦相似度，Jaccard系数，Bregman散度等等。种类丰富，样品奇多，根据不同的需求来选择最适合的邻近度，计算得到相应的邻近度矩阵。

  • 簇与簇之间邻近度的定义

每个簇中的点数不一定相等，如何计算两个不同簇之间的邻近度呢？

常用的有三种方法：单链（MIN准则），全链（MAX准则），组平均技术。

  • 单链：单链（MIN）定义簇的邻近度为不同簇的两个最近的点之间的邻近度。在图中，即不同结点子集中两个节点之间的最短边。下图中的虚线，就是左右两个簇的邻近度。

image-20220414203056206

  • 全链（MAX）定义簇的邻近度为不同簇中两个最远的点之间的邻近度作为簇的邻近度。在图中，即不同结点子集中两个结点之间的最长边。下图中的虚线，就是左右两个簇的邻近度。

image-20220414203214244

  • 组平均是一种基于图的方法。它定义簇邻近度为取自不同簇的所有点对邻近度的平均值（平均长度）。

image-20220414203227635

下图给出了在一个二维数据及上的凝聚聚类过程，指定要寻找三个簇：

image-20220414203239904

最开始，每个点自成一簇，然后在每一步中，相距最近的两个簇被合并。

凝聚算发不能对新数据点做出预测，因此没有predict方法，但是有fit_predict方法。

1
2
3
4
5
6
7
8
9
10
11
12
from sklearn.cluster import AgglomerativeClustering
from sklearn.datasets import make_blobs
import matplotlib.pyplot as plt
import mglearn
X, y = make_blobs(random_state=1)
agg = AgglomerativeClustering(n_clusters=3)
assignment = agg.fit_predict(X)

mglearn.discrete_scatter(X[:,0], X[:,1], assignment)
plt.xlabel("Feature 0")
plt.ylabel("Feature 1")
plt.show()

image-20220414203254162

可以看出，算法完美地完成了聚类。

树状图可视化

层次化的簇分配图如下：

image-20220414203309912

也可采用树状图显示这一过程：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
from scipy.cluster.hierarchy import dendrogram, ward
X, y = make_blobs(random_state=0, n_samples=12)
linkage_array = ward(X)
dendrogram(linkage_array)

ax = plt.gca()
bounds = ax.get_xbound()
ax.plot(bounds, [7.25,7.25], '--', c='k')
ax.plot(bounds, [4,4], '--', c='k')

ax.text(bounds[1], 7.25, ' two clusters', va='center', fontdict={'size': 15})
ax.text(bounds[1], 4, ' three clusters', va='center', fontdict={'size': 15})
plt.xlabel("Sample index")
plt.ylabel("Cluster distance")
plt.show()

image-20220414203327125

树状图在底部显示数据点，然后以这些点（表示单点簇）作为叶节点绘制一棵树，每合并两个簇就添加一个新的父节点。

主要问题

  1. 缺乏全局目标函数

凝聚层次聚类不能视为全局优化目标函数，因为在每一步合并时仅仅局部地确定哪些簇应当合并（或分裂，对于分裂式）。但是避开了解决困难的组合优化问题，并且这样的方法没有局部极小问题或难选择初始点的问题。

  1. 处理不同大小簇的能力

关于如何处理待合并的簇对的相对大小这个问题，解决的方法有两种：一是加权，就是不同簇中的点具有不同的权值；二是非加权，需要考虑每个簇的点数。

比如组平均技术，上面说的组平均实际上是组平均技术的非加权版，全称为“使用算术平均的非加权的对组方法”（UPGMA），Lance-Williams公式的系数也涉及每个被合并簇的大小。而对于加权版本（WPGMA），它的系数是常数：[公式]=1/2，[公式]=1/2，β=0，γ=0。通常非加权方法更可取，除非出于某种原因个体点具有不同的权值，如对象类非均匀的抽样。

  1. 合并决策是最终的

对于合并两个簇，凝聚层次聚类可以使用所有点的逐对相似度信息趋向于作出好的局部决策。然而，一旦做出合并两个簇的决策，以后就不能撤销，这阻碍了局部最优标准变成全局最优标准。一些方法试图克服“合并是最终的”这一限制，如移动树的分支以改善全局目标；使用划分聚类方法（K均值）来创建许多小簇，然后从这些小簇出发进行层次聚类。

3.5.3 DBSCAN

简介

DBSCAN（Density-Based Spatial Clustering of Applications with Noise，具有噪声的基于密度的聚类方法）是一种基于密度的空间聚类算法。 该算法将具有足够密度的区域划分为簇，并在具有噪声的空间数据库中发现任意形状的簇，它将簇定义为密度相连的点的最大集合。

主要优点在于不需要用户先验地设置簇的个数，可以划分具有复杂形状的簇，还可以找出不属于任何簇的点。

算法过程

DBSCAN算法有两个参数：eps和min_samples。

  • eps：epsilon，为以某个点为圆心的圆的半径，默认为0.5
  • min_samples：在这个圆中至少含有数据点的个数

DBSCAN算法将数据点类型分为三类：

  1. 核心点：在半径Eps内含有超过MinPts数目的点。
  2. 边界点：在半径Eps内点的数量小于MinPts,但是落在核心点的邻域内的点。
  3. 噪音点：既不是核心点也不是边界点的点。

算法首先任意选取一个点，然后找到到这个点的距离小于等于eps的所有的点。如果所有点的个数小于min_samples，那么这个点被标记为噪声，不属于任何一个簇；否则这个点标记为核心点，并被分配一个新的簇标签。然后访问该点的所有邻居（在eps范围内），如果他们还没有簇标签，则将核心点的标签赋给他们。如果他们是核心点，那么就依次访问邻居，以此类推。直到在簇的eps距离内没有更多的核心点为止，然后选取另一个尚未被访问的点，并重复相同的过程。

以下展示了两个参数设置对聚类效果的影响：

image-20220414203359665

image-20220414203408588

上图中，较大点为核心点，较小的为边界点，没有颜色的为噪声。

虽然DBSCAN不需要显式地设置簇的个数，但是设置eps可以隐式地控制找到的簇的个数。

在moon数据集上的表现

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
X, y = make_moons(n_samples=200, noise=0.05, random_state=0)

# 数据预处理
# 将数据缩放成平均值为0，方差为1
scaler = StandardScaler()
scaler.fit(X)
X_scaled = scaler.transform(X)

dbscan = DBSCAN()
clusters = dbscan.fit_predict(X_scaled)
# 绘制簇分配
plt.scatter(X_scaled[:,0], X_scaled[:,1], c=clusters, cmap=mglearn.cm2, s=60)
plt.xlabel("Feature 0")
plt.ylabel("Feature 1")
plt.show()

image-20220414203420906

采用默认eps可以达到理想的效果，如果减小eps，可能会出现更多的簇，以下是去eps=0.2时的情况：

image-20220414203430239

3.5.4 聚类算法的对比与评估

用真实值评估聚类

有一些指标可以用于评估聚类算法相对于真是聚类的结果。其中重要的有：

  • 调整兰德指数（adjusted rand index，ARI）： 取值在[-1,1]之间，负数代表结果不好，越接近于1越好；对任意数量的聚类中心和样本数，随机聚类的ARI都非常接近于0；可用于聚类算法之间的比较。缺点在于ARI需要真实的标签信息对聚类结果进行评价。
  • 归一化互信息（normalized mutual information，NMI）：将互信息放在[0,1]之间，容易评价算法的好坏。

下面是利用ARI来比较三种聚类算法：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
from numpy.core.fromnumeric import size
from numpy.random import seed
from sklearn import cluster
from kmeans_proc import K_means
from sklearn.metrics.cluster import adjusted_rand_score
from sklearn.preprocessing import StandardScaler
from sklearn.datasets import make_moons
from sklearn.cluster import KMeans, AgglomerativeClustering, DBSCAN
import matplotlib.pyplot as plt
import numpy as np
import mglearn
X, y = make_moons(n_samples=200, noise=0.05, random_state=0)

# 将数据缩放成平均值为0，方差为1
scaler = StandardScaler()
scaler.fit(X)
X_scaled = scaler.transform(X)

fix, axes = plt.subplots(1, 4, figsize=(15,3), subplot_kw={'xticks': (), 'yticks': ()})

# 列出要使用的算法
algorithms = [KMeans(n_clusters=2), AgglomerativeClustering(n_clusters=2), DBSCAN()]

# 创建一个随机的簇分配，作为参考
random_state = np.random.RandomState(seed=0)
random_clusters = random_state.randint(low=0, high=2, size=len(X))

# 绘制随机分配
axes[0].scatter(X_scaled[:,0], X_scaled[:,1], c=random_clusters, cmap=mglearn.cm3, s=60)
axes[0].set_title("Random assignment - ARI: {:.2f}".format(adjusted_rand_score(y, random_clusters)))
for ax, algorithm in zip(axes[1:], algorithms):
    # 绘制簇分配和簇中心
    clusters = algorithm.fit_predict(X_scaled)
    ax.scatter(X_scaled[:,0], X_scaled[:,1], c=clusters,cmap=mglearn.cm3, s=60)
    ax.set_title("{} - ARI: {:.2f}".format(algorithm.__class__.__name__,adjusted_rand_score(y, clusters)))
plt.show()

image-20220414203457866

从上图可以看出，调整兰德指数给出了符合直觉的结果。

注意：用这种真实值方式评估聚类时，一个常见的错误就是使用accuracy_score而不是ARI或NMI或其他聚类指标来评估聚类算法，使用精度的问题在于，他要求分配的簇标签与真实值完全匹配，但簇标签本身毫无意义——唯一重要的是哪些点位于同一个簇中。

1
2
3
4
5
6
7
8
from sklearn.metrics import accuracy_score
# 这两种点标签对应于相同的聚类
clusters1 = [0,0,1,1,0]
clusters2 = [1,1,0,0,1]
# 精度为0，因为两者标签完全不同
print("Accuracy: {:.2f}".format(accuracy_score(clusters1, clusters2)))
# 调整兰德指数为1，因为两者的聚类完全相同
print("ARI: {:.2f}".format(adjusted_rand_score(clusters1, clusters2)))

image-20220414203507543

在没有真实值的情况下评估聚类

  • 轮廓系数：轮廓分数计算一个簇的紧致度，其值越大越好，最高为1，但是对于复杂形状效果不好。

算法：假设我们已经通过一定算法，将待分类数据进行了聚类。常用的比如使用K-means ，将待分类数据分为了 k 个簇 。对于簇中的每个向量。分别计算它们的轮廓系数。

对于其中的一个点 i 来说：

计算 a(i) = average(i向量到所有它属于的簇中其它点的距离)

计算 b(i) = min (i向量到与它相邻最近的一簇内的所有点的平均距离)

那么 i 向量轮廓系数就为：

image-20220414203519813

可见轮廓系数的值是介于 [-1,1] ，越趋近于1代表内聚度和分离度都相对较优。将所有点的轮廓系数求平均，就是该聚类结果总的轮廓系数。a(i) ：i向量到同一簇内其他点不相似程度的平均值；b(i) ：i向量到其他簇的平均不相似程度的最小值。

以下是利用轮廓系数对三种算法的评估：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
from sklearn.metrics.cluster import silhouette_score
X, y = make_moons(n_samples=200, noise=0.05, random_state=0)

# 将数据缩放成平均值为0，方差为1
scaler = StandardScaler()
scaler.fit(X)
X_scaled = scaler.transform(X)

fix, axes = plt.subplots(1, 4, figsize=(15,3), subplot_kw={'xticks': (), 'yticks': ()})

# 列出要使用的算法
algorithms = [KMeans(n_clusters=2), AgglomerativeClustering(n_clusters=2), DBSCAN()]

# 创建一个随机的簇分配，作为参考
random_state = np.random.RandomState(seed=0)
random_clusters = random_state.randint(low=0, high=2, size=len(X))

# 绘制随机分配
axes[0].scatter(X_scaled[:,0], X_scaled[:,1], c=random_clusters, cmap=mglearn.cm3, s=60)
axes[0].set_title("Random assignment: {:.2f}".format(silhouette_score(X_scaled, random_clusters)))
for ax, algorithm in zip(axes[1:], algorithms):
    # 绘制簇分配和簇中心
    clusters = algorithm.fit_predict(X_scaled)
    ax.scatter(X_scaled[:,0], X_scaled[:,1], c=clusters,cmap=mglearn.cm3, s=60)
    ax.set_title("{} : {:.2f}".format(algorithm.__class__.__name__,silhouette_score(X_scaled, clusters)))
plt.show()

image-20220414203547178

可以看出，k均值分数最高，但实际聚类效果不如分数较低的DBSCAN算法。

3.5.5 聚类方法小结

聚类的应用于评估是一个非常定性的过程，通常在数据分析的探索阶段很有帮助。

k均值可以用簇的平均值来表示簇，他还可以看作一种分解方法，每个数据点都有其簇中心表示。

DBSCAN可以检测到没有分配任何簇的噪声点，还可以帮助自动判断簇的数量，它允许簇有复杂的形状。

凝聚聚类可以提供数据的可能划分的整个层次结构，可以通过树状图轻松查看。

第四章 数据表示与特征工程

  • 连续特征：数据由浮点数组成。

  • 离散特征（或分类特征）：这种特征通常不是数值，不以连续的方式变化，例如产品的品牌，颜色，部门等。

数据表示方式都会对机器学习模型的性能产生巨大影响。对于某个特定应用，如何找到最佳数据表示，这个问题被称为特征工程。

4.1 分类变量

adult数据集：任务是预测一名工人的收入是高于5w，还是低于5w，属于分类任务，在这个数据集中age，hours-per-week属于连续特征，但是workclass，education，gender等属于离散特征，他们都来自于一系列固定的可能取值（而不是一个范围），表示的是定性属性（而不是数量）。

我们需要换一种方式来表示数据。

4.1.1 One-Hot编码（虚拟变量）

表示分类变量最常用的方法就是使用one-hot编码（或N取1编码），也叫虚拟变量。虚拟变量背后的思想是将一个分类变量替换为一个或多个新特征，新特征取值为0或1。

将数据转换为分类变量的one-hot编码常常使用pandas。

1
2
3
4
5
6
7
8
9
10
11
import pandas as pd
from IPython.display import display

data = pd.read_csv("data/adult.data", header=None, index_col=False, names=[
    'age','workclass','fnlwgt','education','education-num','martial-status','occupation',
    'relationship','race','gender','capital-gain','capital-loss','hours-per-week','native-country','income'
])
# 
data = data[['age','workclass','education','gender','hours-per-week','occupation','income']]

display(data.head())

image-20220414203614356

用pandas编码数据可以使用get_dummies函数，它将自动变换所有具有对象类型（例如字符串）的列或所有分类的列。

1
2
3
print("Original features:\n", list(data.columns), "\n")
data_dummies = pd.get_dummies(data)
print("Features after get_dummies:\n", list(data_dummies.columns))

image-20220414203630326

可以看出，连续特征age和hours-per-week没有发生变化，而分类特征的每个可能取值都被扩展为一个新特征：

image-20220414203640001

训练模型

下面使用values属性将data_dummies数据框转换为NumPy数组，然后在其上训练一个机器学习模型。在训练模型之前，要把两个目标变量（现在被编码为两个income列）从数据中分离出来。

常见错误：将输出变量或输出变量的一些导出属性包含在特征表示中，是构建监督学习模型的一个常见错误。

1
2
3
4
5
features = data_dummies.loc[:,'age':'occupation_ Transport-moving']
# 提取NumPy数组
X = features.values
y = data_dummies['income_ >50K'].values
print("X.shape: {} y.shape: {}".format(X.shape, y.shape))

image-20220414203648795

现在的数据表示就可以被scikit-learn处理：

1
2
3
4
5
6
from sklearn.linear_model import LogisticRegression
from sklearn.model_selection import train_test_split
X_train ,X_test, y_train, y_test = train_test_split(X, y, random_state=0)
logreg = LogisticRegression()
logreg.fit(X_train, y_train)
print("Test score: {:.2f}".format(logreg.score(X_test, y_test)))

image-20220414203656231

4.1.2 数字可以编码分类变量

分类特征通常用整数进行编码，它们是数字并不意味着它们必须被视为连续特征。一个整数特征被视为连续的还是离散的，有时并不明确。

pandas的get_dummies函数将所有数字看作是连续特征，不会为其创建虚拟变量。可以将数据框中的数值列转换成字符串，再进行处理：

1
2
demo_df = pd.DataFrame({'Integer Feature': [0,1,2,1], 'Categorical Feature': ['socks','fox','socks','bix']})
display(demo_df)

image-20220414203703932

以上是包含分类字符串和整数特征的数据框。

使用get_dummies函数只会编码字符串特征，不会改变整数特征：

1
display(pd.get_dummies(demo_df))

image-20220414203712867

如果想为Interger Feature这一列创建虚拟变量，可以使用columns参数显式地给出想要编码的列，于是两个特征都会被当做分类处理特征处理：

1
2
3
demo_df['Integer Feature'] = demo_df['Integer Feature'].astype(str)
demo_hot = pd.get_dummies(demo_df, columns=['Integer Feature', 'Categorical Feature'])
display(demo_hot)

image-20220414203721756

4.2 分箱，离散化，线性模型与树

数据表示的最佳方法不仅取决于数据的语义，还取决于所使用的的模型种类。以下是在wave数据集上比较线性回归和决策树的图：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
import numpy as np
from sklearn.linear_model import LinearRegression
from sklearn.tree import DecisionTreeRegressor
from mglearn.datasets import make_wave
import matplotlib.pyplot as plt
import numpy as np
X ,y = make_wave(n_samples=100)
line = np.linspace(-3, 3, 1000, endpoint=False).reshape(-1, 1)

reg = DecisionTreeRegressor(min_samples_split=3).fit(X, y)
plt.plot(line, reg.predict(line), label="decision tree")

reg = LinearRegression().fit(X, y)
plt.plot(line, reg.predict(line), label="linear regression")

plt.plot(X[:,0], y, 'o', c='k')
plt.ylabel("Regression output")
plt.xlabel("Input feature")
plt.legend(loc="best")
plt.show()

image-20220414203732221

有一种方法可以让线性模型在连续数据上变得更加强大，就是利用特征分箱（或者离散化）将其划分为多个特征。

将特征的输入范围划分成固定个数的箱子（bin），比如10个，那么数据点就可以利用它们所在的箱子来表示。

1
2
3
# 在-3,3的区间上划分十个段（箱子）
bins = np.linspace(-3, 3, 11)
print("bins: {}".format(bins))

image-20220414203744087

接下来记录每个数据点所属的箱子。

1
2
3
which_bin = np.digitize(X, bins=bins)
print("\nData points:\n", X[:5])
print("\nBin membership for data points:\n", which_bin[:5])

image-20220414203752427

这里做的就是将wave数据集中单个连续输入特征变换为一个分类特征（离散特征）。

下面将这个离散特征用one-hot编码表示：

1
2
3
4
5
6
7
8
from sklearn.preprocessing import OneHotEncoder
# 使用OneHotEncoder进行变换
encoder = OneHotEncoder(sparse=False)
# encoder.fit找到which_bin中的唯一值
encoder.fit(which_bin)
# transform创建one-hot编码
X_binned = encoder.transform(which_bin)
print(X_binned[:5])

image-20220414203800634

由于我们指定了10个箱子，因此变换后的X_binned数据集现在包含10个特征。

现在在X_binned数据集上进行构建模型：

1
2
3
4
5
6
7
8
9
10
11
12
line_binned = encoder.transform(np.digitize(line, bins=bins))

reg = DecisionTreeRegressor(min_samples_split=3).fit(X_binned, y)
plt.plot(line, reg.predict(line_binned), label="decision tree")
reg = LinearRegression().fit(X_binned, y)
plt.plot(line, reg.predict(line_binned), label="linear regression")

plt.plot(X[:,0], y, 'o', c='k')
plt.ylabel("Regression output")
plt.xlabel("Input feature")
plt.legend(loc="best")
plt.show()

image-20220414203809875

线性模型和决策树的线段完全重合，说明两者做出了相同的预测。对于每个箱子，二者都预测了同一个常数值。实际上每个线段的拐点就是两个箱子的边界。可以看出分箱，让线性模型的表现力得到了一定的提高。

适用范围

对于特定的数据集，如果有理由选择使用线性模型，例如数据集很大，维度很高，但有些特征与输出的关系是非线性的，那么分箱是提高建模能力的好方法。

4.3 交互特征与多项式特征

想要丰富特征表示，特别是对于线性模型而言，另一种方式是添加原始特征的交互特征和多项式特征。

交互特征

在上面的例子中，加入原始特征（即x轴），这样会得到11维的数据集：

1
2
3
X_binned = encoder.transform(which_bin)
X_combined = np.hstack([X, X_binned])
print(X_combined.shape)

image-20220414203819830

现在利用X_combined数据集构建模型：

1
2
3
4
5
6
7
8
9
10
reg = LinearRegression().fit(X_combined, y)
line_combined = np.hstack([line, line_binned])
plt.plot(line, reg.predict(line_combined), label="linear regression combined")
for bin in bins:
    plt.plot([bin,bin], [-3,3], ':', c='k')
plt.ylabel("Regression output")
plt.xlabel("Input feature")
plt.legend(loc="best")
plt.plot(X[:,0], y, 'o', c='k')
plt.show()

image-20220414203828991

可以看出，模型学习到了一个向下的斜率，每个箱子都是相同的，现在我们要让每个箱子学习到不同的斜率，我们可以添加一个交互特征，这个特征是箱子指示符与原始特征的乘积。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
# 将X_binned数据集变换，添加交互特征
X_product = np.hstack([X_binned, X * X_binned])
reg = LinearRegression().fit(X_product, y)

line_product = np.hstack([line_binned, line * line_binned])

plt.plot(line, reg.predict(line_product), label="linear regression product")

for bin in bins:
    plt.plot([bin,bin], [-3,3], ':', c='k')
plt.ylabel("Regression output")
plt.xlabel("Input feature")
plt.legend(loc="best")
plt.plot(X[:,0], y, 'o', c='k')
plt.show()

image-20220414203838724

现在模型中每个箱子都有自己的偏移和斜率。

多项式特征

另一种方法是使用原始特征的多项式特征。对于给定特征x，可以考虑x乘2，乘3，平方等等。

1
2
3
4
5
6
from sklearn.preprocessing import PolynomialFeatures
# 包含直到x**10的多项式
poly = PolynomialFeatures(degree=10, include_bias=False)
poly.fit(X)
X_poly = poly.transform(X)
print("X_poly.shape: {}".format(X_poly.shape))

image-20220414203846840

此时，经过变换后的X_poly数据集的特征为x，一直乘到10，共有10个特征。

将多项式特征与线性回归模型一起使用，可以得到经典的多项式回归模型

1
2
3
4
5
6
7
8
9
reg = LinearRegression().fit(X_poly, y)
line_poly = poly.transform(line)

plt.plot(line, reg.predict(line_poly), label="polynomial linear regression")
plt.plot(X[:,0], y, 'o', c='k')
plt.ylabel("Regression output")
plt.xlabel("Input feature")
plt.legend(loc="best")
plt.show()

image-20220414203859194

多项式特征在这个一维数据上得到了非常平滑的拟合，但高次多项式在边界上或数据很少的区域上可能会有极端的表现。

作为对比，下面是在原始数据上学到的核SVM模型，没有做任何变换：

1
2
3
4
5
6
7
8
9
from sklearn.svm import SVR
for gamma in [1,10]:
    svr = SVR(gamma=gamma).fit(X, y)
    plt.plot(line, svr.predict(line), label='SVR gamma={}'.format(gamma))
plt.plot(X[:,0], y, 'o', c='k')
plt.ylabel("Regression output")
plt.xlabel("Input feature")
plt.legend(loc="best")
plt.show()

image-20220414203908378

多项式特征的构造方式和作用

先使用boston数据集进行说明：首先加载数据，然后利用MinMaxScaler将其缩放到0到1之间：

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
from numpy.core.fromnumeric import transpose
from sklearn.datasets import load_boston
from sklearn.model_selection import train_test_split
from sklearn.preprocessing import MinMaxScaler
from sklearn.preprocessing import PolynomialFeatures
boston = load_boston()
X_train, X_test, y_train, y_test = train_test_split(boston.data, boston.target, random_state=0)

# 缩放数据
scaler = MinMaxScaler()
X_train_scaled = scaler.fit_transform(X_train)
X_test_scaled = scaler.fit_transform(X_test)

# 下面提取交互特征和多项式特征，次数最高为2
poly = PolynomialFeatures(degree=2).fit(X_train_scaled)
X_trian_poly = poly.transform(X_train_scaled)
X_test_poly = poly.transform(X_test_scaled)
print("X_train.shape: {}".format(X_train.shape))
print("X_train_poly.shape: {}".format(X_trian_poly.shape))
print("Polynomial feature names:\n{}".format(poly.get_feature_names()))

image-20220414203917243

原始数据有13个特征，现在被扩展到105个特征，具体如上面所示。第一个特征为常数特征，为‘1’，接下来13个特征为原始特征，然后是为第一个特征的平方，以及它与其他特征的组合。

下面利用岭回归Ridge在有交互特征的数据上和没有交互特征的数据上进行性能对比：

1
2
3
4
5
6
from sklearn.linear_model import Ridge

ridge = Ridge().fit(X_train_scaled, y_train)
print("Score without interactions: {:.3f}".format(ridge.score(X_test_scaled, y_test)))
ridge = Ridge().fit(X_trian_poly, y_train)
print("Score without interactions: {:.3f}".format(ridge.score(X_test_poly, y_test)))

image-20220414203925771

显然，在使用岭回归时，交互特征和多项式特征对性能提升有很大的帮助，但如果使用更加复杂的模型（比如随机森林），情况会稍有不同：

1
2
3
4
5
from sklearn.ensemble import RandomForestRegressor
rf = RandomForestRegressor(n_estimators=100).fit(X_train_scaled, y_train)
print("Score without interactions: {:.3f}".format(rf.score(X_test_scaled, y_test)))
rf = RandomForestRegressor(n_estimators=100).fit(X_trian_poly, y_train)
print("Score without interactions: {:.3f}".format(rf.score(X_test_poly, y_test)))

image-20220414203933263

没有额外的特征，随机森林的性能也要优于岭回归，但是可以看出，添加了交互特征和多项式特征后，实际上会略微降低随机森林的性能。

4.4 单变量非线性变换

其它变换通常对某些特征也有作用，特别是应用数学函数，例如log，exp或sin。大部分模型都在每个特征大致遵循高斯分布时表现最好，下面的例子：

1
2
3
4
5
6
7
8
9
10
11
12
rnd = np.random.RandomState(0)
X_org = rnd.normal(size=(1000,3))
w = rnd.normal(size=3)

X = rnd.poisson(10 * np.exp(X_org))
y = np.dot(X_org, w)

bins = np.bincount(X[:,0])
plt.bar(range(len(bins)), bins, color='r')
plt.ylabel("Number of appearnces")
plt.xlabel("Value")
plt.show()

image-20220414203941259

利用这个数据集拟合一个岭回归模型：

1
2
3
X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=0)
score = Ridge().fit(X_train, y_train).score(X_test, y_test)
print("Test score: {:.3f}".format(score))

image-20220414203949356

现在对特征应用log变换：变换之后，数据分布的不对称性变小，也没有非常大的异常值。

1
2
3
4
5
6
X_train_log = np.log(X_train + 1)
X_test_log = np.log(X_test + 1)
plt.hist(X_train_log[:,0], bins=25, color='gray')
plt.ylabel("Number of appearnces")
plt.xlabel("Value")
plt.show()

image-20220414203957711

然后拟合一个岭回归模型：

1
2
score = Ridge().fit(X_train_log, y_train).score(X_test_log, y_test)
print("Test score: {:.3f}".format(score))

在所举的例子中，所有特征都具有相同的性质。通常来说，只有一部分特征应该进行变换，有时每个特征的变换方式也各不相同。

4.5 自动化特征选择

在添加新特征或处理一般的高维数据集时，最好将特征的数量减少到只包含最有用的那些特征，并且删除其余特征，这样会得到泛化能力更好，更简单的模型。如何判断每个特征的作用，方法主要有三种：

  • 单变量统计（univariate statistics）
  • 基于模型的选择（model-based selection）
  • 迭代选择（iterative selection）

4.5.1 单变量统计

在单变量统计中，我们计算每个特征和目标值之间的关系是否存在统计显著性，然后选择具有最高置信度的特征。这个方法只单独考虑每个特征，如果一个特征只有在与另一个特征合并时才具有信息量，那么这个特征将被舍弃。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
rom sklearn.datasets import load_breast_cancer
from sklearn.feature_selection import SelectPercentile
from sklearn.model_selection import train_test_split
import matplotlib.pyplot as plt
import numpy as np
cancer = load_breast_cancer()

# 获得确定性的随机数
rng = np.random.RandomState(42)
noise = rng.normal(size=(len(cancer.data), 50))

# 添加噪声特征
# 前30个来自数据集，后50个是噪声
X_w_noise = np.hstack([cancer.data, noise])

X_train, X_test, y_train, y_test = train_test_split(X_w_noise, cancer.target, random_state=0, test_size=.5)

# 使用f_classif(默认值)和SelectPercentile来选择50%的特征
select = SelectPercentile(percentile=50)
select.fit(X_train, y_train)

# 对训练集进行变换
X_train_selected = select.transform(X_train)

print("X_train.shape: {}".format(X_train.shape))
print("X_train_selected.shape: {}".format(X_train_selected.shape))

image-20220414204016423

特征的数量从80减少到40个。布尔遮罩如下：

1
2
3
4
5
6
mask = select.get_support()
print(mask)
# 将遮罩可视化，黑色为被选中，白色为舍弃
plt.matshow(mask.reshape(1,-1), cmap='gray_r')
plt.xlabel("Sample index")
plt.show()

image-20220414204027348

可以看出大多数所选择的特征都是原始特征，并且大多数噪声特征都已被删除。下面是Logistic回归的拟合情况：

1
2
3
4
5
6
7
8
9
10
from sklearn.linear_model import LogisticRegression

# 对测试集进行变换
X_test_selected = select.transform(X_test)

lr = LogisticRegression()
lr.fit(X_train, y_train)
print("Score with all features: {:.3f}".format(lr.score(X_test, y_test)))
lr.fit(X_train_selected, y_train)
print("Score with only selected features: {:.3f}".format(lr.score(X_train_selected, y_test)))

image-20220414204049970

4.5.2 基于模型的特征选择

基于模型的特征选择使用一个监督机器学习模型来判断每个特征的重要性，并且仅保留最重要的特征。用于特征选择的模型不需要与用于最终监督建模的模型相同。与单变量选择不同，基于模型的选择同时考虑所有特征。

下面的例子使用包含100棵树的随机森林分类器来计算特征重要性。

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
from sklearn.feature_selection import SelectFromModel
from sklearn.ensemble import RandomForestClassifier
from sklearn.datasets import load_breast_cancer
from sklearn.model_selection import train_test_split
import matplotlib.pyplot as plt
import numpy as np
cancer = load_breast_cancer()
# 获得确定性的随机数
rng = np.random.RandomState(42)
noise = rng.normal(size=(len(cancer.data), 50))
# 添加噪声特征
# 前30个来自数据集，后50个是噪声
X_w_noise = np.hstack([cancer.data, noise])
X_train, X_test, y_train, y_test = train_test_split(X_w_noise, cancer.target, random_state=0, test_size=.5)
# 采用中位数作为阈值，这样可以得到一半的特征
select = SelectFromModel(RandomForestClassifier(n_estimators=100, random_state=42),threshold="median")
select.fit(X_train, y_train)
X_train_l1 = select.transform(X_train)
print("X_train.shape: {}".format(X_train.shape))
print("X_train_l1.shape: {}".format(X_train_l1.shape))
mask = select.get_support()
# 将遮罩可视化，黑色为被选中，白色为舍弃
plt.matshow(mask.reshape(1,-1), cmap='gray_r')
plt.xlabel("Sample index")
plt.show()

image-20220414204108072

image-20220414204118389

1
2
3
X_test_l1 = select.transform(X_test)
score = LogisticRegression(max_iter=5000).fit(X_train_l1, y_train).score(X_test_l1, y_test)
print("Test score: {:.3f}".format(score))

image-20220414204126541

4.5.3 迭代特征选择

在迭代特征选择中，将会构建一系列模型，每个模型都使用不同数量的特征。有两种基本方法：

  • 开始时没有特征，然后逐个添加特征，直到满足某个终止条件。
  • 从所有特征开始，然后逐个删除特征，直到满足某个终止条件。

其中一种特殊方法是递归特征消除（RFE），它从所有特征开始构建模型，并根据模型舍弃最不重要的特征，然后使用除被舍弃特征之外的所有特征来构建一个新的模型，如此继续，直到剩下预设数量的特征。

1
2
3
4
5
6
7
from sklearn.feature_selection import RFE
select = RFE(RandomForestClassifier(n_estimators=100, random_state=42),n_features_to_select=40)
select.fit(X_train, y_train)
mask = select.get_support()
plt.matshow(mask.reshape(1,-1), cmap='gray_r')
plt.xlabel("Sample index")
plt.show()

image-20220414204137686

对一个随机森林训练了40次，每训练每一次就删除一个特征。

1
2
3
4
5
X_train_rfe = select.transform(X_train)
X_test_rfe = select.transform(X_test)

score = LogisticRegression(max_iter=1000).fit(X_train_rfe, y_train).score(X_test_rfe, y_test)
print("Test score: {:.3f}".format(score))

image-20220414204145980

4.6 小结

对于线性模型，可能会从分箱，添加多项式和交互项而生成的新特征中大大受益，对于更加复杂的非线性模型（例如随机森林和SVM），在无需显示扩展特征空间的前提下就可以学习更加复杂的任务。在实践中，所使用的特征（以及特征与方法之间的匹配）通常是使机器学习方法表现良好的重要因素。