机器学习的数据清理和探索-四-机器学习的数据清理和探索（四）第十一章：第十一章：决策树和随机森林分类决策树和随机森林

机器学习的数据清理和探索（四）

原文：annas-archive.org/md5/25ad0fee8d118820a9d79ad1484952bd

译者：飞龙

协议：CC BY-NC-SA 4.0

第十一章：第十一章：决策树和随机森林分类

决策树和随机森林是非常流行的分类模型。这部分的理由是它们易于训练和解释。它们也非常灵活。我们可以建模复杂性，而无需必然增加特征空间或转换特征。我们甚至不需要对多类问题应用算法做任何特殊处理，这是我们在逻辑回归中必须做的。

另一方面，决策树可能不如其他分类模型稳定，对训练数据中的微小变化相当敏感。当存在显著的类别不平衡（一个类别的观测值比另一个类别多得多）时，决策树也可能存在偏差。幸运的是，这些问题可以通过诸如袋装法来减少方差和过采样来处理不平衡的技术来解决。我们将在本章中探讨这些技术。

在本章中，我们将涵盖以下主题：

关键概念
决策树模型
实现随机森林
实现梯度提升

技术要求

除了我们迄今为止一直在使用的 scikit-learn 模块之外，我们还将使用来自 Imbalanced-learn 的 SMOTENC。我们将使用 SMOTENC 来解决类别不平衡问题。可以使用pip install -U imbalanced-learn命令安装 Imbalanced-learn 库。本章中的所有代码都使用 scikit-learn 版本 0.24.2 和 1.0.2 进行了测试。

关键概念

决策树是一个非常有用的机器学习工具。它们是非参数的，易于解释，并且可以处理各种类型的数据。不假设特征与目标之间关系的线性，也不假设误差项的正态性。甚至不需要对数据进行缩放。决策树还经常能够很好地捕捉预测变量和目标之间的复杂关系。

决策树算法的灵活性和其建模数据中复杂和未预见的关联的能力，归功于用于分割数据的递归分割过程。决策树根据其特征值对观测值进行分组。这是通过一系列二进制决策来完成的，从根节点处的初始分割开始，以每个分组的叶子节点结束。每个分割都是基于提供关于目标最多信息的特征和特征值。更精确地说，分割的选择基于它是否产生最低的 Gini 不纯度得分。我们将在稍后更详细地讨论 Gini 不纯度。

沿着从根节点到叶子的分支，具有相同值或相同值范围的全部新观测值，都会得到相同的预测目标值。当目标是分类时，这就是该叶子节点训练观测值中目标的最频繁值。

以下图表提供了一个相当直接的决策树示例，其中包含虚构的数据和大学完成情况模型的结果。对于这个决策树，通过高中学业成绩将那些成绩在 3.0 或以下和成绩高于 3.0 的人分开，发现与其他可用特征以及其他阈值相比，这种初始分割会导致最低的不纯度。因此，高中学业成绩是我们的根节点，也称为深度 0：

图 11.1 – 完成大学教育的决策树

根节点处的二分分割导致树左侧的观测值占 45%，右侧的占 55%。在深度 1，两侧都有基于父母收入的二分分割，尽管阈值不同；左侧为 $80k，右侧为$ 45k。对于高中学业成绩大于 3 且父母收入高于$80k 的情况，没有更多的分割。在这里，我们得到了毕业的预测。这是一个叶子节点。

我们可以从每个叶子节点向上导航树，描述树是如何分割数据的，就像我们对父母收入高于 $80k 和高中学业成绩高于 3.0 的个人所做的那样。例如，决策树预测那些从学校获得低水平支持、父母收入超过$ 45k 和高中学业成绩低于或等于 3 的个人不会毕业。

那么，决策树算法是如何施展这种魔法的？它是如何选择特征和阈值或类值的？为什么是父母收入大于 $80k 或$ 45k？为什么对于父母收入低于或等于$45k 的情况，在深度 2（分割 3）处收到补助，而对于其他叶子节点的深度 2 则是学生支持水平？甚至有一个叶子节点在深度 2 都没有进一步的分割。

一种衡量二分分割对类信息提供的信息的方法是它帮助我们区分类内和类外成员资格的程度。我们经常使用基尼不纯度计算来进行这种评估，尽管有时也使用熵。基尼不纯度统计量告诉我们每个节点上类成员资格分割得有多好。这可以在以下公式中看到：

在这里，是属于k类和m个类的概率，而m是类的数量。如果一个节点上的类成员资格相等，那么基尼不纯度为 0.5。当完全纯净时，它为 0。

尝试手动计算 Gini 不纯度可能会有所帮助，以更好地理解它是如何工作的。我们可以为以下图中显示的非常简单的决策树做这件事。这里只有两个叶节点——一个是为高中 GPA 大于 3 的个人，另一个是为高中 GPA 小于或等于 3 的个人。（再次强调，这些计数是为了说明目的而编造的。我们假设图 11.1中使用的百分比是基于 100 人中的计数。）请看以下图表：

图 11.2 - 具有一个分割和 Gini 不纯度计算的决策树

根据这个模型，对于 GPA 高的个人，Graduated会被预测，因为其中大多数，45 人中的 40 人，都毕业了。Gini 不纯度相对较低，这是好的。我们可以使用前面的公式计算该节点的 Gini 不纯度：

我们的模型会预测高中 GPA 小于或等于 3 的个人未毕业，因为这些人中的大多数都没有从大学毕业。然而，这里的纯度要低得多。该节点的 Gini 不纯度值如下：

决策树算法计算从给定点开始的所有可能分割的 Gini 不纯度值的加权总和，并选择得分最低的分割。如果使用熵而不是 Gini 不纯度，算法将遵循类似的过程。在本章中，我们将使用 scikit-learn 的分类和回归树（CART）算法来构建决策树。该工具默认使用 Gini 不纯度，尽管我们可以让它使用熵。

决策树是我们所说的贪婪学习器。算法选择在当前级别给出最佳 Gini 不纯度或熵得分的分割。它不会检查该选择如何影响随后可用的分割，也不会基于该信息重新考虑当前级别的选择。这使得算法比其他情况下更有效率，但它可能不会提供全局最优解。

决策树的主要缺点是它们的方差高。它们可能会过度拟合训练数据中的异常观测值，因此在新数据上表现不佳。根据我们数据的特点，我们每次拟合决策树时都可能得到一个非常不同的模型。我们可以使用集成方法，如 bagging 或随机森林，来解决这个问题。

使用随机森林进行分类

随机森林，可能不会令人惊讶，是一系列决策树的集合。但这并不能区分随机森林和自助聚合（通常称为 bagging）。Bagging 通常用于减少具有高方差的机器学习算法（如决策树）的方差。使用 bagging，我们从数据集中生成随机样本，比如说 100 个。然后，我们在每个样本上运行我们的模型，例如决策树分类器，并对预测进行平均。

然而，使用 bagging 生成的样本可能存在相关性，并且产生的决策树可能有很多相似之处。这种情况在只有少数特征可以解释大部分变化时更为可能。随机森林通过限制每个分割可以选定的特征数量来解决这一问题。对于决策树分类模型的一个好的经验法则是取可用特征数的平方根来确定要使用的特征数。例如，如果有 25 个特征，我们会在每个分割中使用 5 个。

让我们更精确地描述构建随机森林所涉及的步骤：

从训练数据中随机采样实例（样本具有与原始数据集相同的观测数。）
随机选择样本中的特征（每次选择的好数量是可用特征总数的平方根。）
从步骤 2中随机选择的特征中确定一个特征，分割会导致具有最大纯度的节点。
内部循环：重复执行步骤 2和步骤 3，直到构建出一个决策树。
外部循环：重复所有步骤，包括内部循环，直到创建出所需数量的树。所有树的结果由投票决定；也就是说，基于所有树对于给定特征值的最高频率类别标签来预测类别。

这个过程的另一个有趣副作用是它为我们生成了测试数据。所谓的自助抽样过程——即带替换的抽样——会导致许多实例被排除在一个或多个树之外，通常多达三分之一。这些实例，被称为袋外样本，可以用来评估模型。

基于许多不相关的决策树进行分类预测，对方差（降低它！）有积极的影响，这是你可以预期的。随机森林模型通常比决策树模型更具泛化能力。它们不太容易过拟合，也不太可能被异常数据所影响。但这也带来了一定的代价。构建一百个或更多的决策树比只构建一个需要更多的系统资源。我们也失去了决策树易于解释的优点；解释每个特征的重要性变得更加困难。

使用梯度提升决策树

从概念上讲，梯度提升决策树与随机森林相似。它们依赖于多个决策树来提高模型性能。但它们是顺序执行的，每个树都从前一个树学习。每个新的树都从前一个迭代的残差开始工作。

梯度提升决策树的学习速率由ɑ超参数决定。你可能想知道为什么我们不希望我们的模型尽可能快地学习。更快的学习速率更有效率，对系统资源的消耗也更少。然而，我们可以通过降低学习速率构建一个更具泛化能力的模型。过度拟合的风险更小。最佳学习速率最终是一个经验问题。我们需要进行一些超参数调整来找到它。我们将在本章的最后部分进行这项工作。

正如随机森林的情况一样，我们可以通过进行二元分类问题的步骤来提高我们对梯度提升的直觉：

根据样本中目标变量的平均值对目标进行初步预测。对于二元目标，这是一个比例。将此预测分配给所有观测值。（我们在这里使用类别成员概率的对数。）
计算每个实例的残差，对于类内实例，将是 1 减去初始预测，对于类外实例，是 0 减去预测，或-预测。
构建一个决策树来预测残差。
基于决策树模型生成新的预测。
根据新的预测（按学习率缩放）调整每个实例的先前预测。如前所述，我们使用学习率是因为我们不希望预测移动得太快。
如果未达到最大树的数量或残差非常小，则回退到步骤 3。

虽然这是梯度提升工作原理的简化解释，但它确实为我们提供了算法所做工作的良好感觉。希望这也有助于你理解为什么梯度提升变得如此受欢迎。该算法反复调整以适应先前错误，但这样做相对高效，并且比单独的决策树有更小的过度拟合风险。

在本章的其余部分，我们将讨论决策树、随机森林和梯度提升的示例。我们将讨论如何调整超参数以及如何评估这些模型。我们还将讨论每种方法的优缺点。

决策树模型

在本章中，我们将再次使用心脏病数据。这将是一个很好的方法来比较我们的逻辑回归模型的结果与决策树等非参数模型的结果。按照以下步骤进行：

首先，我们加载到目前为止一直在使用的相同库。新的模块是来自 scikit-learn 的DecisionTreeClassifier和来自 Imbalance Learn 的SMOTENC，这将帮助我们处理不平衡数据：

import pandas as pd
import numpy as np
from sklearn.model_selection import train_test_split
from sklearn.preprocessing import OneHotEncoder
from imblearn.pipeline import make_pipeline
from sklearn.compose import ColumnTransformer
from sklearn.model_selection import RandomizedSearchCV
from imblearn.over_sampling import SMOTENC
from sklearn.tree import DecisionTreeClassifier, plot_tree
from scipy.stats import randint
import sklearn.metrics as skmet
import os
import sys
sys.path.append(os.getcwd() + "/helperfunctions")
from preprocfunc import MakeOrdinal,\
  ReplaceVals

让我们加载心脏病数据。我们将把我们的目标heartdisease转换为0和1的整数。我们不会重复之前章节中从数据中生成的频率和描述性统计信息。如果你没有完成第十章，逻辑回归，快速查看它们可能是有帮助的：

healthinfo = pd.read_csv("data/healthinfosample.csv")
healthinfo.set_index("personid", inplace=True)
healthinfo.heartdisease.value_counts()
No     27467
Yes     2533
Name: heartdisease, dtype: int64
healthinfo['heartdisease'] = \
  np.where(healthinfo.heartdisease=='No',0,1).\
  astype('int')
healthinfo.heartdisease.value_counts()
0    27467
1     2533
Name: heartdisease, dtype: int64

注意类别不平衡。不到 10%的观察结果患有心脏病。我们将在模型中处理这个问题。

让我们再看看年龄类别特征的值，因为它有点不寻常。它包含年龄范围的字符数据。我们将使用我们加载的MakeOrdinal类将其转换为有序特征。例如，18-24的值将给我们一个0的值，而50-54的值在转换后将变为6：

healthinfo.agecategory.value_counts().\
  sort_index().reset_index()
          index  agecategory
0         18-24         1973
1         25-29         1637
2         30-34         1688
3         35-39         1938
4         40-44         2007
5         45-49         2109
6         50-54         2402
7         55-59         2789
8         60-64         3122
9         65-69         3191
10        70-74         2953
11        75-79         2004
12  80 or older         2187

我们应该按数据类型组织我们的特征，这样将使一些任务更容易。我们还将设置一个字典来重新编码genhealth和diabetic特征：

num_cols = ['bmi','physicalhealthbaddays',
   'mentalhealthbaddays','sleeptimenightly']
binary_cols = ['smoking','alcoholdrinkingheavy',
  'stroke','walkingdifficult','physicalactivity',
  'asthma','kidneydisease','skincancer']
cat_cols = ['gender','ethnicity']
spec_cols1 = ['agecategory']
spec_cols2 = ['genhealth','diabetic']
rep_dict = {
  'genhealth': {'Poor':0,'Fair':1,'Good':2,
    'Very good':3,'Excellent':4},
  'diabetic': {'No':0,
    'No, borderline diabetes':0,'Yes':1,
    'Yes (during pregnancy)':1}           
}

现在，让我们创建训练和测试数据框：

X_train, X_test, y_train, y_test =  \
  train_test_split(healthinfo[num_cols + 
    binary_cols + cat_cols + spec_cols1 +
    spec_cols2],\
  healthinfo[['heartdisease']], test_size=0.2,
    random_state=0)

接下来，我们将设置列转换。我们将使用我们的自定义类将agecategory特征编码为有序，并将genhealth和diabetic的特征的字符值替换为数值。

我们不会转换数值列，因为在使用决策树时通常不需要缩放这些特征。我们也不会担心异常值，因为决策树对它们的敏感性低于逻辑回归。我们将设置remainder为passthrough，以便转换器将剩余的列（数值列）原样通过：

ohe = OneHotEncoder(drop='first', sparse=False)
spectrans1 = make_pipeline(MakeOrdinal())
spectrans2 = make_pipeline(ReplaceVals(rep_dict))
bintrans = make_pipeline(ohe)
cattrans = make_pipeline(ohe)
coltrans = ColumnTransformer(
  transformers=[
    ("bin", bintrans, binary_cols),
    ("cat", cattrans, cat_cols),
    ("spec1", spectrans1, spec_cols1),
    ("spec2", spectrans2, spec_cols2),
  ],
    remainder = 'passthrough'
)

在运行模型之前，我们需要做一些工作。正如你将在下一步看到的，我们需要知道单变量编码器将返回多少个特征。同时，我们也应该获取新的特征名称。我们稍后还需要它们。（我们只需要对数据的一个小随机样本进行列转换器拟合即可。）看看下面的代码：
```
coltrans.fit(X_train.sample(1000))
new_binary_cols = \
  coltrans.\
  named_transformers_['bin'].\
  named_steps['onehotencoder'].\
  get_feature_names(binary_cols)
new_cat_cols = \
  coltrans.\
  named_transformers_['cat'].\
  named_steps['onehotencoder'].\
  get_feature_names(cat_cols)
```

让我们查看特征名称：

new_cols = np.concatenate((new_binary_cols, 
  new_cat_cols, np.array(spec_cols1 + spec_cols2 +
  num_cols)))
