2.2　训练一个传统的机器学习模型

运用机器学习方法分析数据、建立预测模型本身是一个非常复杂的过程，很难在较短的篇幅内完全说清楚，所以这里我们将会结合一个简单的实战案例，将基本概念结合代码实现出来。同时，本节也将穿插介绍如何使用Python的数据科学套装，如表2-2所示。

表2-2　使用Python数据的科学套装

如果使用了第1章中指定的开发环境，这里直接使用就好，否则需要注意一下，版本号的偏差可能会造成程序报告警告，说某一用法是以前的，未来版本会抛弃，也有可能会直接报错。如果本章接下来的程序提示错误，请大家首先确认使用的是否为指定的开发环境以及正确的版本号。

2.2.1　第一步，观察数据

我们这里使用sklearn官方提供的鸢尾花分类数据集。这个数据集最初是埃德加·安德森从加拿大加斯帕半岛上的鸢尾属花朵中提取的地理变异数据，包含150个样本，属于鸢尾属下的三个亚属，即山鸢尾（setosa）、变色鸢尾（versicolor）和维吉尼亚鸢尾（virginica）。四个特征被用作样本的定量分析，分别是花萼（Sepal）和花瓣（Petal）的长度和宽度。这个案例的目的是通过建立一种数学模型，尝试使用四个特征去预测某一鸢尾花属于哪一个亚种。具体如何建立这种模型，接下来我们将会逐步讲解。

要建立模型，首先需要观察数据，而观察数据的第一步是要阅读数据相关的说明文档，弄清楚我们拿到的数据是哪一种具体格式。虽然我们知道了数据有多少个样本、多少种特征，但是目前还不清楚数据是以什么格式给我们的，是文件还是数据库，或者是一个python对象，所以需要去sklearn官网查阅说明文档，地址是http://scikit-learn.org/stable/modules/generated/sklearn.datasets.load_iris.html。

文档提到返回的对象如下：

文档信息写到，这是一个Python对象（object），并且具有字典（dictionary）的特征。具体的调用方法如下：

    from sklearn.datasets import load_iris
    data = load_iris()

2.2.2　第二步，预览数据

既然上一步提到这是一个具有字典特征的对象，我们就可以用Python字典调用的方法对数据有一个总体的预览。

这里首先简单介绍Python的字典。Python主要包括列表、元组以及字典这几种较为高级的数据结构，如果读者熟悉C++，其实就是C++标准库里的vector、set以及unordered_map。字典结构其实就是一个键值对（unordered_map），但这里的键值对相比C++而言，用起来相对容易一些，可以方便地指向字符串，甚至列表、数组和其他字典：

运行结果：

          #out:
          0
          [1, 2, 3] 2
          4

可以使用for循环遍历整个字典：

    for k in map_abc:
          print("key:",k, "\nvalue:", map_abc[k])

运行结果：

          #out:
          key: b
          value: [1, 2, 3]
          key: a
          value: 0
          key: c
          value: {'cc': 4}

同样，也来遍历我们的数据：

data.data就是需要的150×4的输入特征矩阵，四个特征存在data.feature_names里。分类标签在data.target中，其中的0、1、2分别代表data.target_names里面的setosa、versicolor以及virginica。由此，我们完成了数据的基本预览。

实际上，在预览阶段还有很多问题需要思考：

• 不同分类的样本分类是否均匀？
• 数据是否受到极值、缺失值的影响？
• 能否不建立模型就看出三者的分类关系？
• 能否对样本的分类情况进行简单的预览？

这些问题可以用很简单的程序回答。其实我们想更快地回答这些问题，因此接下来将引入Python数据分析套装，用pandas对数据进行类似SQL的合理结构化，然后基于可视化分析库matplotlib与seaborn，通过图形回答这些问题。

对矩阵格式的数据进行结构化处理，转换为pandas的dataframe。

    df = pd.DataFrame(data.data)
    df.columns = data.feature_names
    df['species'] = [ data['target_names'][x] for x in data.target ]
    df.head()

不同分类的样本分类是否均匀？

对species分组计数：

    df_cnt = df['species'].value_counts().reset_index()
    df_cnt

结果如下：（注意这里不加reset_index的话，将不会返回数据框的形式，不方便制图。）

对结果做图：

    sns.barplot(data=df_cnt, x='index', y='species')

其结果如图2-2所示。

图2-2　分组计数结果做图效果

数据是否受到极值、缺失值的影响？

极值、缺失值会影响对数据的整体认识，对模型产生干扰，因此，我们在正式分析数据之前，首先应该关注它们会产生什么样的影响。这里的极值是指相比其他数字而言，非常大或者非常小的值，又称为离群值，既可能是由于收集阶段的错误造成的，也可能是由于一些意外因素造成的。极值会对数据分析结果造成一定干扰，典型的例子如同国家统计局发布人均年收入数据时，很多人发现自己“被平均”，因此对于极值，会根据实际需求选择是否剔除。

缺失值同样也会有影响。比如小明给自己量体温，量了一周的体温值：