new_cols
array(['smoking_Yes', 'alcoholdrinkingheavy_Yes',
       'stroke_Yes', 'walkingdifficult_Yes',
       'physicalactivity_Yes', 'asthma_Yes',
       'kidneydisease_Yes', 'skincancer_Yes',
       'gender_Male', 'ethnicity_Asian',
       'ethnicity_Black', 'ethnicity_Hispanic',
       'ethnicity_Other', 'ethnicity_White',
       'agecategory', 'genhealth', 'diabetic', 'bmi',
       'physicalhealthbaddays', 'mentalhealthbaddays',
       'sleeptimenightly'], dtype=object)

在我们拟合决策树之前，我们需要处理我们的不平衡数据集。我们可以使用来自 Imbalanced-learn 的SMOTENC模块来对心脏病类别进行过采样。这将生成足够的心脏病类别的代表性实例，以平衡类别成员资格。

接下来，我们必须实例化一个决策树分类器，并指出叶子节点需要至少有五个观察结果，并且树的深度不能超过两个。然后，我们将使用列转换器转换训练数据并拟合模型。

我们将在本节稍后进行一些超参数调整。现在，我们只想生成一个易于解释和可视化的决策树。

SMOTENC需要知道哪些列是分类的。当我们设置列转换器时，我们首先编码二进制列，然后是分类列。因此，二进制列数加上分类列数给出了那些列索引的终点。然后，我们必须传递一个范围，从 0 开始，到分类列数结束，到SMOTENC的categorical_features参数。

现在，我们可以创建一个包含列转换、过采样和决策树分类器的管道，并对其进行拟合：

catcolscnt = new_binary_cols.shape[0] + \
  new_cat_cols.shape[0]
smotenc = \
  SMOTENC(categorical_features=np.arange(0,catcolscnt),
  random_state=0)
dtc_example = DecisionTreeClassifier(
  min_samples_leaf=5, max_depth=2)
pipe0 = make_pipeline(coltrans, smotenc, dtc_example)
pipe0.fit(X_train, y_train.values.ravel())

注意

当我们担心我们的模型在捕捉一个类别的变化方面做得不好，因为我们有太多少的该类实例，相对于一个或多个其他类别时，过采样可以是一个好的选择。过采样会复制该类别的实例。

合成少数过采样技术（SMOTE）是一个使用 KNN 来复制实例的算法。本章中 SMOTE 的实现来自 Imbalanced-learn，特别是 SMOTENC，它可以处理分类数据。

当类不平衡比这个数据集更严重时，通常会进行过采样，比如 100 到 1。尽管如此，我认为在本章中演示如何使用 SMOTE 和类似工具是有帮助的。

运行拟合后，我们可以查看哪些特征被识别为重要。agecategory、genhealth和diabetic是此简单模型中的重要特征：

feature_imp = \
  pipe0.named_steps['decisiontreeclassifier'].\
  tree_.compute_feature_importances(normalize=False)
feature_impgt0 = feature_imp>0
feature_implabs = np.column_stack((feature_imp.\
  ravel(), new_cols))
feature_implabs[feature_impgt0]
array([[0.10241844433036575, 'agecategory'],
       [0.04956947743193013, 'genhealth'],
       [0.012777650193266089, 'diabetic']],
      dtype=object)

接下来，我们可以生成决策树的图形：

plot_tree(pipe0.named_steps['decisiontreeclassifier'],
  feature_names=new_cols, 
  class_names=['No Disease','Disease'], fontsize=10)

这会产生以下图形：

图 11.3 – 以心脏病为目标的心脏病决策树示例

初始的二分分割，在根节点（也称为深度 0），是基于agecategory是否小于或等于6。 (回想一下，agecategory最初是一个字符特征。编码后的初始值50-54得到6的值。) 如果根节点语句为真，它将导致下一级节点向左。如果语句为假，它将导致下一级节点向右。样本数是该节点到达的观察数。因此，diabetic<=0.001（即非糖尿病患者）节点上的样本值12576反映了从父节点中得到的语句为真的实例数；也就是说，有12576个实例的年龄类别值小于或等于6。

每个节点内的value列表给我们提供了训练数据中每个类别的实例。在这种情况下，第一个值是无疾病观察值的计数。第二个值是有心脏病观察值的计数。例如，在diabetic<=0.001节点上，有10781个无疾病观察值和1795个疾病观察值。

决策树在叶子模式下预测最频繁的类别。因此，这个模型会预测 54 岁或以下且不是糖尿病患者的个体没有疾病（agecategory<=6和diabetic<=0.001）。在训练数据中，该组有10142个无疾病观察结果，990个疾病观察结果。这给我们一个具有非常好的基尼不纯度0.162的叶子节点。

如果一个人有糖尿病，即使他们 54 岁或以下，我们的模型也会预测心脏病。然而，这种预测并不那么确定。基尼不纯度为0.493。相比之下，预测 54 岁以上（agecategory<=6.001为假）且健康状况不佳（genhealth<=3.0）的个体的疾病，其基尼不纯度显著较低。

我们的模型预测 54 岁以上且一般健康状况等于4的人没有疾病。（当我们进行列转换时，我们将一般健康状况值Excellent编码为4。）然而，较差的基尼不纯度分数表明，我们的模型并不完全有信心做出那个预测。

让我们看看这个模型的某些指标。该模型有不错的敏感性，但不是很好，大约 70%的时间预测有心脏病。当我们做出积极预测时，精确度相当低。只有 19%的时间我们做出的积极预测是正确的：

pred = pipe0.predict(X_test)
print("accuracy: %.2f, sensitivity: %.2f, specificity: %.2f, precision: %.2f"  %
  (skmet.accuracy_score(y_test.values.ravel(), pred),
  skmet.recall_score(y_test.values.ravel(), pred),
  skmet.recall_score(y_test.values.ravel(), pred,
    pos_label=0),
  skmet.precision_score(y_test.values.ravel(), pred)))
accuracy: 0.72, sensitivity: 0.70, specificity: 0.73, precision: 0.19

这个模型在超参数方面做出了一些相当随意的决定，将最大深度设置为2，并将叶子的最小样本数设置为5。我们应该探索这些超参数的替代值，以获得性能更好的模型。让我们进行随机网格搜索以找到这些值。

让我们设置一个包含列转换、过采样和决策树分类器的管道。我们还将创建一个包含最小叶子大小和最大树深度超参数范围的字典。请注意，对于每个字典键，decisiontreeclassifier后面有两个下划线：
```
dtc = DecisionTreeClassifier(random_state=0)
pipe1 = make_pipeline(coltrans, smotenc, dtc)
dtc_params = {
 'decisiontreeclassifier__min_samples_leaf': randint(100, 1200),
 'decisiontreeclassifier__max_depth': randint(2, 11)
}
```

现在，我们可以运行随机网格搜索。我们将运行 20 次迭代以测试我们超参数的多个值：

rs = RandomizedSearchCV(pipe1, dtc_params, cv=5,
  n_iter=20, scoring="roc_auc")
rs.fit(X_train, y_train.values.ravel())
rs.best_params_
{'decisiontreeclassifier__max_depth': 9,
 'decisiontreeclassifier__min_samples_leaf': 954}
rs.best_score_
0.7964540832005679

让我们看看每次迭代的得分：

results = \
  pd.DataFrame(rs.cv_results_['mean_test_score'], \
    columns=['meanscore']).\
  join(pd.DataFrame(rs.cv_results_['params'])).\
  sort_values(['meanscore'], ascending=False).\
  rename(columns=\
    {'decisiontreeclassifier__max_depth':'maxdepth',
     'decisiontreeclassifier__min_samples_leaf':\
     'samples'})

这产生了以下输出。表现最好的模型与我们之前构建的模型相比，max_depth有显著增加。我们的模型在叶子中每个叶子的最小实例数也更高：

    meanscore  maxdepth  samples
15      0.796         9      954
13      0.796         8      988
4       0.795         7      439
19      0.795         9      919
12      0.794         9      856
3       0.794         9      510
2       0.794         9     1038
5       0.793         8      575
0       0.793        10     1152
10      0.793         7     1080
6       0.793         6     1013
8       0.793        10      431
17      0.793         6      896
14      0.792         6      545
16      0.784         5      180
1       0.778         4      366
11      0.775         4      286
9       0.773         4      138
18      0.768         3      358
7       0.765         3      907

让我们生成一个混淆矩阵：

pred2 = rs.predict(X_test)
cm = skmet.confusion_matrix(y_test, pred2)
cmplot = \
  skmet.ConfusionMatrixDisplay(confusion_matrix=cm,
  display_labels=['Negative', 'Positive'])
cmplot.plot()
cmplot.ax_.\
  set(title='Heart Disease Prediction Confusion Matrix', 
  xlabel='Predicted Value', ylabel='Actual Value')

这产生了以下图表：

图 11.4 – 决策树模型的冠心病混淆矩阵

这张图你可能立刻注意到的一个问题是我们的精确度有多低。绝大多数时间我们预测为阳性，我们都是错误的。

让我们查看准确率、灵敏度、特异性和精确度分数，看看模型在没有超参数调整的情况下，这些指标是否有很大改进。我们在灵敏度上表现明显更差，现在为 63%，而之前为 70%。然而，我们在特异性上做得稍微好一些。我们现在在测试数据上正确预测负面的概率为 79%，而之前的模型为 73%：

print("accuracy: %.2f, sensitivity: %.2f, specificity: %.2f, precision: %.2f"  %
  (skmet.accuracy_score(y_test.values.ravel(), pred2),
   skmet.recall_score(y_test.values.ravel(), pred2),
   skmet.recall_score(y_test.values.ravel(), pred2,
     pos_label=0),
   skmet.precision_score(y_test.values.ravel(), 
     pred2)))
accuracy: 0.77, sensitivity: 0.63, specificity: 0.79, precision: 0.21

决策树是分类模型的良好起点。它们对底层数据几乎没有假设，并且不需要太多的预处理。注意，在这个例子中我们没有进行任何缩放或异常值检测，因为在决策树中这通常不是必要的。我们还得到一个相当容易理解或解释的模型，只要我们限制深度数量。

我们通常可以通过随机森林来提高我们的决策树模型的性能，原因我们在本章开头已经讨论过。一个关键的原因是，当使用随机森林而不是决策树时，方差会降低。

在下一节中，我们将探讨随机森林。

实现随机森林

让我们尝试使用随机森林来提高我们的心脏病模型：

首先，让我们加载与上一节相同的库，但这次我们将导入随机森林分类器：

import pandas as pd
import numpy as np
from imblearn.pipeline import make_pipeline
from sklearn.ensemble import RandomForestClassifier
from sklearn.model_selection import RandomizedSearchCV
from scipy.stats import randint
import sklearn.metrics as skmet
import os
import sys
sys.path.append(os.getcwd() + "/helperfunctions")
import healthinfo as hi

我们还加载了healthinfo模块；它加载健康信息数据并执行我们的预处理。这里没有太多花哨的东西。我们之前步骤中使用的预处理代码只是复制到了当前工作目录的helperfunctions子文件夹中。

现在，让我们获取由healthinfo模块处理过的数据，以便我们可以使用它来为我们的随机森林分类器：
```
X_train = hi.X_train
X_test = hi.X_test
y_train = hi.y_train
y_test = hi.y_test
```

让我们实例化一个随机森林分类器并为网格搜索创建一个管道。我们还将创建一个用于搜索的超参数字典。除了我们用于决策树的max_depth和min_samples_leaf超参数外，随机森林的一个重要超参数是n_estimators。这表示用于该搜索迭代的树的数量。我们还将添加entropy作为标准，除了我们迄今为止使用的gini：

rfc = RandomForestClassifier(random_state=0)
pipe1 = make_pipeline(hi.coltrans, hi.smotenc, rfc)
rfc_params = {
 'randomforestclassifier__min_samples_leaf':
    randint(100, 1200),
 'randomforestclassifier__max_depth': 
    randint(2, 11),
 'randomforestclassifier__n_estimators': 
    randint(100, 3000),
 'randomforestclassifier__criterion': 
    ['gini','entropy']
}
rs = RandomizedSearchCV(pipe1, rfc_params, cv=5, 
  n_iter=20, scoring="roc_auc")
rs.fit(X_train, y_train.values.ravel())

我们可以使用随机网格搜索对象的best_params_和best_score_属性来找到最佳参数和相应的分数。最佳模型有1023棵树和最大深度为9。

在这里，我们可以看到roc_auc分数相对于上一节中的决策树模型有所提高：

rs.best_params_
{'randomforestclassifier__criterion': 'gini',
 'randomforestclassifier__max_depth': 9,
 'randomforestclassifier__min_samples_leaf': 667,
 'randomforestclassifier__n_estimators': 1023}
rs.best_score_
0.8210934290375318

随机森林的结果比单个决策树的结果更难以解释，但一个好的起点是查看特征重要性。前三个特征与我们在决策树中看到的是相同的——即agecategory、genhealth和diabetic：

feature_imp = \
  rs.best_estimator_['randomforestclassifier'].\
  feature_importances_
feature_implabs = np.column_stack((feature_imp.\
  ravel(), hi.new_cols))
pd.DataFrame(feature_implabs,
  columns=['importance','feature']).\
  sort_values(['importance'], ascending=False)
         importance       feature
14       0.321            agecategory
15       0.269            genhealth
16       0.159            diabetic
13       0.058            ethnicity_White
0        0.053            smoking_Yes
18       0.033            physicalhealthbaddays
8        0.027            gender_Male
3        0.024            walkingdifficult_Yes
20       0.019            sleeptimenightly
11       0.010            ethnicity_Hispanic
17       0.007            bmi
19       0.007            mentalhealthbaddays
1        0.007            alcoholdrinkingheavy_Yes
5        0.003            asthma_Yes
10       0.002            ethnicity_Black
4        0.001            physicalactivity_Yes
7        0.001            skincancer_Yes
2        0.000            stroke_Yes
6        0.000            kidneydisease_Yes
9        0.000            ethnicity_Asian
12       0.000            ethnicity_Other

让我们看看一些指标。与先前的模型相比，敏感性有所提高，但在其他任何指标上都没有太大变化。我们保持了相同的相对良好的特异性分数。总体而言，这是我们迄今为止最好的模型：

print("accuracy: %.2f, sensitivity: %.2f, specificity: %.2f, precision: %.2f"  %
  (skmet.accuracy_score(y_test.values.ravel(), pred),
  skmet.recall_score(y_test.values.ravel(), pred),
  skmet.recall_score(y_test.values.ravel(), pred,
    pos_label=0),
  skmet.precision_score(y_test.values.ravel(), pred)))
accuracy: 0.77, sensitivity: 0.69, specificity: 0.78, precision: 0.22

我们可能仍然能够提高我们模型性能指标。我们至少应该尝试一些梯度提升。正如我们在本章的使用梯度提升决策树部分所讨论的，梯度提升决策树有时可能比随机森林产生更好的模型。这是因为每个树都是从先前树的错误中学习的。

实现梯度提升

在本节中，我们将尝试使用梯度提升来改进我们的随机森林模型。我们必须注意的一个问题是过拟合，这在梯度提升决策树中可能比在随机森林中更成问题。这是因为随机森林的树不会从其他树中学习，而梯度提升中，每个树都是基于先前树的学习的。我们在这里选择的超参数至关重要。让我们开始吧：

我们将首先导入必要的库。我们将使用与随机森林相同的模块，但我们将从ensemble导入GradientBoostingClassifier而不是RandomForestClassifier：

import pandas as pd
import numpy as np
from imblearn.pipeline import make_pipeline
from sklearn.model_selection import RandomizedSearchCV
from sklearn.ensemble import GradientBoostingClassifier
import sklearn.metrics as skmet
from scipy.stats import uniform
from scipy.stats import randint
import matplotlib.pyplot as plt
import os
import sys
sys.path.append(os.getcwd() + "/helperfunctions")
import healthinfo as hi

现在，让我们获取由healthinfo模块处理过的数据，用于我们的随机森林分类器：
```
X_train = hi.X_train
X_test = hi.X_test
y_train = hi.y_train
y_test = hi.y_test
```
接下来，我们将实例化一个梯度提升分类器实例，并将其添加到一个管道中，包括我们用于预处理健康数据的步骤（这些步骤在我们导入的模块中）。

我们将创建一个包含梯度提升分类器超参数的字典。这些包括熟悉的每个叶子的最小样本数、最大深度以及估计器的数量超参数。我们还将添加用于检查学习率的值：

gbc = GradientBoostingClassifier(random_state=0)
pipe1 = make_pipeline(hi.coltrans, hi.smotenc, gbc)
gbc_params = {
 'gradientboostingclassifier__min_samples_leaf':
     randint(100, 1200),
 'gradientboostingclassifier__max_depth':
     randint(2, 20),
 'gradientboostingclassifier__learning_rate':
     uniform(loc=0.02, scale=0.25),
 'gradientboostingclassifier__n_estimators':
     randint(100, 1200)
}

现在，我们准备进行网格搜索。（注意，这可能在您的机器上运行需要一些时间。）在我们运行的七次迭代中，最佳模型的学习率为0.25，估计器或树的数量为308。它有一个相当不错的roc_auc分数为0.82：

rs = RandomizedSearchCV(pipe1, gbc_params, cv=5, 
  n_iter=7, scoring="roc_auc")
rs.fit(X_train, y_train.values.ravel())
rs.best_params_
{'gradientboostingclassifier__learning_rate': 0.2528,
 'gradientboostingclassifier__max_depth': 3,
 'gradientboostingclassifier__min_samples_leaf': 565,
 'gradientboostingclassifier__n_estimators': 308}
rs.best_score_
0.8162378796382679

让我们看看特征重要性：

feature_imp = pd.Series(rs.\
  best_estimator_['gradientboostingclassifier'].\
  feature_importances_, index=hi.new_cols)
feature_imp.loc[feature_imp>0.01].\
    plot(kind='barh')
plt.tight_layout()    
plt.title('Gradient Boosting Feature Importance')

这会产生以下图表：

图 11.5 – 梯度提升特征重要性

让我们看看一些指标。有趣的是，我们得到了出色的准确性和特异性，但敏感性却非常糟糕。这可能是因为过拟合：

pred = rs.predict(X_test)
print("accuracy: %.2f, sensitivity: %.2f, specificity: %.2f, precision: %.2f"  %
  (skmet.accuracy_score(y_test.values.ravel(), pred),
  skmet.recall_score(y_test.values.ravel(), pred),
  skmet.recall_score(y_test.values.ravel(), pred,
    pos_label=0),
  skmet.precision_score(y_test.values.ravel(), pred)))
accuracy: 0.91, sensitivity: 0.19, specificity: 0.97, precision: 0.40

虽然我们的结果参差不齐，但梯度提升决策树通常是一个很好的分类模型选择。这尤其适用于我们正在建模特征与目标之间的复杂关系时。当决策树是一个合适的选择，但我们担心决策树模型的高方差时，梯度提升决策树至少与随机森林一样是一个好的选择。

现在，让我们总结一下本章学到的内容。

摘要

本章探讨了如何使用决策树来解决分类问题。尽管本章中的所有示例都涉及二元目标，但我们使用的算法也可以处理多类问题。与从逻辑回归到多项式逻辑回归的转变不同，当我们的目标值超过两个时，为了有效地使用这些算法，我们通常不需要做出太多改变。

我们探讨了两种处理决策树高方差的方法。一种方法是使用随机森林，这是一种袋装方法。这将减少我们预测中的方差。另一种方法是使用梯度提升决策树。提升可以帮助我们捕捉数据中的非常复杂的关系，但存在非平凡的过拟合风险。考虑到这一点，调整我们的超参数尤为重要。

在下一章中，我们将探讨另一个著名的分类算法：K 最近邻算法。

第十二章：第十二章：用于分类的 K-最近邻

K-最近邻 (KNN) 是当观察或特征不多，且预测类别成员不需要非常高效时，用于分类模型的一个很好的选择。它是一个懒惰学习器，因此比其他分类算法更快地拟合，但在对新观察进行分类时则慢得多。它也可能在极端情况下产生不太准确的预测，但通过适当地调整 k 可以改进这一点。我们将在本章开发的模型中仔细考虑这些选择。

KNN 可能是我们能选择的最为直接的非参数算法之一，使其成为一个良好的诊断工具。不需要对特征的分布或特征与目标之间的关系做出任何假设。没有很多超参数需要调整，而且两个关键的超参数——最近邻和距离度量——都很容易理解。

KNN 可以成功用于二元和多类问题，而无需对算法进行任何扩展。

在本章中，我们将涵盖以下主题：

KNN 的关键概念
二元分类的 KNN
多类分类的 KNN

技术要求

除了常用的 scikit-learn 库之外，我们还需要 imblearn（不平衡学习）库来运行本章中的代码。这个库帮助我们处理显著的类别不平衡。imblearn 可以通过 pip install imbalanced-learn 安装，或者如果你使用 Anaconda，可以通过 conda install -c conda-forge imbalanced-learn 安装。所有代码都已使用 scikit-learn 版本 0.24.2 和 1.0.2 进行测试。

KNN 的关键概念

KNN 可能是我们将在本书中讨论的最直观的算法。其思想是找到 k 个属性最相似的实例，其中这种相似性对目标很重要。最后一个条款是一个重要但可能显然的限定条件。我们关注与目标值相关的属性之间的相似性。

对于每个需要预测目标的观察，KNN 会找到与该观察的特征最相似的 k 个训练观察。当目标是分类时，KNN 会选择 k 个训练观察中目标的最频繁值。（我们通常选择奇数个 k 以避免分类问题中的平局。）

我所说的通过训练观察到的，是指那些具有已知目标值的观察。KNN 由于是一个懒惰学习器，所以不需要进行真正的训练。我将在本节中更详细地讨论这一点。

以下图表说明了使用 KNN 进行分类，其中 k 的值为 1 和 3。当 k=1 时，我们会预测新的观察值 X 将属于圆形类别。当 k=3 时，它将被分配到正方形类别：

图 12.1 – k 的值为 1 和 3 的 KNN

但我们所说的相似或最近的实例是什么意思呢？有几种方法可以衡量相似性，但最常用的度量是欧几里得距离。欧几里得距离是两点之间平方差的和。这可能会让你想起勾股定理。从点 a 到点 b 的欧几里得距离如下：

欧几里得距离的一个合理的替代方案是曼哈顿距离。从点 a 到点 b 的曼哈顿距离如下：

scikit-learn 中的默认距离度量是闵可夫斯基距离。从点 a 到点 b 的闵可夫斯基距离如下：

注意到当 p 为 1 时，它与曼哈顿距离相同。当 p 为 2 时，它与欧几里得距离相同。

曼哈顿距离有时被称为出租车距离。这是因为它反映了两个点在网格路径上的距离。以下图表说明了曼哈顿距离并将其与欧几里得距离进行了比较：

图 12.2 – 欧几里得和曼哈顿距离度量

使用曼哈顿距离可以在特征类型或尺度差异很大时产生更好的结果。然而，我们可以将距离度量的选择视为一个经验问题；也就是说，我们可以尝试两者（或其他的距离度量）并看看哪个给我们带来性能最好的模型。我们将在下一节通过网格搜索来演示这一点。

如你所怀疑的那样，KNN 模型对 k 的选择很敏感。较低的 k 值会导致一个试图识别观察之间细微差异的模型。在非常低的 k 值时，存在过度拟合的实质性风险。但在 k 值较高时，我们的模型可能不够灵活。我们再次面临方差-偏差权衡。较低的 k 值导致偏差较少而方差较多，而较高的值则相反。

对于 k 的选择没有明确的答案。但一个好的经验法则是使用观察数的平方根。然而，就像我们对距离度量所做的那样，我们应该测试模型在不同 k 值下的性能。KNN 是一种非参数算法。不对底层数据的属性做出假设，例如线性或正态分布的特征。这使得 KNN 非常灵活。它可以用来模拟特征与目标之间的各种关系。

如前所述，KNN 是一种懒惰学习算法。在训练时间不进行任何计算。学习仅在测试时发生。这有其优点和缺点。当数据中有许多实例或维度时，它可能不是一个好的选择，而且预测速度很重要。KNN 也往往在稀疏数据上表现不佳，例如包含许多 0 值的数据集。

在下一节中，我们将使用 KNN 构建一个二分类模型，然后在下一节构建几个多分类模型。

KNN 用于二分类

KNN 算法与决策树算法有一些相同的优点。不需要满足关于特征或残差的分布的先验假设。它是我们试图在前两章中构建的心脏病模型的一个合适的算法。数据集不是很大（30,000 个观测值）并且没有太多特征。

注意

心脏病数据集可在www.kaggle.com/datasets/kamilpytlak/personal-key-indicators-of-heart-disease公开下载。它来源于 2020 年美国疾病控制与预防中心对超过 40 万人的调查数据。我已经从这个数据集中随机抽取了 30,000 个观测值用于本节的分析。数据列包括受访者是否曾经患有心脏病、体重指数、吸烟史、大量饮酒、年龄、糖尿病和肾病。

让我们开始构建我们的模型：

首先，我们必须加载我们在过去几章中使用的一些相同的库。我们还将加载KneighborsClassifier：

import pandas as pd
import numpy as np
from imblearn.pipeline import make_pipeline
from sklearn.model_selection import RandomizedSearchCV,\
  RepeatedStratifiedKFold
from sklearn.neighbors import KNeighborsClassifier
from sklearn.feature_selection import SelectKBest, chi2
from scipy.stats import randint
import sklearn.metrics as skmet
from sklearn.model_selection import cross_validate
import os
import sys
sys.path.append(os.getcwd() + "/helperfunctions")
import healthinfo as hi

healthinfo模块包含了我们在第十章，逻辑回归中使用的所有代码，用于加载健康信息数据并进行预处理。这里没有必要重复这些步骤。如果你还没有阅读第十章，逻辑回归，至少浏览一下该章节的第二部分的代码可能会有所帮助。这将让你更好地了解特征。

现在，让我们获取由healthinfo模块处理过的数据并显示特征名称：

X_train = hi.X_train
X_test = hi.X_test
y_train = hi.y_train
y_test = hi.y_test
new_cols = hi.new_cols
new_cols
array(['smoking_Yes', 'alcoholdrinkingheavy_Yes',
       'stroke_Yes', 'walkingdifficult_Yes',
       'physicalactivity_Yes', 'asthma_Yes',
       'kidneydisease_Yes', 'skincancer_Yes',
       'gender_Male', 'ethnicity_Asian',
       'ethnicity_Black', 'ethnicity_Hispanic',
       'ethnicity_Other', 'ethnicity_White',
       'agecategory', 'genhealth', 'diabetic', 'bmi',
       'physicalhealthbaddays', 'mentalhealthbaddays',
       'sleeptimenightly'], dtype=object)

我们可以使用 K 折交叉验证来评估这个模型。我们已经在第六章，准备模型评估中讨论了 K 折交叉验证。我们将指定我们想要重复 10 次 10 个分割。默认值分别是5和10。

我们模型的精确度，即我们预测心脏病时的正确率，异常低，为0.17。灵敏度，即存在心脏病时预测心脏病的比率，也较低，为0.56：

knn_example = KNeighborsClassifier(n_neighbors=5, n_jobs=-1)
kf = RepeatedStratifiedKFold(n_splits=10, n_repeats=10, random_state=0)
pipe0 = make_pipeline(hi.coltrans, hi.smotenc, knn_example)
scores = cross_validate(pipe0, X_train,
  y_train.values.ravel(), \
  scoring=['accuracy','precision','recall','f1'], \
  cv=kf, n_jobs=-1)
print("accuracy: %.2f, sensitivity: %.2f, precision: %.2f, f1: %.2f"  %
  (np.mean(scores['test_accuracy']),\
  np.mean(scores['test_recall']),\
  np.mean(scores['test_precision']),\
  np.mean(scores['test_f1'])))
accuracy: 0.73, sensitivity: 0.56, precision: 0.17, f1: 0.26

我们可以通过一些超参数调整来提高我们模型的性能。让我们为几个邻居和距离度量创建一个字典。我们还将尝试使用我们的filter方法选择不同数量的特征：

knn = KNeighborsClassifier(n_jobs=-1)
pipe1 = make_pipeline(hi.coltrans, hi.smotenc,
   SelectKBest(score_func=chi2), knn)
knn_params = {
 'selectkbest__k':
    randint(1, len(new_cols)),
 'kneighborsclassifier__n_neighbors':
    randint(5, 300),
 'kneighborsclassifier__metric':
    ['euclidean','manhattan','minkowski']
}

我们将在网格搜索中的评分将基于接收者操作特征曲线（ROC 曲线）下的面积。我们在第六章，准备模型评估中介绍了 ROC 曲线：
```
rs = RandomizedSearchCV(pipe1, knn_params, cv=5, scoring="roc_auc")
rs.fit(X_train, y_train.values.ravel())
```

我们可以使用随机网格搜索的最佳估计器属性从selectkbest获取选定的特征：

selected = rs.best_estimator_['selectkbest'].\
  get_support()
selected.sum()
11
new_cols[selected]
array(['smoking_Yes', 'alcoholdrinkingheavy_Yes',
       'walkingdifficult_Yes', 'ethnicity_Black',
       'ethnicity_Hispanic', 'agecategory',
       'genhealth', 'diabetic', 'bmi',
       'physicalhealthbaddays','mentalhealthbaddays'],
      dtype=object)

我们还可以查看最佳参数和最佳得分。11 个特征（17 个特征中的 11 个）被选中，正如我们在上一步中看到的。一个k（n_neighbors）为254和曼哈顿距离度量是得分最高的模型的另一个超参数：
```
rs.best_params_
{'kneighborsclassifier__metric': 'manhattan',
 'kneighborsclassifier__n_neighbors': 251,
 'selectkbest__k': 11}
rs.best_score_
0.8030553205304845
```

让我们看看这个模型的更多指标。我们在敏感性方面做得很好，但其他指标并不好：

pred = rs.predict(X_test)
print("accuracy: %.2f, sensitivity: %.2f, specificity: %.2f, precision: %.2f"  %
  (skmet.accuracy_score(y_test.values.ravel(), pred),
  skmet.recall_score(y_test.values.ravel(), pred),
  skmet.recall_score(y_test.values.ravel(), pred, 
    pos_label=0),
  skmet.precision_score(y_test.values.ravel(), pred)))
accuracy: 0.67, sensitivity: 0.82, specificity: 0.66, precision: 0.18

我们还应该绘制混淆矩阵。为此，我们可以查看相对较好的敏感性。在这里，我们正确地将大多数实际阳性识别为阳性。然而，这是以许多假阳性为代价的。我们可以从上一步的精确度得分中看到这一点。大多数时候我们预测阳性，我们都是错误的：
```
cm = skmet.confusion_matrix(y_test, pred)
cmplot = skmet.ConfusionMatrixDisplay(
  confusion_matrix=cm,
  display_labels=['Negative', 'Positive'])
cmplot.plot()
cmplot.ax_.set(title='Heart Disease Prediction Confusion Matrix', 
  xlabel='Predicted Value', ylabel='Actual Value')
```