问小明这一周的平均体温是多少？这种情况下，真实的平均体温我们是不知道的，因为不知道小明周四的体温是多少，所以整周的平均值也是不知道的：

    import numpy as np
    np_temp = [36.2, 36.3, 36.4, np.nan, 36.3, 36.2, 36.4]
    np.mean(np_temp)

运行结果：

          # out:
            nan

遇到这种情况，实际上是对数据的一种浪费，因为周四没有收集到体温数据，这一周其他6天的数据就失去了作用。在大规模的数据中，我们很难保证所有数据都被合理地收集到，此时如果有几万个数据，因为个别的缺失就全部扔掉，这种浪费就更加明显了。对于这种情况，通常的做法是在计算的情况下不考虑未知量，或者是用平均值、中位数去填补未知值。

方法1　不考虑周四

    np_temp = [36.2, 36.3, 36.4, np.nan, 36.3, 36.2, 36.4]
    np.nanmean(np_temp)

运行结果：

          # out
          36.3...

方法2　用其他数字平均值填补缺失值后计算

    np_temp[3] = np.nanmean(np_temp)
    np.mean(np_temp)

运行结果：

    # out
    36.3...

方法3　用数字0填补缺失值后计算（大错特错，不要这样）

    np_temp[3] = 0
    np.mean(np_temp)

之所以把这一条写上来，是因为现实中真的有人会用0直接填补缺失值。对这种情况，不是说不行，很多时候其实可以用0填补，比如统计某人口稀少地区的各种汽车销量，如果某一汽车销量缺失，我们是可以用0填补的，因为这个数字可能真的很小，可能真的是0，以至于统计人员收集数据时直接忽略了这种汽车。但小明周三的体温，实在不太像是0度，因此这里不可以用0来填补缺失值。

介绍完极值与缺失值将会带来的问题后，我们想知道如何用最快的速度判断数据里面是否存在这两种情况，最简单的办法就是对数据框执行.describe()操作：

    df.describe()

运行结果：

      # out:

这里的极值（min/max）看起来并不明显，平均值、标准差范围也比较接近，并且不存在数据缺失。如果有缺失，最后一行会单独显示缺失了几个值，继而我们进一步确认这四个特征是否为正态分布。这种检验通常使用QQ-Plot，即横坐标是标准正态分布的Quantile（如图2-3所示的x轴），纵坐标是样本实际分布的Quantile，如果实际分布是正态分布，两者之间就会存在线性的关系。

图2-3　横坐标是标准正态分布的Quantile

对我们的四个特征进行正态分布检验：

运行结果如图2-4所示。

这里的前两个特征可以认为来自同一个正态分布，如此很多统计假设（如t检验）等就可以成立，继而在后续的建模阶段方便我们使用很多有效的算法。后两个特征更像是两条线中间断开了，似乎是来自两个正态分布。这种情况未必是坏事，如果某个分类种类完全来自其中一个正态分布，其他的来自另一个正态分布，那么这个特征就成了一个非常具有区分度的特征，需要我们重点关注。

图2-4　四个特征的正态分布

因此下一步的目标就是求具体每个分类的平均值、方差，如何做呢？这里当然可以提取不同分类对应的子矩阵，然后计算子矩阵的平均值、方差，读者可以自行尝试。我们这里提供另一种基于透视表（pivot_table）的方法。首先，对于150×4的二维矩阵特征，将其变为600×1的一维矩阵特征：

    pd.melt(df, id_vars=['species'])

运行结果：

          # out:

然后将这个一维矩阵特征通过透视表操作，针对每个特征值，以分类结果为行名称，求它们的平均值以及方差：

运行结果：

      # out:

我们关注petal length和petal width两个特征，发现二者在setosa这个种类中要小于其他两个种类，并且setosa与versicolor之间的距离（4.260-1.464）也远高于三倍方差值（0.0301×3+0.2208×3），正好对应在图2-4中的两个间隔较远的正态分布。

继而我们以第三个特征为例，对其不同组分别进行正态分布检验：

其结果如图2-5所示。

图2-5　第三个特征的正态分布

由此可见，该特征在不同组内部是符合正态分布的。

能不能不建立模型，就看出三者的分类关系？

首先，这里说一下“不建立模型”的逻辑是什么。不建立模型就看出分类关系，是因为有时候分类的定义很简单，用复杂模型多少有点“杀鸡用牛刀”的意味。比如判断一个人是否酒驾，就一个指标，每百毫升血液里酒精含量是否大于20mg，大于就是酒驾，否则不是，就是一个简单的if…else判断关系，不需要复杂的模型。实际工作中，如果找到了这种简单的if…else标准，其实是件好事，首先我们不需要费时费力地挖特征、训练模型，其次得出的结论也很简单，不涉及复杂的组合，也利于被人接受。

我们的鸢尾花数据集是否存在这种简单标准？可以借助简单的可视化技术来看一眼：

    sns.pairplot(df, hue="species")

运行结果如图2-6所示。

图2-6　鸢尾花数据集的可视化

这个图包含了两层信息：第一层是单一值是否足以分开不同的分类结果，就是对角线上的四个分布图；第二层是数据间的两两组合是否足以分开不同的分类结果，就是非对角线上的散点图，其中左下角与右上角部分x - y坐标轴是对称的。