这会产生以下图表：

图 12.3 – 超参数调整后的心脏病混淆矩阵

在本节中，你学习了如何使用具有二进制目标的 KNN。我们可以遵循非常相似的步骤来使用 KNN 进行多类分类问题。

KNN 多类分类

构建 KNN 多类模型相当简单，因为它不需要对算法进行特殊扩展，例如将逻辑回归应用于具有两个以上值的目标所需的扩展。我们可以通过使用我们在第十章，逻辑回归中的多项式逻辑回归部分使用的相同机器故障数据来看到这一点。

注意

这份关于机器故障的数据集可在www.kaggle.com/datasets/shivamb/machine-predictive-maintenance-classification公开使用。有 10,000 个观测值，12 个特征，以及两个可能的靶标。一个是二进制靶标，指定机器是否故障。另一个包含故障类型。该数据集中的实例是合成的，由一个旨在模仿机器故障率和原因的过程生成。

让我们构建我们的机器故障类型模型：

首先，让我们加载现在熟悉的模块：

import pandas as pd
import numpy as np
from sklearn.model_selection import train_test_split
from sklearn.preprocessing import OneHotEncoder, MinMaxScaler
from imblearn.pipeline import make_pipeline
from sklearn.impute import SimpleImputer
from sklearn.compose import ColumnTransformer
from sklearn.model_selection import RandomizedSearchCV
from sklearn.neighbors import KNeighborsClassifier
from imblearn.over_sampling import SMOTENC
from sklearn.feature_selection import SelectKBest, chi2
import sklearn.metrics as skmet
import os
import sys
sys.path.append(os.getcwd() + "/helperfunctions")
from preprocfunc import OutlierTrans

让我们加载机器故障数据并查看其结构。共有 10,000 个观测值，没有缺失数据。数据包括分类数据和数值数据的组合：

machinefailuretype = pd.read_csv("data/machinefailuretype.csv")
machinefailuretype.info()
<class 'pandas.core.frame.DataFrame'>
RangeIndex: 10000 entries, 0 to 9999
Data columns (total 10 columns):
 #   Column               Non-Null Count  Dtype
---  ------               --------------  -----  
 0   udi                  10000 non-null  int64
 1   product              10000 non-null  object 
 2   machinetype          10000 non-null  object 
 3   airtemp              10000 non-null  float64
 4   processtemperature   10000 non-null  float64
 5   rotationalspeed      10000 non-null  int64
 6   torque               10000 non-null  float64
 7   toolwear             10000 non-null  int64
 8   fail                 10000 non-null  int64
 9   failtype             10000 non-null  object 
dtypes: float64(3), int64(4), object(3)
memory usage: 781.4+ KB

让我们也看看一些观测值：

machinefailuretype.head()
   udi product machinetype airtemp processtemperature\
0   1   M14860      M         298         309 
1   2   L47181      L         298         309 
2   3   L47182      L         298         308 
3   4   L47183      L         298         309 
4   5   L47184      L         298         309 
 Rotationalspeed  Torque  toolwear  fail   failtype
0      1551         43         0      0    No Failure
1      1408         46         3      0    No Failure
2      1498         49         5      0    No Failure
3      1433         40         7      0    No Failure
4      1408         40         9      0    No Failure

我们还应该对分类特征进行一些频率分析。绝大多数观测值，97%，没有出现故障。这种相当明显的类别不平衡可能很难建模。有三种机器类型——高质量、低质量和中等质量：

machinefailuretype.failtype.value_counts(dropna=False).sort_index()
Heat Dissipation Failure  112
No Failure                9652
Overstrain Failure        78
Power Failure             95
Random Failures           18
Tool Wear Failure         45
Name: failtype, dtype: int64
machinefailuretype.machinetype.\
  value_counts(dropna=False).sort_index()
H      1003
L      6000
M      2997
Name: machinetype, dtype: int64

让我们合并一些failtype值并检查我们的工作。首先，我们将定义一个函数setcode，将故障类型文本映射到故障类型代码。我们将随机分配故障和工具磨损故障到代码5，用于其他故障：

def setcode(typetext):
  if (typetext=="No Failure"):
    typecode = 1
  elif (typetext=="Heat Dissipation Failure"):
    typecode = 2
  elif (typetext=="Power Failure"):
    typecode = 3
  elif (typetext=="Overstrain Failure"):
    typecode = 4
  else:
    typecode = 5
  return typecode
machinefailuretype["failtypecode"] = \
  machinefailuretype.apply(lambda x: setcode(x.failtype), axis=1)
machinefailuretype.groupby(['failtypecode','failtype']).size().\
  reset_index()
   failtypecode  failtype                   0
0   1            No Failure                 9652
1   2            Heat Dissipation Failure   112
2   3            Power Failure              95
3   4            Overstrain Failure         78
4   5            Random Failures            18
5   5            Tool Wear Failure          45

我们应该查看我们数值特征的描述性统计：

num_cols = ['airtemp', 'processtemperature',
  'rotationalspeed', 'torque', 'toolwear']
cat_cols = ['machinetype']
machinefailuretype[num_cols].agg(['min','median','max']).T
                       min       median       max
airtemp                295       300          304
processtemperature     306       310          314
rotationalspeed        1,168     1,503        2,886
torque                 4         40           77
toolwear               0         108          253

现在，我们准备创建训练和测试数据框。我们将使用我们刚刚创建的故障类型代码作为我们的目标：

X_train, X_test, y_train, y_test =  \
  train_test_split(machinefailuretype[num_cols + cat_cols],\
  machinefailuretype[['failtypecode']], test_size=0.2, random_state=0)

现在，让我们设置列转换。对于数值特征，我们将将异常值设置为中位数，然后缩放数据。我们将使用最小-最大缩放，这将返回从 0 到 1 的值（MinMaxScaler的默认值）。我们使用这个缩放器，而不是标准缩放器，以避免负值。我们稍后将使用的特征选择方法selectkbest不能与负值一起使用：
```
ohe = OneHotEncoder(drop='first', sparse=False)
cattrans = make_pipeline(ohe)
standtrans = make_pipeline(
  OutlierTrans(3),SimpleImputer(strategy="median"),
  MinMaxScaler())
coltrans = ColumnTransformer(
  transformers=[
    ("cat", cattrans, cat_cols),
    ("stand", standtrans, num_cols),
  ]
)
```

让我们也看看编码后的列。我们需要在过度采样之前做这件事，因为SMOTENC模块需要分类特征的列索引。我们进行过度采样是为了处理显著的类别不平衡。我们已在第十一章中更详细地讨论了这一点，决策树和随机森林分类：

coltrans.fit(X_train.sample(1000))
new_cat_cols = \
  coltrans.\
  named_transformers_['cat'].\
  named_steps['onehotencoder'].\
  get_feature_names(cat_cols)
new_cols = np.concatenate((new_cat_cols, np.array(num_cols)))
print(new_cols)
['machinetype_L' 'machinetype_M' 'airtemp' 
'processtemperature' 'rotationalspeed' 'torque' 
'toolwear']

接下来，我们将为我们的模型设置一个管道。该管道将执行列转换，使用SMOTENC进行过度采样，使用selectkbest进行特征选择，然后运行 KNN 模型。记住，我们必须将分类特征的列索引传递给SMOTENC，以便它能够正确运行：
```
catcolscnt = new_cat_cols.shape[0]
smotenc = SMOTENC(categorical_features=np.arange(0,catcolscnt), random_state=0)
knn = KNeighborsClassifier(n_jobs=-1)
pipe1 = make_pipeline(coltrans, smotenc, SelectKBest(score_func=chi2), knn)
```

现在，我们准备拟合我们的模型。我们将进行随机网格搜索，以确定 KNN 的最佳值和距离度量。我们还将搜索特征选择的最佳k值：

knn_params = {
 'selectkbest__k': np.arange(1, len(new_cols)),
 'kneighborsclassifier__n_neighbors': np.arange(5, 175, 2),
 'kneighborsclassifier__metric': ['euclidean','manhattan','minkowski']
}
rs = RandomizedSearchCV(pipe1, knn_params, cv=5, scoring="roc_auc_ovr_weighted")
rs.fit(X_train, y_train.values.ravel())

让我们看看网格搜索发现了什么。除了processtemperature之外的所有特征都值得保留在模型中。KNN 的最佳值和距离度量分别是125和minkowski。基于 ROC 曲线下的面积的最佳分数是0.9：

selected = rs.best_estimator_['selectkbest'].get_support()
selected.sum()
6
new_cols[selected]
array(['machinetype_L', 'machinetype_M', 'airtemp', 
       'rotationalspeed', 'torque', 'toolwear'], 
       dtype=object)
rs.best_params_
{'selectkbest__k': 6,
 'kneighborsclassifier__n_neighbors': 125,
 'kneighborsclassifier__metric': 'minkowski'}
rs.best_score_
0.899426752716227

让我们看看一个混淆矩阵。查看第一行，我们可以看到当没有发生故障时，发现了相当数量的故障。然而，我们的模型正确地识别了大多数实际的热量、功率和过载故障。这可能不是一个可怕的精确度和敏感度权衡。根据问题，我们可能接受大量假阳性，以在我们的模型中获得可接受的敏感度水平：
```
pred = rs.predict(X_test)
cm = skmet.confusion_matrix(y_test, pred)
cmplot = skmet.ConfusionMatrixDisplay(confusion_matrix=cm,
   display_labels=['None', 'Heat','Power','Overstrain','Other'])
cmplot.plot()
cmplot.ax_.set(title='Machine Failure Type Confusion Matrix', 
  xlabel='Predicted Value', ylabel='Actual Value')
```

这产生了以下图表：

图 12.4 – 超参数调整后机器故障类型的混淆矩阵

我们还应该查看一个分类报告。你可能还记得从第六章，为模型评估做准备，宏平均是简单地在类别间取平均。在这里，我们更感兴趣的是加权平均。加权 F1 分数为0.81并不差。记住，F1 是精确率和敏感度的调和平均数：

print(skmet.classification_report(y_test, pred,
  target_names=\
  ['None', 'Heat','Power','Overstrain','Other']))
              Precision  recall  f1-score  support
        None  0.99       0.71    0.83      1927
        Heat  0.11       0.90    0.20      21
       Power  0.15       0.61    0.24      18
  Overstrain  0.36       0.76    0.49      21
       Other  0.01       0.31    0.02      13
    accuracy                     0.71      2000
   macro avg  0.33       0.66    0.36      2000
weighted avg  0.96       0.71    0.81      2000

机器故障类型的类别不平衡使得建模特别困难。尽管如此，我们的 KNN 模型表现相对较好，假设大量假阳性不是问题。在这种情况下，一个假阳性可能不像一个假阴性那样成问题。它可能只是需要对看似有故障风险的机器进行更多检查。如果我们将其与实际机器故障的惊讶相比，偏向于敏感度而不是精确度可能是合适的。

让我们在另一个多类问题上尝试 KNN。

字母识别的 KNN

我们可以采取与预测机器故障时使用字母识别相同的策略。只要我们有能够很好地区分字母的特征，KNN 就是该模型的合理选择。我们将在本节尝试这种方法。

注意

在本节中，我们将使用字母识别数据。这些数据可在archive-beta.ics.uci.edu/ml/datasets/letter+recognition公开使用。有 26 个字母（全部为大写）和 20 种不同的字体。16 个不同的特征捕捉每个字母的不同属性。

让我们构建模型：

首先，我们将加载我们已经使用过的相同库：

import pandas as pd
import numpy as np
from sklearn.model_selection import train_test_split
from sklearn.model_selection import StratifiedKFold, \
  GridSearchCV
from sklearn.neighbors import KNeighborsClassifier
import sklearn.metrics as skmet
from scipy.stats import randint

现在，我们将加载数据并查看前几个实例。有 20,000 个观测值和 17 列。letter是我们的目标：

letterrecognition = pd.read_csv("data/letterrecognition.csv")
letterrecognition.shape
(20000, 17)
letterrecognition.head().T
            0     1     2     3     4
letter      T     I     D     N     G
xbox        2     5     4     7     2
ybox        8     12    11    11    1
width       3     3     6     6     3
height      5     7     8     6     1
onpixels    1     2     6     3     1
xbar        8     10    10    5     8
ybar        13    5     6     9     6
x2bar       0     5     2     4     6
y2bar       6     4     6     6     6
xybar       6     13    10    4     6
x2ybar      10     3    3     4     5
xy2bar      8     9     7     10    9
x-ege       0     2     3     6     1
xegvy       8     8     7     10    7
y-ege       0     4     3     2     5
yegvx       8     10    9     8     10

现在，让我们创建训练和测试数据框：

X_train, X_test, y_train, y_test =  \
  train_test_split(letterrecognition.iloc[:,1:],\
  letterrecognition.iloc[:,0:1], test_size=0.2, 
  random_state=0)

接下来，让我们实例化一个 KNN 实例。我们还将设置分层 K 折交叉验证和超参数的字典。我们将寻找k（n_neighbors）和距离度量的最佳超参数：

knn = KNeighborsClassifier(n_jobs=-1)
kf = StratifiedKFold(n_splits=5, shuffle=True, random_state=0)
knn_params = {
  'n_neighbors': np.arange(3, 41, 2),
  'metric': ['euclidean','manhattan','minkowski']
}

现在，我们已经准备好进行彻底的网格搜索。在这里我们进行彻底搜索是因为我们没有很多超参数需要检查。表现最好的距离度量是欧几里得距离。最近邻的k值是3。这个模型使我们几乎达到 95%的准确率：
```
gs = GridSearchCV(knn, knn_params, cv=kf, scoring='accuracy')
gs.fit(X_train, y_train.values.ravel())
gs.best_params_
{'metric': 'euclidean', 'n_neighbors': 3}
gs.best_score_
0.9470625
```

让我们生成预测并绘制一个混淆矩阵：

pred = gs.best_estimator_.predict(X_test)
letters = np.sort(letterrecognition.letter.unique())
cm = skmet.confusion_matrix(y_test, pred)
cmplot = \
  skmet.ConfusionMatrixDisplay(confusion_matrix=cm,
  display_labels=letters)
cmplot.plot()
cmplot.ax_.set(title='Letters', 
  xlabel='Predicted Value', ylabel='Actual Value')

这会产生以下图表：

图 12.5 – 字母预测的混淆矩阵

让我们快速总结本章所学的内容。

摘要

本章展示了使用 KNN 进行二分类或多分类分类是多么容易。由于 KNN 不对正态性或线性做出假设，因此它可以在逻辑回归可能不会产生最佳结果的情况下使用。这种灵活性确实带来了过拟合的真实风险，因此在选择k时必须谨慎。我们还在本章探讨了如何调整二分类和多分类模型的超参数。最后，当我们在乎预测速度或处理大型数据集时，KNN 并不是一个很好的选择。在上一章中我们探讨的决策树或随机森林分类，在这些情况下通常是一个更好的选择。

另一个非常好的选择是支持向量分类。我们将在下一章探讨支持向量分类。

第十三章：第十三章：支持向量机分类

支持向量分类模型和 k 最近邻模型之间有一些相似之处。它们都是直观且灵活的。然而，由于算法的性质，支持向量分类比 k 最近邻有更好的可扩展性。与逻辑回归不同，它可以很容易地处理非线性模型。使用支持向量机进行分类的策略和问题与我们讨论的相似，见第八章，支持向量回归，当时我们使用支持向量机进行回归。

支持向量分类（SVC）的一个关键优势是它赋予我们减少模型复杂度的能力，同时不增加我们的特征空间。但它也提供了多个我们可以调整的杠杆，以限制过拟合的可能性。我们可以选择线性模型，或从几个非线性核中选择。我们可以使用正则化参数，就像我们在逻辑回归中所做的那样。通过扩展，我们还可以使用这些相同的技巧来构建多类模型。

在本章中，我们将探讨以下主题：

SVC 的关键概念
线性 SVC 模型
非线性 SVM 分类模型
多类分类的 SVM

技术要求

在本章中，我们将坚持使用 pandas、NumPy 和 scikit-learn 库。本章中的所有代码都使用 scikit-learn 版本 0.24.2 和 1.0.2 进行了测试。显示决策边界的代码需要 scikit-learn 版本 1.1.1 或更高版本。

SVC 的关键概念

我们可以使用支持向量机（SVMs）来找到一条线或曲线，通过类别来分离实例。当类别可以通过一条线来区分时，它们被称为线性可分。

然而，正如我们在图 13.1中可以看到的，可能有许多可能的线性分类器。每条线都成功地使用两个特征 x1 和 x2 来区分由点表示的两个类别，关键的区别在于这些线如何对新的实例进行分类，这些新实例由透明矩形表示。使用离正方形最近的线会将透明矩形分类为点。使用其他两条线中的任何一条都会将其分类为正方形。

图 13.1 – 三种可能的线性分类器

当线性判别器非常接近训练实例时，就像图 13.2中的两条线一样，就有更大的风险将新实例分类错误。我们希望得到一个能够给出最大类别间间隔的线；一个离每个类别的边界数据点最远的线。那就是图 13.1中的中间线，但在图 13.2中可以看得更清楚：

图 13.2 – SVM 分类和最大间隔

粗体线分割了最大间隔，被称为决策边界。每个类别的边界数据点被称为支持向量。

我们使用支持向量机来寻找类别之间具有最大间隔的线性判别式。它是通过找到一个可以最大化的间隔的方程来实现的，其中间隔是数据点到分离超平面的距离。在具有两个特征的情况下，如图 13.2所示，该超平面只是一条线。然而，这可以推广到具有更多维度的特征空间。

对于像图 13.2中的数据点，我们可以使用所谓的硬间隔分类而不会出现问题；也就是说，我们可以对每个类别的所有观察值在决策边界的正确一侧非常严格。但如果我们数据点的样子像图 13.3中的那些呢？在这里，有一个方形非常接近点。硬间隔分类器是左侧的线，给我们非常小的间隔。

图 13.3 – 带有硬间隔和软间隔的支持向量机

如果我们使用软间隔分类，则得到右侧的线。软间隔分类放宽了所有实例都必须正确分离的约束。正如图 13.3中的数据所示，允许训练数据中有少量错误分类可以给我们更大的间隔。我们忽略偏离的方形，得到由软间隔线表示的决策边界。

约束放宽的程度由C超参数决定。C的值越大，对间隔违规的惩罚就越大。不出所料，具有较大C值的模型更容易过拟合。图 13.4说明了间隔如何随着C值的改变而变化。在C = 1时，错误分类的惩罚较低，给我们比C为 100 时更大的间隔。然而，即使在C为 100 的情况下，仍然会发生一些间隔违规。

图 13.4 – 不同 C 值下的软间隔

在实际操作中，我们几乎总是用软间隔构建我们的 SVC 模型。scikit-learn 中C的默认值是 1。

非线性支持向量机和核技巧

我们尚未完全解决 SVC 的线性可分性问题。为了简单起见，回到一个涉及两个特征的分类问题是有帮助的。假设两个特征与分类目标的关系图看起来像图 13.5中的插图。目标有两个可能的值，由点和正方形表示。x1 和 x2 是数值，具有负值。

图 13.5 – 使用两个特征无法线性分离的类别标签

在这种情况下，我们如何识别类之间的边界？通常情况下，在更高的维度中可以识别出边界。在这个例子中，我们可以使用图 13.6 中所示的多项式变换：

图 13.6 – 使用多项式变换建立边界

现在有一个第三维度，它是 x1 和 x2 平方和的总和。点都高于平方。这与我们使用多项式变换进行线性回归的方式相似。

这种方法的缺点之一是我们可能会迅速拥有太多特征，以至于模型无法良好地执行。这就是核技巧大显身手的地方。SVC 可以使用核函数隐式地扩展特征空间，而不实际创建更多特征。这是通过创建一个可以用来拟合非线性边界的值向量来实现的。

虽然这允许我们拟合一个类似于图 13.6 中假设的假设多项式变换，但 SVC 中最常用的核函数是径向基函数（RBF）。RBF 之所以受欢迎，是因为它比其他常见的核函数更快，并且可以使用伽马超参数提供额外的灵活性。RBF 核的方程如下：

在这里，和是数据点。伽马值，，决定了每个点的影响力大小。当伽马值较高时，点必须非常接近才能被分组在一起。在伽马值非常高的情况下，我们开始看到点的岛屿。

当然，伽马值或C的高值取决于我们的数据。一个好的方法是，在大量建模之前，创建不同伽马值和C值的决策边界可视化。这将让我们了解在不同的超参数值下，我们是否过度拟合或欠拟合。在本章中，我们将绘制不同伽马值和C值的决策边界。

SVC 的多类分类

到目前为止，我们关于支持向量机（SVC）的所有讨论都集中在二元分类上。幸运的是，适用于二元分类支持向量机的所有关键概念也适用于我们的目标值超过两个的可能值时的分类。我们将多类问题建模为一对一或一对余问题，从而将其转化为二元分类问题。

图 13.7 – 多类 SVC 选项

在三类示例中，一对一分类很容易说明，如图 13.7 的左侧所示。每个类别与每个其他类别之间估计一个决策边界。例如，虚线是点类与正方形类之间的决策边界。实线是点与椭圆形之间的决策边界。

在一对一分类中，每个类别与不属于该类别的实例之间构建一个决策边界。这如图 13.7 的右侧所示。实线是点与不是点（即正方形或椭圆形）的实例之间的决策边界。虚线和双线分别是正方形与剩余实例和椭圆形与剩余实例之间的决策边界。

我们可以使用一对一或一对一分类来构建线性和非线性 SVC 模型。我们还可以指定C的值来构建软边界。然而，使用这些技术中的每一个构建更多的决策边界需要比二分类 SVC 更多的计算资源。如果我们有大量的观察结果、许多特征和多个参数需要调整，我们可能需要非常好的系统资源才能及时获得结果。

三类示例隐藏了一个关于一对一和一对一分类器不同之处。对于三个类别，它们使用相同数量的分类器（三个），但随着一对一的增加，分类器的数量相对迅速地增加。在一对一分类中，分类器的数量始终等于类别值的数量，而在一对一分类中，它等于以下：

在这里，S是分类器的数量，N是目标的目标值（类别值）的基数。因此，基数是 4 时，一对一分类需要 4 个分类器，而一对一分类使用 6 个。

我们在本章的最后部分探讨了多类 SVC 模型，但让我们从一个相对简单的线性模型开始，看看 SVC 的实际应用。在为 SVC 模型进行预处理时，有两个需要注意的事项。首先，SVC 对特征的规模很敏感，因此在我们拟合模型之前需要解决这一点。其次，如果我们使用硬边界或高C值，异常值可能会对我们的模型产生很大的影响。

线性 SVC 模型

我们通常可以通过使用线性 SVC 模型获得良好的结果。当我们有超过两个特征时，没有简单的方法来可视化我们的数据是否线性可分。我们通常根据超参数调整来决定是线性还是非线性。在本节中，我们将假设我们可以通过线性模型和软边界获得良好的性能。

在本节中，我们将处理关于美国职业篮球联赛（NBA）比赛的数据。数据集包含了从 2017/2018 赛季到 2020/2021 赛季每场 NBA 比赛的统计数据。这包括主队、主队是否获胜、客队、客队和主队的投篮命中率、失误、篮板和助攻，以及许多其他指标。

注意

NBA 比赛数据可在www.kaggle.com/datasets/wyattowalsh/basketball供公众下载。该数据集从 1946/1947 赛季的 NBA 赛季开始。它使用nba_api从nba.com获取统计数据。该 API 可在github.com/swar/nba_api找到。

让我们构建一个线性 SVC 模型：

我们首先加载熟悉的库。唯一的新模块是LinearSVC和DecisionBoundaryDisplay。我们将使用DecisionBoundaryDisplay来显示线性模型的边界：

import pandas as pd
import numpy as np
from sklearn.model_selection import train_test_split
from sklearn.preprocessing import OneHotEncoder, StandardScaler
from sklearn.svm import LinearSVC
from scipy.stats import uniform
from sklearn.impute import SimpleImputer
from sklearn.pipeline import make_pipeline
from sklearn.compose import ColumnTransformer
from sklearn.feature_selection import RFECV
from sklearn.inspection import DecisionBoundaryDisplay
from sklearn.model_selection import cross_validate, \
  RandomizedSearchCV, RepeatedStratifiedKFold
import sklearn.metrics as skmet
import seaborn as sns
import os
import sys
sys.path.append(os.getcwd() + "/helperfunctions")
from preprocfunc import OutlierTrans

我们已经准备好加载 NBA 比赛数据。我们只需做一些清理工作。少数观测值的目标WL_HOME（主队是否获胜）有缺失值。我们移除这些观测值。我们将WL_HOME特征转换为0和1特征。

在这里，类别不平衡的问题并不大。这将在以后为我们节省一些时间：

nbagames = pd.read_csv("data/nbagames2017plus.csv", parse_dates=['GAME_DATE'])
nbagames = \
  nbagames.loc[nbagames.WL_HOME.isin(['W','L'])]
nbagames.shape
(4568, 149)
nbagames['WL_HOME'] = \
  np.where(nbagames.WL_HOME=='L',0,1).astype('int')
nbagames.WL_HOME.value_counts(dropna=False)
1    2586
0    1982
Name: WL_HOME, dtype: int64

让我们按数据类型组织我们的特征：

num_cols = ['FG_PCT_HOME','FTA_HOME','FG3_PCT_HOME',
  'FTM_HOME','FT_PCT_HOME','OREB_HOME','DREB_HOME',
  'REB_HOME','AST_HOME','STL_HOME','BLK_HOME',
  'TOV_HOME','FG_PCT_AWAY','FTA_AWAY','FG3_PCT_AWAY',
  'FT_PCT_AWAY','OREB_AWAY','DREB_AWAY','REB_AWAY',
  'AST_AWAY','STL_AWAY','BLK_AWAY','TOV_AWAY']
cat_cols = ['SEASON']

让我们看看一些描述性统计。（为了节省空间，我已经从打印输出中省略了一些特征。）我们需要缩放这些特征，因为它们的范围差异很大。没有缺失值，但当我们为极端值分配缺失值时，我们将生成一些缺失值：

nbagames[['WL_HOME'] + num_cols].agg(['count','min','median','max']).T
                  count     min     median  max
WL_HOME           4,568     0.00    1.00    1.00
FG_PCT_HOME       4,568     0.27    0.47    0.65
FTA_HOME          4,568     1.00    22.00   64.00
FG3_PCT_HOME      4,568     0.06    0.36    0.84
FTM_HOME          4,568     1.00    17.00   44.00
FT_PCT_HOME       4,568     0.14    0.78    1.00
OREB_HOME         4,568     1.00    10.00   25.00
DREB_HOME         4,568     18.00   35.00   55.00
REB_HOME          4,568     22.00   45.00   70.00
AST_HOME          4,568     10.00   24.00   50.00
.........
FT_PCT_AWAY       4,568     0.26    0.78    1.00
OREB_AWAY         4,568     0.00    10.00   26.00
DREB_AWAY         4,568     18.00   34.00   56.00
REB_AWAY          4,568     22.00   44.00   71.00
AST_AWAY          4,568     9.00    24.00   46.00
STL_AWAY          4,568     0.00    8.00    19.00
BLK_AWAY          4,568     0.00    5.00    15.00
TOV_AWAY          4,568     3.00    14.00   30.00

我们还应该回顾一下特征的相关性：

corrmatrix = nbagames[['WL_HOME'] + \
  num_cols].corr(method="pearson")
sns.heatmap(corrmatrix, 
  xticklabels=corrmatrix.columns,
  yticklabels=corrmatrix.columns, cmap="coolwarm")
plt.title('Heat Map of Correlation Matrix')
plt.tight_layout()
plt.show()

这产生了以下图表：

图 13.8 – NBA 比赛统计数据相关性热图

一些特征与目标相关，包括主队的投篮命中率（FG_PCT_HOME）和主队的防守篮板球（DREB_HOME）。

特征之间也存在相关性。例如，主队的投篮命中率（FG_PCT_HOME）和主队的 3 分投篮命中率（FG3_PCT_HOME）呈正相关，这并不令人意外。此外，主队的篮板球（REB_HOME）和防守篮板球（DREB_HOME）可能过于紧密地相关，以至于任何模型都无法分离它们的影响。

接下来，我们创建训练和测试数据框：

X_train, X_test, y_train, y_test =  \
  train_test_split(nbagames[num_cols + cat_cols],\
  nbagames[['WL_HOME']], test_size=0.2, random_state=0)

我们需要设置列转换。对于数值列，我们检查异常值并缩放数据。我们将一个分类特征SEASON进行独热编码。我们将在网格搜索中使用这些转换：

ohe = OneHotEncoder(drop='first', sparse=False)
cattrans = make_pipeline(ohe)
standtrans = make_pipeline(OutlierTrans(2),
  SimpleImputer(strategy="median"), StandardScaler())
coltrans = ColumnTransformer(
  transformers=[
    ("cat", cattrans, cat_cols),
    ("stand", standtrans, num_cols)
  ]
)

在构建我们的模型之前，让我们看看一个线性 SVC 模型的决策边界。我们基于与目标相关的两个特征来设置边界：主队的投篮命中率(FG_PCT_HOME)和主队的防守篮板(DREB_HOME)。

我们创建了一个函数dispbound，它将使用DecisionBoundaryDisplay模块来显示边界。这个模块在 scikit-learn 版本 1.1.1 或更高版本中可用。DecisionBoundaryDisplay需要一个模型来拟合，两个特征和目标值：

pipe0 = make_pipeline(OutlierTrans(2),
  SimpleImputer(strategy="median"), StandardScaler())
X_train_enc = pipe0.\
  fit_transform(X_train[['FG_PCT_HOME','DREB_HOME']])
def dispbound(model, X, xvarnames, y, title):
  dispfit = model.fit(X,y)
  disp = DecisionBoundaryDisplay.from_estimator(
    dispfit, X, response_method="predict",
    xlabel=xvarnames[0], ylabel=xvarnames[1],
    alpha=0.5,
  )
  scatter = disp.ax_.scatter(X[:,0], X[:,1],
    c=y, edgecolor="k")

  disp.ax_.set_title(title)
  legend1 = disp.ax_.legend(*scatter.legend_elements(),
    loc="lower left", title="Home Win")
  disp.ax_.add_artist(legend1)
dispbound(LinearSVC(max_iter=1000000,loss='hinge'),
  X_train_enc, ['FG_PCT_HOME','DREB_HOME'],
  y_train.values.ravel(),
  'Linear SVC Decision Boundary')

这产生了以下图表：

图 13.9 – 双特征线性 SVC 模型的决策边界

我们只使用两个特征就得到了一个相当不错的线性边界。这是个好消息，但让我们构建一个更精心设计的模型。

为了构建我们的模型，我们首先实例化一个线性 SVC 对象并设置递归特征消除。然后我们将列转换、特征选择和线性 SVC 添加到管道中并拟合它：

svc = LinearSVC(max_iter=1000000, loss='hinge',
   random_state=0)