对于第一层信息，我们注意到第三行、第三列的点图中setosa与其他两类在petal length这个特征上是可以完全区分开的，即petal length小于2的是setosa、大于2的是其他两种分类。另外两个种类区分得并不好，有很多的重叠部分。

对于第二层信息，根据点图判断的话，所有第三行、第三列的点图中，setosa与其他两个种类也是可以完全区分开的。因为这些图是特征的两两组合，这些组合都使用第一层信息就可以明确区分setosa的特征——petal length，自然在两两组合中也可以区分。通过两两组合，另外两种也并非完全分开，但是重叠部分的交集有多有少——sepal的长宽重合得多，petal的长宽重合就少得多。

这时候读者会有新的问题了——既然两个两个组合，看起来后两类重叠的面积少了很多，那么特征如果进行三个组合，画三维立体图重叠的是否会更少？这里不太好画三维的点图，读者可以画四个三维的图看一眼。这里之所以不画并非是笔者懒，而是既然可以三个特征组合，就可以四个一起组合，这种情况下就没法画图了，毕竟我们生活在一个三维的宏观世界，要表示一个四维的坐标系难度其实挺大的。

因此，我们这里的结论是，如果想找出setosa，那么不用建模也可以得出结论，花瓣（petal）长度小于2cm的就是。如果想找出其他两种，结论就不是肯定的了。

当然，对于这种直观方法的局限性，读者应该也发现了，就是特征两两组合会造成维度灾难。我们这里是四个特征，两两组合就会得到4×3÷2 = 6种不同的组合方式。如果是几百个特征，就是上万种组合方式，几万特征就是上亿种组合方式。这种情况下，要是再画这种组合图，人工找出特征组合就不现实了。我们需要计算机帮助组合特征，而计算机组合特征的方法就是借助机器学习模型。

对于各种机器学习模型，本书的篇幅主要是介绍监督学习方法，就是之前提到的应试教育。在此之前，我们先简单介绍下非监督学习找组合特征的方法，并且基于非监督学习的特征组合，对数据进行一个分类结果的预览。

能否对样本的分类情况进行简单的预览？

上文提到，直接将特征以两两组合的方式画在二维平面上，在特征过多时，会由于组合过多使得这种方法难以发挥作用。这时能否直接在二维平面上看到四维的分布情况呢？可以通过降维做到。

提到“降维”，最经典的一段文学描述莫过于《三体》这部科幻小说，在最后阶段描写的外星人使用“二向箔”，对太阳系进行的一场降维打击：

在二维化的太空艇上，可以看到二维展开后的三维构造，可以分辨出座舱和聚变发动机等部位，还有座舱中那个卷曲的人体。在另一个二维化的人体上，可以清楚地分辨出骨骼和脉络，也可以认出身体的各个部位。在二维化的过程中，三维物体上的每个点都按照精确的几何规则投射到二维平面上，以至于这个二维体成为三维太空艇和三维人体的两张最完整最精确的图纸，其所有的内部结构都在平面上排列出来，没有任何隐藏。

这是一段非常具有启发意义的描述。因为其实降维打击可以很容易，外星人完全可以用一个巨大的苍蝇拍直接拍扁整个太阳系，不过这样就做不到“其所有的内部结构都在平面上排列出来，没有任何隐藏”。可见降维十分容易，降维同时尽可能地保留原有的结构难，例如前面的seaborn.pairplot实际上就是进行了12次四维到二维的降维，每次都用了其中的两组特征，然后扔掉了剩下的两组。于是我们就思考如何尽可能多地在二维平面保留高维特征背后的信息。

主成分分析（PCA）是其中的一种思路。PCA有两种通俗易懂的解释，一是最大化投影后数据的方差，即让数据在投影到的平面上更分散；二是最小化投影造成的损失，即让数据到投影平面的垂直距离最小。更通俗地说，就是PCA实际上是在寻找一个能把苍蝇在墙上拍得最扁、展开面积最大的一种角度。这种降维方法与二向箔给太阳系降维时使用的展示所有细节的黑科技降维方法当然没法比，但也非常有用，可以将其运用在鸢尾花的数据集上。

这种算法运用在鸢尾花数据集上的效果如下。需要读者的注意是，为了方便展示，这里用了鸢尾花数据集的前三个特征进行主成分分析，实际运用中需要展示全部特征。

其运行结果如图2-7所示。

图2-7　前三个特征进行主成分分析

如果希望知道投影平面是什么，可以通过协方差矩阵得到。

其运行结果如图2-8所示。

图2-8　通过协方差矩阵得到投影平面

由此可见，这里找到的投影平面是数据展开效果非常好的一个平面，符合PCA定义的预期。

最后，请读者思考PCA的定义。首先，投影平面是否必须是一个平面？能否是一个曲面，比如广义相对论里那种扭曲的空间。其次，这里“最大化投影后数据的方差”的目标，是否可以修改？答案当然是肯定的，这部分的详细介绍请进一步学习“流形学习”（manifold learning）相关的内容。这些方法在我们的数据集上的简单使用效果如下：

其运行结果如图2-9所示。

图2-9　运行结果