rfecv = RFECV(estimator=svc, cv=5)
pipe1 = make_pipeline(coltrans, rfecv, svc)
pipe1.fit(X_train, y_train.values.ravel())

让我们看看从我们的递归特征消除中选择了哪些特征。我们首先需要获取一元编码后的列名。然后我们可以使用rfecv对象的get_support方法来获取选定的特征。（如果你使用的是 scikit-learn 版本 1 或更高版本，你会得到一个关于get_feature_names的弃用警告。不过，你可以使用get_feature_names_out代替，尽管这不会与 scikit-learn 的早期版本兼容。）

new_cat_cols = \
  pipe1.named_steps['columntransformer'].\
  named_transformers_['cat'].\
  named_steps['onehotencoder'].\
  get_feature_names(cat_cols)
new_cols = np.concatenate((new_cat_cols, np.array(num_cols)))
sel_cols = new_cols[pipe1['rfecv'].get_support()]
np.set_printoptions(linewidth=55)
sel_cols
array(['SEASON_2018', 'SEASON_2019', 'SEASON_2020',
       'FG_PCT_HOME', 'FTA_HOME', 'FG3_PCT_HOME',
       'FTM_HOME', 'FT_PCT_HOME', 'OREB_HOME',
       'DREB_HOME', 'REB_HOME', 'AST_HOME',
       'TOV_HOME', 'FG_PCT_AWAY', 'FTA_AWAY',
       'FG3_PCT_AWAY', 'FT_PCT_AWAY', 'OREB_AWAY',
       'DREB_AWAY', 'REB_AWAY', 'AST_AWAY',
       'BLK_AWAY', 'TOV_AWAY'], dtype=object)

我们应该看看系数。对于每个选定的列的系数可以通过linearsvc对象的coef_属性来访问。也许并不令人惊讶，主队的投篮命中率(FG_PCT_HOME)和客队的投篮命中率(FG_PCT_AWAY)是主队获胜的最重要正负预测因子。接下来最重要的特征是客队和主队的失误次数：

pd.Series(pipe1['linearsvc'].\
  coef_[0], index=sel_cols).\
  sort_values(ascending=False)
FG_PCT_HOME     2.21
TOV_AWAY        1.20
REB_HOME        1.19
FTM_HOME        0.95
FG3_PCT_HOME    0.94
FT_PCT_HOME     0.31
AST_HOME        0.25
OREB_HOME       0.18
DREB_AWAY       0.11
SEASON_2018     0.10
FTA_HOME       -0.05
BLK_AWAY       -0.07
SEASON_2019    -0.11
SEASON_2020    -0.19
AST_AWAY       -0.44
OREB_AWAY      -0.47
DREB_HOME      -0.49
FT_PCT_AWAY    -0.53
REB_AWAY       -0.63
FG3_PCT_AWAY   -0.80
FTA_AWAY       -0.81
TOV_HOME       -1.19
FG_PCT_AWAY    -1.91
dtype: float64

让我们看看预测结果。我们的模型在预测主队获胜方面做得很好：

pred = pipe1.predict(X_test)
print("accuracy: %.2f, sensitivity: %.2f, specificity: %.2f, precision: %.2f"  %
  (skmet.accuracy_score(y_test.values.ravel(), pred),
  skmet.recall_score(y_test.values.ravel(), pred),
  skmet.recall_score(y_test.values.ravel(), pred, pos_label=0),
  skmet.precision_score(y_test.values.ravel(), pred)))
accuracy: 0.93, sensitivity: 0.95, specificity: 0.92, precision: 0.93

我们应该通过进行交叉验证来确认这些指标不是偶然的。我们使用重复的有放回分层 k 折来进行验证，这意味着我们想要 7 个折和 10 次迭代。我们得到的结果与之前步骤中的结果几乎相同：

kf = RepeatedStratifiedKFold(n_splits=7,n_repeats=10,\
  random_state=0)
scores = cross_validate(pipe1, X_train, \
  y_train.values.ravel(), \
  scoring=['accuracy','precision','recall','f1'], \
  cv=kf, n_jobs=-1)
print("accuracy: %.2f, precision: %.2f, sensitivity: %.2f, f1: %.2f"  %
  (np.mean(scores['test_accuracy']),\
  np.mean(scores['test_precision']),\
  np.mean(scores['test_recall']),\
  np.mean(scores['test_f1'])))
accuracy: 0.93, precision: 0.93, sensitivity: 0.95, f1: 0.94

到目前为止，我们一直在使用C的默认值1。我们可以尝试使用随机网格搜索来识别一个更好的C值：

svc_params = {
 'linearsvc__C': uniform(loc=0, scale=100)
}
rs = RandomizedSearchCV(pipe1, svc_params, cv=10, 
  scoring='accuracy', n_iter=20, random_state=0)
rs.fit(X_train, y_train.values.ravel())
rs.best_params_
{'linearsvc__C': 54.88135039273247}
rs.best_score_
0.9315809566584325

最佳的C值是 2.02，最佳的准确度得分是 0.9316。

让我们仔细看看 20 次网格搜索中每次的得分。每个得分是 10 个折的准确度得分的平均值。实际上，无论C值如何，我们得到的分数都差不多：

results = \
  pd.DataFrame(rs.cv_results_['mean_test_score'], \
    columns=['meanscore']).\
  join(pd.DataFrame(rs.cv_results_['params'])).\
  sort_values(['meanscore'], ascending=False)
results
              meanscore        linearsvc__C
0             0.93             54.88
8             0.93             96.37
18            0.93             77.82
17            0.93             83.26
13            0.93             92.56
12            0.93             56.80
11            0.93             52.89
1             0.93             71.52
10            0.93             79.17
7             0.93             89.18
6             0.93             43.76
5             0.93             64.59
3             0.93             54.49
2             0.93             60.28
19            0.93             87.00
9             0.93             38.34
4             0.93             42.37
14            0.93             7.10
15            0.93             8.71
16            0.93             2.02

让我们现在看看一些预测结果。我们的模型在各个方面都做得很好，但并没有比初始模型做得更好：

pred = rs.predict(X_test)
print("accuracy: %.2f, sensitivity: %.2f, specificity: %.2f, precision: %.2f"  %
  (skmet.accuracy_score(y_test.values.ravel(), pred),
  skmet.recall_score(y_test.values.ravel(), pred),
  skmet.recall_score(y_test.values.ravel(), pred, pos_label=0),
  skmet.precision_score(y_test.values.ravel(), pred)))
accuracy: 0.93, sensitivity: 0.95, specificity: 0.92, precision: 0.93

让我们也看看一个混淆矩阵：

cm = skmet.confusion_matrix(y_test, pred)
cmplot = \
  skmet.ConfusionMatrixDisplay(confusion_matrix=cm, display_labels=['Loss', 'Won'])
cmplot.plot()
cmplot.ax_.set(title='Home Team Win Confusion Matrix', 
  xlabel='Predicted Value', ylabel='Actual Value')

这产生了以下图表：

图 13.10 – 主队胜负的混淆矩阵

我们的模型在很大程度上正确预测了主队的胜负。调整C的值并没有带来太大的变化，因为我们几乎无论C值如何都能获得相同的准确率。

注意

你可能已经注意到，我们在处理 NBA 比赛数据时比在之前章节中处理的心脏病和机器故障数据时更频繁地使用准确率指标。我们更关注那个数据的敏感性。这有两个原因。首先，当类别几乎平衡时，准确率是一个更有说服力的度量标准，正如我们在第六章“准备模型评估”中详细讨论的那样。其次，在预测心脏病和机器功率故障时，我们倾向于敏感性，因为那些领域中的假阴性成本高于假阳性。而对于预测 NBA 比赛，则没有这样的偏见。

线性 SVC 模型的一个优点是它们很容易解释。我们能够查看系数，这有助于我们理解模型并与其他人沟通我们预测的基础。尽管如此，确认我们不会使用非线性模型获得更好的结果也是有帮助的。我们将在下一节中这样做。

非线性 SVM 分类模型

虽然非线性 SVC 在概念上比线性 SVC 更复杂，正如我们在本章第一节中看到的，使用 scikit-learn 运行非线性模型相对简单。与线性模型的主要区别是我们需要进行相当多的超参数调整。我们必须指定C、gamma的值以及我们想要使用的核函数。

虽然有一些理论上的理由可以假设对于特定的建模挑战，某些超参数值可能比其他值更有效，但我们通常通过经验方法（即超参数调整）来解决这个问题。我们将在本节中使用与上一节相同的 NBA 比赛数据来尝试这样做：

我们加载了上一节中使用的相同库。我们还导入了LogisticRegression模块。我们稍后将会使用该模块与特征选择包装器方法结合：

import pandas as pd
import numpy as np
from sklearn.preprocessing import MinMaxScaler
from sklearn.pipeline import make_pipeline
from sklearn.svm import SVC
from sklearn.linear_model import LogisticRegression
from scipy.stats import uniform
from sklearn.feature_selection import RFECV
from sklearn.impute import SimpleImputer
from scipy.stats import randint
from sklearn.model_selection import RandomizedSearchCV
import sklearn.metrics as skmet
import os
import sys
sys.path.append(os.getcwd() + "/helperfunctions")
from preprocfunc import OutlierTrans

我们导入了nbagames模块，其中包含加载和预处理 NBA 比赛数据的代码。这仅仅是我们在上一节中运行以准备建模数据的代码的副本。在这里没有必要重复那些步骤。

我们还导入了上一节中使用的dispbound函数来显示决策边界。我们将那段代码复制到了当前目录下的helperfunctions子目录中，文件名为displayfunc.py：

import nbagames as ng
from displayfunc import dispbound

我们使用nbagames模块来获取训练和测试数据：

X_train = ng.X_train
X_test = ng.X_test
y_train = ng.y_train
y_test = ng.y_test

在构建模型之前，让我们看看具有两个特征（主队的投篮命中率 FG_PCT_HOME 和主队的防守篮板 DREB_HOME）的几个不同核的决策边界。我们首先使用 rbf 核，并使用不同的 gamma 和 C 值：

pipe0 = make_pipeline(OutlierTrans(2),
  SimpleImputer(strategy="median"),
  StandardScaler())
X_train_enc = \
  pipe0.fit_transform(X_train[['FG_PCT_HOME',
   'DREB_HOME']])
dispbound(SVC(kernel='rbf', gamma=30, C=1),
  X_train_enc,['FG_PCT_HOME','DREB_HOME'],
  y_train.values.ravel(),
  "SVC with rbf kernel-gamma=30, C=1")

以几种不同的方式运行此操作会产生以下图表：

图 13.11 – 使用 rbf 核和不同 gamma 和 C 值的决策边界

在 gamma 和 C 的值接近默认值时，我们看到决策边界有一些弯曲，以适应损失类中的几个偏离的点。这些是主队尽管有很高的防守篮板总数却输掉比赛的情况。使用 rbf 核，其中两个这样的实例现在被正确分类。还有一些主队投篮命中率很高但防守篮板很低的实例，现在也被正确分类。然而，与上一节中的线性模型相比，我们的预测整体上并没有太大变化。

但如果我们增加 C 或 gamma 的值，这种变化会显著。回想一下，C 的较高值会增加误分类的惩罚。这导致边界围绕实例旋转。

将 gamma 增加到 30 会导致严重的过度拟合。gamma 的高值意味着数据点必须非常接近才能被分组在一起。这导致决策边界紧密地与少数实例相关联，有时甚至只有一个实例。

我们还可以展示多项式核的边界。我们将保持默认的 C 值，以关注改变度数的影响：

dispbound(SVC(kernel='poly', degree=7),
  X_train_enc, ['FG_PCT_HOME','DREB_HOME'],
  y_train.values.ravel(),
  "SVC with polynomial kernel - degree=7")

以几种不同的方式运行此操作会产生以下图表：

图 13.12 – 使用多项式核和不同度数的决策边界

我们可以看到在较高度数级别上决策边界的某些弯曲，以处理几个不寻常的实例。这里并没有过度拟合，但我们的预测也没有真正得到很大改善。

这至少暗示了当我们构建模型时可以期待什么。我们应该尝试一些非线性模型，但有很大可能性它们不会比我们在上一节中使用的线性模型带来太多改进。

现在，我们已准备好设置我们将用于非线性 SVC 的管道。我们的管道将执行列转换和递归特征消除。我们使用逻辑回归进行特征选择：
```
rfecv = RFECV(estimator=LogisticRegression())
svc = SVC()
pipe1 = make_pipeline(ng.coltrans, rfecv, svc)
```
我们创建一个字典用于我们的超参数调整。这个字典的结构与我们用于此目的的其他字典略有不同。这是因为某些超参数只能与某些其他超参数一起使用。例如，gamma不能与线性核一起使用：
```
svc_params = [
  {
    'svc__kernel': ['rbf'],
    'svc__C': uniform(loc=0, scale=20),
    'svc__gamma': uniform(loc=0, scale=100)
  },
  {
    'svc__kernel': ['poly'],
    'svc__degree': randint(1, 5),
    'svc__C': uniform(loc=0, scale=20),
    'svc__gamma': uniform(loc=0, scale=100)
  },
  {
    'svc__kernel': ['linear','sigmoid'],
    'svc__C': uniform(loc=0, scale=20)
  }
]
```
注意

你可能已经注意到我们将使用的一个核是线性的，并想知道这与我们在上一节中使用的线性 SVC 模块有何不同。LinearSVC通常会更快地收敛，尤其是在大型数据集上。它不使用核技巧。我们可能也会得到不同的结果，因为优化在几个方面都不同。

现在我们已经准备好拟合一个 SVC 模型。最佳模型实际上是一个线性核的模型：

rs = RandomizedSearchCV(pipe1, svc_params, cv=5, 
  scoring='accuracy', n_iter=10, n_jobs=-1,
  verbose=5, random_state=0)
rs.fit(X_train, y_train.values.ravel())
rs.best_params_
{'svc__C': 1.1342595463488636, 'svc__kernel': 'linear'}
rs.best_score_
0.9299405955437289

让我们更仔细地看看选定的超参数和相关的准确率分数。我们可以从网格对象的cv_results_字典中获取params列表中的 20 个随机选择的超参数组合。我们也可以从同一个字典中获取平均测试分数。

我们按准确率分数降序排序。线性核优于多项式核和rbf核，尽管在3、4和5度上并不比多项式核显著更好。rbf核的表现尤其糟糕：

results = \
  pd.DataFrame(rs.cv_results_['mean_test_score'], \
    columns=['meanscore']).\
  join(pd.json_normalize(rs.cv_results_['params'])).\
  sort_values(['meanscore'], ascending=False)
results
            C          gamma     kernel     degree
meanscore                              
0.93        1.13       NaN       linear     NaN
0.89        1.42       64.82     poly       3.00
0.89        9.55       NaN       sigmoid    NaN
0.89        11.36      NaN       sigmoid    NaN
0.89        2.87       75.86     poly       5.00
0.64        12.47      43.76     poly       4.00
0.64        15.61      72.06     poly       4.00
0.57        11.86      84.43     rbf        NaN
0.57        16.65      77.82     rbf        NaN
0.57        19.57      79.92     rbf        NaN

注意

我们使用 pandas 的json_normalize方法来处理我们从params列表中提取的有些混乱的超参数组合。这是因为不同的超参数取决于所使用的核。这意味着params列表中的 20 个字典将具有不同的键。例如，多项式核将具有度数的值。线性核和rbf核则没有。

我们可以通过best_estimator_属性访问支持向量。有 625 个支持向量支撑着决策边界：
```
rs.best_estimator_['svc'].\
  support_vectors_.shape
(625, 18)
```

最后，我们可以看一下预测结果。不出所料，我们没有比上一节中运行的线性 SVC 模型得到更好的结果。我说不出所料，因为最佳模型被发现是一个线性核的模型：

pred = rs.predict(X_test)
print("accuracy: %.2f, sensitivity: %.2f, specificity: %.2f, precision: %.2f"  %
  (skmet.accuracy_score(y_test.values.ravel(), pred),
  skmet.recall_score(y_test.values.ravel(), pred),
  skmet.recall_score(y_test.values.ravel(), pred, 
    pos_label=0),
  skmet.precision_score(y_test.values.ravel(), pred)))
accuracy: 0.93, sensitivity: 0.94, specificity: 0.91, precision: 0.93

虽然我们没有改进上一节中的模型，但尝试一些非线性模型仍然是一项值得的练习。事实上，我们通常就是这样发现我们是否有可以成功线性分离的数据。这通常很难可视化，所以我们依赖于超参数调整来告诉我们哪个核最适合我们的数据。

本节和上一节展示了使用 SVM 进行二类分类的关键技术。我们迄今为止所做的大部分内容也适用于多类分类。在下一节中，当我们的目标值超过两个时，我们将探讨 SVC 建模策略。

多类分类的 SVM

当我们进行多类分类时，所有我们在使用 SVC 进行二类分类时遇到的问题都适用。我们需要确定类别是否线性可分，如果不是，哪个核将产生最佳结果。正如本章第一节所讨论的，我们还需要决定这种分类是否最好建模为一对一或一对多。一对一找到将每个类别与其他每个类别分开的决策边界。一对多找到将每个类别与所有其他实例区分开的决策边界。我们在本节中尝试这两种方法。

我们将使用我们在前几章中使用过的机器故障数据。

注意

这个关于机器故障的数据集可以在www.kaggle.com/datasets/shivamb/machine-predictive-maintenance-classification上公开使用。有 10,000 个观测值，12 个特征，以及两个可能的目标。一个是二元的：机器故障或未故障。另一个是故障类型。这个数据集中的实例是合成的，由一个旨在模拟机器故障率和原因的过程生成。

让我们构建一个多类 SVC 模型：

我们首先加载本章中一直在使用的相同库：

import pandas as pd
from sklearn.model_selection import train_test_split
from sklearn.preprocessing import OneHotEncoder, MinMaxScaler
from sklearn.pipeline import make_pipeline
from sklearn.svm import SVC
from scipy.stats import uniform
from sklearn.impute import SimpleImputer
from sklearn.compose import ColumnTransformer
from sklearn.model_selection import RandomizedSearchCV
import sklearn.metrics as skmet
import os
import sys
sys.path.append(os.getcwd() + "/helperfunctions")
from preprocfunc import OutlierTrans

我们将加载机器故障类型数据集并查看其结构。这里有字符和数值数据的混合。没有缺失值：

machinefailuretype = pd.read_csv("data/machinefailuretype.csv")
machinefailuretype.info()
<class 'pandas.core.frame.DataFrame'>
RangeIndex: 10000 entries, 0 to 9999
Data columns (total 10 columns):
 #   Column                 Non-Null Count     Dtype
---  ------                 --------------     -----  
 0   udi                    10000 non-null     int64
 1   product                10000 non-null     object 
 2   machinetype            10000 non-null     object 
 3   airtemp                10000 non-null     float64
 4   processtemperature     10000 non-null     float64
 5   rotationalspeed        10000 non-null     int64
 6   torque                 10000 non-null     float64
 7   toolwear               10000 non-null     int64
 8   fail                   10000 non-null     int64
 9   failtype               10000 non-null     object 
dtypes: float64(3), int64(4), object(3)
memory usage: 781.4+ KB

让我们看看一些观测值：

machinefailuretype.head()
   udi product machinetype airtemp processtemperature\
0  1   M14860       M         298        309 
1  2   L47181       L         298        309 
2  3   L47182       L         298        308 
3  4   L47183       L         298        309 
4  5   L47184       L         298        309 
   rotationalspeed  torque  toolwear  fail  failtype  
0        1551         43        0       0   No Failure  
1        1408         46        3       0   No Failure  
2        1498         49        5       0   No Failure  
3        1433         40        7       0   No Failure  
4        1408         40        9       0   No Failure

让我们也看看目标值的分布。我们有显著的类别不平衡，所以我们需要以某种方式处理这个问题：

machinefailuretype.failtype.\
  value_counts(dropna=False).sort_index()
Heat Dissipation Failure     112
No Failure                   9652
Overstrain Failure           78
Power Failure                95
Random Failures              18
Tool Wear Failure            45
Name: failtype, dtype: int64

我们可以通过为故障类型创建一个数字代码来节省一些麻烦，我们将使用这个数字代码而不是字符值。我们不需要将其放入管道中，因为我们没有在转换中引入任何数据泄露：

def setcode(typetext):
  if (typetext=="No Failure"):
    typecode = 1
  elif (typetext=="Heat Dissipation Failure"):
    typecode = 2
  elif (typetext=="Power Failure"):
    typecode = 3
  elif (typetext=="Overstrain Failure"):
    typecode = 4
  else:
    typecode = 5
  return typecode
machinefailuretype["failtypecode"] = \
  machinefailuretype.apply(lambda x: setcode(x.failtype), axis=1)

我们还应该查看一些描述性统计。我们需要对特征进行缩放：

num_cols = ['airtemp','processtemperature',
  'rotationalspeed','torque','toolwear']
cat_cols = ['machinetype']
machinefailuretype[num_cols].agg(['min','median','max']).T
                       min        median     max
airtemp                295.30     300.10     304.50
processtemperature     305.70     310.10     313.80
rotationalspeed        1,168.00   1,503.00   2,886.00
torque                 3.80       40.10      76.60
toolwear               0.00       108.00     253.00

现在让我们创建训练和测试数据框。我们还应该使用stratify参数来确保训练和测试数据中目标值的分布均匀：

X_train, X_test, y_train, y_test =  \
  train_test_split(machinefailuretype[num_cols + cat_cols],\
  machinefailuretype[['failtypecode']],\
  stratify=machinefailuretype[['failtypecode']], \
  test_size=0.2, random_state=0)

我们设置了需要运行的列转换。对于数值列，我们将异常值设置为中位数，然后缩放值。我们对一个分类特征machinetype进行了一元编码。它有H、M和L值，分别代表高质量、中质量和低质量：

ohe = OneHotEncoder(drop='first', sparse=False)
cattrans = make_pipeline(ohe)
standtrans = make_pipeline(OutlierTrans(3),
  SimpleImputer(strategy="median"),
  MinMaxScaler())
coltrans = ColumnTransformer(
  transformers=[
    ("cat", cattrans, cat_cols),
    ("stand", standtrans, num_cols),
  ]
)

接下来，我们设置一个包含列转换和 SVC 实例的管道。我们将class_weight参数设置为balanced以处理类别不平衡。这会应用一个与目标类别频率成反比的权重：
```
svc = SVC(class_weight='balanced', probability=True)
pipe1 = make_pipeline(coltrans, svc)
```

在这种情况下，我们只有少量特征，所以我们不会担心特征选择。（我们可能仍然会担心高度相关的特征，但在这个数据集中这不是一个问题。）

我们创建了一个字典，包含用于网格搜索的超参数组合。这基本上与我们在上一节中使用的字典相同，只是我们添加了一个决策函数形状键。这将导致网格搜索尝试一对一（ovo）和一对多（ovr）分类：

svc_params = [
  {
    'svc__kernel': ['rbf'],
    'svc__C': uniform(loc=0, scale=20),
    'svc__gamma': uniform(loc=0, scale=100),
    'svc__decision_function_shape': ['ovr','ovo']
  },
  {
    'svc__kernel': ['poly'],
    'svc__degree': np.arange(0,6),
    'svc__C': uniform(loc=0, scale=20),
    'svc__gamma': uniform(loc=0, scale=100),
    'svc__decision_function_shape': ['ovr','ovo']
  },
  {
    'svc__kernel': ['linear','sigmoid'],
    'svc__C': uniform(loc=0, scale=20),
    'svc__decision_function_shape': ['ovr','ovo']
  }
]

现在我们已经准备好运行随机网格搜索。我们将基于 ROC 曲线下的面积来评分。最佳超参数包括一对一决策函数和rbf核：

rs = RandomizedSearchCV(pipe1, svc_params, cv=7, scoring="roc_auc_ovr", n_iter=10)
rs.fit(X_train, y_train.values.ravel())
rs.best_params_
{'svc__C': 5.609789456747942,
 'svc__decision_function_shape': 'ovo',
 'svc__gamma': 27.73459801111866,
 'svc__kernel': 'rbf'}
rs.best_score_
0.9187636814475847

让我们看看每次迭代的分数。除了我们在上一步中看到的最佳模型外，还有几个其他超参数组合的分数几乎一样高。使用线性核的一对多几乎与表现最好的模型一样好：

results = \
  pd.DataFrame(rs.cv_results_['mean_test_score'], \
    columns=['meanscore']).\
  join(pd.json_normalize(rs.cv_results_['params'])).\
  sort_values(['meanscore'], ascending=False)
results
meanscore  svc__C svc__decision_function_shape  svc__gamma svc__kernel
7     0.92     5.61     ovo     27.73     rbf
5     0.91     9.43     ovr     NaN       linear
3     0.91     5.40     ovr     NaN       linear
0     0.90     19.84    ovr     28.70     rbf
8     0.87     5.34     ovo     93.87     rbf
6     0.86     8.05     ovr     80.57     rbf
9     0.86     4.41     ovo     66.66     rbf
1     0.86     3.21     ovr     85.35     rbf
4     0.85     0.01     ovo     38.24     rbf
2     0.66     7.61     ovr     NaN       sigmoid

我们应该看一下混淆矩阵：

pred = rs.predict(X_test)
cm = skmet.confusion_matrix(y_test, pred)
cmplot = skmet.ConfusionMatrixDisplay(confusion_matrix=cm,
   display_labels=['None', 'Heat','Power','Overstrain','Other'])
cmplot.plot()
cmplot.ax_.set(title='Machine Failure Type Confusion Matrix', 
  xlabel='Predicted Value', ylabel='Actual Value')

这产生了以下图表：

图 13.13 – 预测机器故障类型的混淆矩阵

让我们也做一个分类报告。尽管我们的模型在预测热和过载故障方面做得相当不错，但对于大多数类别，我们并没有获得很高的敏感性分数：

print(skmet.classification_report(y_test, pred,
  target_names=['None', 'Heat','Power', 'Overstrain', 'Other']))
              precision    recall  f1-score   support
        None       0.99      0.97      0.98      1930
        Heat       0.50      0.91      0.65        22
       Power       0.60      0.47      0.53        19
  Overstrain       0.65      0.81      0.72        16
       Other       0.06      0.15      0.09        13
    accuracy                           0.96      2000
   macro avg       0.56      0.66      0.59      2000
weighted avg       0.97      0.96      0.96      2000

当建模目标，如具有高类别不平衡的机器故障类型时，我们通常更关心除了准确性之外的指标。这部分取决于我们的领域知识。避免假阴性可能比避免假阳性更重要。过早地对机器进行彻底检查肯定比过晚进行更好。

96%到 97%的加权精确度、召回率（敏感性）和 f1 分数并不能很好地反映我们模型的表现。它们主要反映了类别不平衡很大，以及预测没有机器故障非常容易的事实。远低于宏观平均值（这些只是类别间的简单平均值）表明，我们的模型在预测某些类型的机器故障方面存在困难。

这个例子说明了将 SVC 扩展到具有多于两个值的目标的模型相对容易。我们可以指定是否想要使用一对一或一对多分类。当类别数量超过三个时，一对一方法可能会更快，因为将训练更少的分类器。

摘要

在本章中，我们探讨了实现 SVC 的不同策略。我们使用了线性 SVC（不使用核），当我们的类别是线性可分时，它可以表现得非常好。然后我们检查了如何使用核技巧将 SVC 扩展到类别不可分的情况。最后，我们使用一对一和一对多分类来处理多于两个值的目标。

SVC 是一种特别有用的二分类和多分类技术。它能够处理特征与目标之间的简单和复杂关系。对于几乎所有的监督学习问题，SVMs 至少都应该被考虑。然而，它在处理非常大的数据集时效率并不高。

在下一章中，我们将探讨另一种流行且灵活的分类算法，朴素贝叶斯。

第十四章：第十四章：朴素贝叶斯分类

本章中，我们将探讨在哪些情况下朴素贝叶斯可能比我们迄今为止考察的某些分类器更有效率。朴素贝叶斯是一个非常直观且易于实现的分类器。假设我们的特征是独立的，我们甚至可能比逻辑回归得到更好的性能，尤其是如果我们不使用后者进行正则化的话。

本章中，我们将讨论朴素贝叶斯的基本假设以及算法如何用于解决我们已探索的一些建模挑战，以及一些新的挑战，如文本分类。我们将考虑何时朴素贝叶斯是一个好的选择，何时不是。我们还将检查朴素贝叶斯模型的解释。

本章我们将涵盖以下主题：

关键概念
朴素贝叶斯分类模型
朴素贝叶斯用于文本分类

技术要求

本章中，我们将主要使用 pandas、NumPy 和 scikit-learn 库。唯一的例外是 imbalanced-learn 库，可以使用 pip install imbalanced-learn 安装。本章中的所有代码都使用 scikit-learn 版本 0.24.2 和 1.0.2 进行了测试。

关键概念

朴素贝叶斯分类器使用贝叶斯定理来预测类别成员资格。贝叶斯定理描述了事件发生的概率与给定新、相关数据的事件发生概率之间的关系。给定新数据的事件的概率称为后验概率。在新的数据之前发生事件的概率适当地称为先验概率。

贝叶斯定理给出了以下方程：

后验概率（给定新数据的事件的概率）等于数据给定事件的概率，乘以事件的先验概率，除以新数据的概率。

稍微不那么口语化地，这通常如下所示：

在这里，A 是一个事件，例如类别成员资格，而 B 是新信息。当应用于分类时，我们得到以下方程：

在这里，是给定实例的特征后实例属于该类别的概率，而是给定类别成员资格的特征概率。P(y) 是类别成员资格的概率，而是特征值的概率。因此，后验概率，，等于给定类别成员资格的特征值概率，乘以类别成员资格的概率，除以特征值的概率。

这里的假设是特征之间相互独立。这就是给这个方法带来朴素这个形容词的原因。然而，作为一个实际问题，特征独立性并不是使用朴素贝叶斯获得良好结果所必需的。

简单贝叶斯可以处理数值或分类特征。当我们主要拥有数值特征时，我们通常使用高斯贝叶斯。正如其名所示，高斯贝叶斯假设特征值的条件概率遵循正态分布。然后可以使用每个类中特征的标准差和均值相对简单地计算出。

当我们的特征是离散的或计数时，我们可以使用多项式贝叶斯。更普遍地说，当特征值的条件概率遵循多项式分布时，它效果很好。多项式贝叶斯的一个常见应用是与使用词袋方法的文本分类。在词袋中，特征是文档中每个词的计数。我们可以应用贝叶斯定理来估计类成员的概率：

在这里，是给定一个词频向量W时属于k类的概率。我们将在本章的最后部分充分利用这一点。

简单贝叶斯适用于相当广泛的文本分类任务。它在情感分析、垃圾邮件检测和新闻故事分类等方面都有应用，仅举几个例子。

简单贝叶斯是一种既适用于训练又适用于预测的高效算法，并且通常表现良好。它相当可扩展，可以很好地处理大量实例和特征。它也非常容易解释。当模型复杂性对于良好的预测不是必需的时候，算法表现最佳。即使简单贝叶斯不太可能是产生最少偏差的方法，它也经常用于诊断目的，或者检查不同算法的结果。

我们在上一章中使用过的 NBA 数据可能是用简单贝叶斯建模的好候选。我们将在下一节中探讨这一点。

简单贝叶斯分类模型

简单贝叶斯的一个吸引力在于，即使数据量很大，你也能快速获得不错的结果。在系统资源上，拟合和预测都相对容易。另一个优点是，可以捕捉相对复杂的关系，而无需转换特征空间或进行大量的超参数调整。我们可以用我们在上一章中使用的 NBA 数据来证明这一点。

在本节中，我们将使用关于国家篮球协会（NBA）比赛的数据。该数据集包含了从 2017/2018 赛季到 2020/2021 赛季每场 NBA 比赛的统计数据。这包括主队；主队是否获胜；客队；客队和主队的投篮命中率；两队的失误、篮板和助攻；以及其他一些指标。

注意

NBA 比赛数据可以通过以下链接由公众下载：www.kaggle.com/datasets/wyattowalsh/basketball。这个数据集包含从 1946/1947 赛季开始的比赛数据。它使用nba_api从nba.com获取统计数据。该 API 可在github.com/swar/nba_api找到。

让我们使用朴素贝叶斯构建一个分类模型：

我们将加载我们在过去几章中使用过的相同库：

import pandas as pd
import numpy as np
from sklearn.model_selection import train_test_split
from sklearn.preprocessing import OneHotEncoder, StandardScaler
from sklearn.impute import SimpleImputer
from sklearn.pipeline import make_pipeline
from sklearn.compose import ColumnTransformer
from sklearn.feature_selection import RFE
from sklearn.naive_bayes import GaussianNB
from sklearn.linear_model import LogisticRegression
from sklearn.model_selection import cross_validate, \
  RandomizedSearchCV, RepeatedStratifiedKFold
import sklearn.metrics as skmet
import os
import sys
sys.path.append(os.getcwd() + "/helperfunctions")
from preprocfunc import OutlierTrans

接下来，我们将加载 NBA 比赛数据。在这里我们需要进行一些数据清理。一些观测值在主队是否获胜（WL_HOME）方面有缺失值。我们将删除这些值，因为那将是我们的目标。我们还将WL_HOME转换为整数。请注意，没有太多的类别不平衡，我们不需要采取激进的措施来处理它：

nbagames = pd.read_csv("data/nbagames2017plus.csv", parse_dates=['GAME_DATE'])
nbagames = nbagames.loc[nbagames.WL_HOME.isin(['W','L'])]
nbagames.shape
(4568, 149)
nbagames['WL_HOME'] = \
  np.where(nbagames.WL_HOME=='L',0,1).astype('int')
nbagames.WL_HOME.value_counts(dropna=False)
1    2586
0    1982
Name: WL_HOME, dtype: int64

现在，让我们创建训练和测试数据框，按数值和分类特征组织它们。我们还应该生成一些描述性统计。由于我们在上一章中已经做了，所以这里不再重复；然而，回顾那些数字可能有助于为建模阶段做好准备：

num_cols = ['FG_PCT_HOME','FTA_HOME','FG3_PCT_HOME',
  'FTM_HOME','FT_PCT_HOME','OREB_HOME','DREB_HOME',
  'REB_HOME','AST_HOME','STL_HOME','BLK_HOME',
  'TOV_HOME', 'FG_PCT_AWAY','FTA_AWAY','FG3_PCT_AWAY',
  'FT_PCT_AWAY','OREB_AWAY','DREB_AWAY','REB_AWAY',
  'AST_AWAY','STL_AWAY','BLK_AWAY','TOV_AWAY']
cat_cols = ['TEAM_ABBREVIATION_HOME','SEASON']
X_train, X_test, y_train, y_test =  \
  train_test_split(nbagames[num_cols + cat_cols],\
  nbagames[['WL_HOME']], test_size=0.2,random_state=0)

现在，我们需要设置列转换。我们将处理数值特征的异常值，将这些值和任何缺失值分配给中位数。然后，我们将使用标准缩放器。我们将为分类特征设置独热编码：

ohe = OneHotEncoder(drop='first', sparse=False)
cattrans = make_pipeline(ohe)
standtrans = make_pipeline(OutlierTrans(2),
  SimpleImputer(strategy="median"), StandardScaler())
coltrans = ColumnTransformer(
  transformers=[
    ("cat", cattrans, cat_cols),
    ("stand", standtrans, num_cols)
  ]
)

现在，我们已经准备好运行一个朴素贝叶斯分类器。我们将在列转换和一些递归特征消除之后，将高斯朴素贝叶斯实例添加到一个管道中：
```
nb = GaussianNB()
rfe = RFE(estimator=LogisticRegression(),
  n_features_to_select=15)
pipe1 = make_pipeline(coltrans, rfe, nb)
```

让我们用 K 折交叉验证来评估这个模型。我们得到了不错的分数，虽然不如上一章中支持向量分类的分数好：

kf = RepeatedStratifiedKFold(n_splits=7,n_repeats=10,\
  random_state=0)
scores = cross_validate(pipe1, X_train, \
  y_train.values.ravel(), \
  scoring=['accuracy','precision','recall','f1'], \
  cv=kf, n_jobs=-1)
print("accuracy: %.2f, precision: %.2f, 
  sensitivity: %.2f, f1: %.2f"  %
  (np.mean(scores['test_accuracy']),\
   np.mean(scores['test_precision']),\
   np.mean(scores['test_recall']),\
   np.mean(scores['test_f1'])))
accuracy: 0.81, precision: 0.84, sensitivity: 0.83, f1: 0.83

对于高斯朴素贝叶斯，我们只有一个超参数需要担心调整。我们可以通过var_smoothing超参数确定使用多少平滑度。我们可以进行随机网格搜索以找出最佳值。

var_smoothing超参数决定了添加到方差中的量，这将导致模型对接近平均值实例的依赖性降低：

nb_params = {
    'gaussiannb__var_smoothing': np.logspace(0,-9, num=100)
}
rs = RandomizedSearchCV(pipe1, nb_params, cv=kf, \
  scoring='accuracy')
rs.fit(X_train, y_train.values.ravel())

我们得到了更好的准确性：

rs.best_params_
{'gaussiannb__var_smoothing': 0.657933224657568}
rs.best_score_
0.8608648056923919

我们还应该查看不同迭代的结果。正如我们所见，较大的平滑值表现更好：

results = \
  pd.DataFrame(rs.cv_results_['mean_test_score'], \
    columns=['meanscore']).\
  join(pd.DataFrame(rs.cv_results_['params'])).\
  sort_values(['meanscore'], ascending=False)
results
         meanscore     gaussiannb__var_smoothing
2        0.86086        0.65793
1        0.85118        0.03511
9        0.81341        0.00152
5        0.81212        0.00043
7        0.81180        0.00019
8        0.81169        0.00002
3        0.81152        0.00000
6        0.81152        0.00000
0        0.81149        0.00000
4        0.81149        0.00000

我们还可以查看每次迭代的平均拟合和评分时间：

print("fit time: %.3f, score time: %.3f"  %
  (np.mean(rs.cv_results_['mean_fit_time']),\
  np.mean(rs.cv_results_['mean_score_time'])))
fit time: 0.660, score time: 0.029

让我们看看最佳模型的预测结果。除了提高准确性外，敏感性也有所提高，从0.83提升到0.92：

pred = rs.predict(X_test)
print("accuracy: %.2f, sensitivity: %.2f, \
  specificity: %.2f, precision: %.2f"  %
  (skmet.accuracy_score(y_test.values.ravel(), pred),
  skmet.recall_score(y_test.values.ravel(), pred),
  skmet.recall_score(y_test.values.ravel(), pred, \
     pos_label=0),
  skmet.precision_score(y_test.values.ravel(), pred)))
accuracy: 0.86, sensitivity: 0.92, specificity: 0.79, precision: 0.83

同时查看一个混淆矩阵来更好地了解模型的表现也是一个好主意：

cm = skmet.confusion_matrix(y_test, pred)
cmplot = skmet.ConfusionMatrixDisplay(
  confusion_matrix=cm, display_labels=['Loss', 'Won'])
cmplot.plot()
cmplot.ax_.set(title='Home Team Win Confusion Matrix', 
  xlabel='Predicted Value', ylabel='Actual Value')

这产生了以下图表：

图 14.1 – 基于高斯朴素贝叶斯模型的主队获胜混淆矩阵

虽然这还不错，但仍然不如我们之前章节中的支持向量模型好。特别是，我们希望在预测损失方面做得更好。这也在我们之前步骤中看到的相对较低的0.79特异性得分中得到了反映。记住，特异性是我们正确预测实际负值的负值比率。

另一方面，拟合和评分运行得相当快。我们也不需要做太多的超参数调整。朴素贝叶斯在建模二元或多类目标时通常是一个好的起点。

朴素贝叶斯已经成为文本分类中更受欢迎的选项。我们将在下一节中使用它。

文本分类的朴素贝叶斯

也许令人惊讶的是，一个基于计算条件概率的算法对文本分类是有用的。但这与一个关键简化假设相当直接。让我们假设我们的文档可以通过文档中每个单词的计数很好地表示，不考虑单词顺序或语法。这被称为词袋。词袋与分类目标之间的关系——比如说，垃圾邮件/非垃圾邮件或正面/负面——可以用多项式朴素贝叶斯成功建模。

在本节中，我们将使用短信数据。我们将使用的数据集包含垃圾邮件和非垃圾邮件的标签。

注意

该短信数据集可以通过www.kaggle.com/datasets/team-ai/spam-text-message-classification供公众下载。它包含两列：短信文本和垃圾邮件或非垃圾邮件（ham）标签。

让我们用朴素贝叶斯进行一些文本分类：

我们将需要几个我们在这本书中尚未使用的模块。我们将导入MultinomialNB，这是我们构建多项式朴素贝叶斯模型所需的。我们还需要CountVectorizer来创建词袋。我们将导入SMOTE模块来处理类别不平衡。请注意，我们将使用一个imbalanced-learn管道而不是一个scikit-learn管道。这是因为我们将在我们的管道中使用SMOTE：
```
import pandas as pd
import numpy as np
from sklearn.model_selection import train_test_split
from imblearn.pipeline import make_pipeline
from imblearn.over_sampling import SMOTE
from sklearn.naive_bayes import MultinomialNB
from sklearn.feature_extraction.text import CountVectorizer
import sklearn.metrics as skmet
```
注意

在本节中，我们使用SMOTE进行过采样；也就是说，我们将复制代表性不足的类别的实例。当我们担心我们的模型在捕捉一个类别的变化方面做得不好，因为我们相对于一个或多个其他类别的实例太少时，过采样可以是一个好的选择。过采样会复制该类别的实例。

接下来，我们将加载短信数据集。我们将把我们的目标转换为整数变量，并确认它按预期工作。注意显著的类别不平衡。让我们查看前几行，以更好地了解数据：

spamtext = pd.read_csv("data/spamtext.csv")
spamtext['spam'] = np.where(spamtext.category=='spam',1,0)
spamtext.groupby(['spam','category']).size()
spam  category
0     ham         4824
1     spam         747
dtype: int64
spamtext.head()
 category  message                                spam
0  ham     Go until jurong point, crazy..         0
1  ham     Ok lar... Joking wif u oni...          0
2  spam    Free entry in 2 a wkly comp to win...  1
3  ham     U dun say so early hor... U c already..0
4  ham     Nah I don't think he goes to usf, ..   0

现在，我们创建训练和测试数据框。我们将使用stratify参数来确保训练和测试数据中目标值的分布相等。

我们还将实例化一个CountVectorizer对象来创建我们后面的词袋。我们指出我们想要忽略一些单词，因为它们不提供有用的信息。我们本可以创建一个停用词列表，但在这里，我们将利用 scikit-learn 的英文停用词列表：

X_train, X_test, y_train, y_test =  \
  train_test_split(spamtext[['message']],\
  spamtext[['spam']], test_size=0.2,\
  stratify=spamtext[['spam']], random_state=0)
countvectorizer = CountVectorizer(analyzer='word', \
  stop_words='english')

让我们看看向量器是如何与我们的数据中的几个观察结果一起工作的。为了便于查看，我们只会从包含少于 50 个字符的消息中提取信息。

使用向量器，我们为每个观察结果中使用的所有非停用词获取计数。例如，like在第一条消息中使用了一次，而在第二条消息中一次也没有使用。这给like在转换数据中的第一个观察结果赋予了一个值为1，而在第二个观察结果中赋予了一个值为0。

我们不会在我们的模型中使用这一步骤中的任何内容。我们这样做只是为了说明目的：

smallsample = \
  X_train.loc[X_train.message.str.len()<50].\
    sample(2, random_state=35)
smallsample
                                        message
2079                I can take you at like noon
5393  I dont know exactly could you ask chechi.
ourvec = \
    pd.DataFrame(countvectorizer.\
    fit_transform(smallsample.values.ravel()).\
    toarray(),\
    columns=countvectorizer.get_feature_names())
ourvec
    ask   chechi  dont   exactly  know  like  noon
0    0    0       0      0        0     1     1
1    1    1       1      1        1     0     0

现在，让我们实例化一个MultinomialNB对象并将其添加到管道中。我们将使用SMOTE进行过采样以处理类别不平衡：

nb = MultinomialNB()
smote = SMOTE(random_state=0)
pipe1 = make_pipeline(countvectorizer, smote, nb)
pipe1.fit(X_train.values.ravel(), 
  y_train.values.ravel())

现在，让我们看看一些预测结果。我们得到了令人印象深刻的0.97准确率和同样好的特异性。这种出色的特异性表明我们没有许多误报。相对较低的反应性表明我们没有捕捉到一些正例（垃圾邮件），尽管我们仍然做得相当不错：

pred = pipe1.predict(X_test.values.ravel())
print("accuracy: %.2f, sensitivity: %.2f, specificity: %.2f, precision: %.2f"  %
  (skmet.accuracy_score(y_test.values.ravel(), pred),
  skmet.recall_score(y_test.values.ravel(), pred),
  skmet.recall_score(y_test.values.ravel(), pred, pos_label=0),
  skmet.precision_score(y_test.values.ravel(), pred)))
accuracy: 0.97, sensitivity: 0.87, specificity: 0.98, precision: 0.87

使用混淆矩阵可视化我们模型的表现是有帮助的：

cm = skmet.confusion_matrix(y_test, pred)
cmplot = skmet.ConfusionMatrixDisplay(
  confusion_matrix=cm, \
  display_labels=['Not Spam', 'Spam'])
cmplot.plot()
cmplot.ax_.set(
  title='Spam Prediction Confusion Matrix', 
  xlabel='Predicted Value', ylabel='Actual Value')

这会产生以下图表：

图 14.2 – 使用多项式朴素贝叶斯进行垃圾邮件预测

朴素贝叶斯在构建文本分类模型时可以产生很好的结果。指标通常非常好且非常高效。这是一个非常直接的二元分类问题。然而，朴素贝叶斯也可以在多类文本分类问题中有效。该算法可以以与我们在这里使用多类目标相同的方式应用。

摘要

朴素贝叶斯是一个很好的算法，可以添加到我们的常规工具包中，用于解决分类问题。它并不总是产生最少偏差的预测方法。然而，另一方面也是真的。过拟合的风险较小，尤其是在处理连续特征时。它也非常高效，能够很好地扩展到大量观察和大量特征空间。

本书接下来的两章将探讨无监督学习算法——那些我们没有预测目标的情况。在下一章中，我们将研究主成分分析，然后在下一章中研究 K-means 聚类。