R-机器学习第二版-二-

R 机器学习第二版（二）

原文：annas-archive.org/md5/e79782597ad902714a518953c5d1e70b

译者：飞龙

协议：CC BY-NC-SA 4.0

第五章. 分治法 – 使用决策树和规则进行分类

在决定是否接受多个薪资和福利不同的工作邀请时，许多人通常通过列出利弊清单，并根据简单的规则排除选项。例如，“如果我需要通勤超过一个小时，我会感到不高兴。”或者，“如果我的收入低于 50k 美元，我将无法养活我的家人。”通过这种方式，预测未来幸福的复杂且艰难的决策可以简化为一系列简单的决策。

本章讲解了决策树和规则学习器——这两种机器学习方法也从一系列简单的选择中做出复杂决策。这些方法然后以逻辑结构的形式呈现其知识，且无需统计学知识也能理解。这一特点使得这些模型在商业战略和流程改进方面特别有用。

在本章结束时，你将学到：

• 树和规则如何“贪婪地”将数据划分为有趣的片段

• 最常见的决策树和分类规则学习器，包括 C5.0、1R 和 RIPPER 算法

• 如何使用这些算法执行现实世界中的分类任务，例如识别风险较高的银行贷款和有毒蘑菇

我们将从研究决策树开始，随后探讨分类规则。然后，我们将通过预览后面的章节来总结所学内容，这些章节讨论了将树和规则作为基础的更先进的机器学习技术。

理解决策树

决策树学习器是强大的分类器，它们利用树结构来建模特征与潜在结果之间的关系。正如下图所示，这种结构得名于其形状类似于字面意义上的树木，它从粗大的树干开始，向上延伸时分成越来越细的枝条。决策树分类器以相同的方式使用分支决策的结构，将示例引导到最终的预测类别值。

为了更好地理解其实际应用，假设我们考虑以下这棵树，它预测一个工作邀请是否应该被接受。一个待考虑的工作邀请从根节点开始，然后通过决策节点传递，根据工作的属性做出选择。这些选择将数据分割到分支上，指示决策的潜在结果，这里表现为“是”或“否”的结果，尽管在某些情况下可能有多个可能性。如果可以做出最终决定，树将通过叶节点（也叫终端节点）终止，叶节点表示一系列决策的结果应该采取的行动。在预测模型中，叶节点提供给定树中一系列事件后的预期结果。

决策树算法的一个巨大优势是，其流程图式的树结构并不一定仅供学习者内部使用。在模型创建之后，许多决策树算法会以人类可读的格式输出生成的结构。这为我们提供了对模型如何以及为何在某项任务中工作或不工作的深刻理解。这也使得决策树特别适用于那些需要因法律原因或为了向他人分享结果以指导未来业务实践而使分类机制透明的应用场景。考虑到这一点，一些潜在的应用包括：

• 信贷评分模型，其中导致申请人被拒绝的标准需要明确记录，并且必须排除偏见。

• 关于客户行为的营销研究，如满意度或流失率，这些研究将与管理层或广告公司共享。

• 基于实验室测量、症状或疾病进展速度的医疗状况诊断。

尽管前述应用展示了树在决策过程中提供价值，但这并不意味着它们的效用就此结束。事实上，决策树可能是最广泛使用的机器学习技术之一，可以应用于几乎任何类型的数据——通常可以获得出色的开箱即用的应用效果。

尽管决策树有广泛的应用，但值得注意的是，在某些情况下，树可能不是理想选择。一个这样的例子是当数据包含大量具有多个级别的名义特征，或者数据包含大量数值特征时。这些情况可能导致决策数量非常庞大，并且树结构过于复杂。它们还可能导致决策树出现过拟合数据的倾向，尽管正如我们很快会看到的，甚至这种弱点也可以通过调整一些简单的参数来克服。

分治法

决策树是使用一种名为递归划分的启发式方法构建的。这种方法通常也被称为分治法，因为它将数据划分成子集，然后将这些子集反复划分成更小的子集，依此类推，直到算法判断子集内的数据足够同质化，或满足其他停止准则为止。

要了解如何将数据集拆分以创建决策树，可以想象一个裸根节点，它将成长为一棵成熟的树。最初，根节点代表整个数据集，因为此时尚未进行拆分。接下来，决策树算法必须选择一个特征来进行拆分；理想情况下，它选择的是最能预测目标类别的特征。然后，根据该特征的不同值，示例数据被分成多个组，树的第一个分支就此形成。

沿着每个分支往下，算法继续分治数据，每次选择最合适的特征来创建另一个决策节点，直到达到停止标准。分治法可能会在某个节点停止，情形如下：

• 节点处的所有（或几乎所有）示例都属于同一类

• 没有剩余的特征可以用来区分各个示例

• 该树已经增长到预设的大小限制

为了说明树形结构的构建过程，让我们考虑一个简单的例子。假设你在一家好莱坞电影公司工作，负责决定公司是否应该继续制作那些有潜力的新人作家所提交的电影剧本。度假回来后，你的办公桌上堆满了提案。由于没有时间逐一阅读每份提案，你决定开发一个决策树算法，用来预测一部潜在电影是否会落入以下三类之一：关键成功、主流热片或票房失败。

为了构建决策树，你查阅了公司档案，分析了影响公司最近发布的 30 部电影成功与失败的因素。你很快注意到，电影的预估拍摄预算、一线明星的出演数量和电影的成功程度之间存在某种关系。对这一发现感到兴奋，你制作了一个散点图来展示这种模式：

使用分治策略，我们可以从这些数据中构建一个简单的决策树。首先，为了创建树的根节点，我们根据明星数量这一特征进行划分，将电影分为有无大量一线明星的两组：

接下来，在具有较多明星的电影组中，我们可以再次进行划分，区分有无高预算的电影：

到这个阶段，我们已经将数据划分为三组。图表左上角的组完全由获得好评的电影组成。这个组的特点是有大量的明星出演，且预算相对较低。在右上角，绝大多数电影都是票房热片，具有高预算和大量的明星阵容。最后一个组虽然没有太多明星，但预算从小到大不等，包含了票房失败的电影。

如果我们愿意，我们可以继续通过基于越来越具体的预算和名人数量范围来划分数据，直到每个当前错误分类的值都位于自己的小分区中，并且被正确分类。然而，不建议以这种方式过度拟合决策树。虽然没有什么可以阻止我们无限制地划分数据，但过于具体的决策并不总是能够更广泛地泛化。我们将通过在此停止算法来避免过拟合的问题，因为每个组中超过 80%的示例都来自同一类。这构成了我们停止标准的基础。

提示

你可能已经注意到，斜线可能会更加干净地划分数据。这是决策树知识表示的一个局限性，它使用的是轴对齐划分。每次划分只考虑一个特征，这使得决策树无法形成更复杂的决策边界。例如，可以通过一个决策来创建一条斜线，询问：“名人数量是否大于预估预算？”如果是，那么“它将是一个关键性的成功”。

我们用于预测电影未来成功的模型可以用一个简单的树表示，如下图所示。为了评估剧本，按照每个决策的分支，直到预测出剧本的成功或失败。很快，你将能够从一堆积压的剧本中识别出最有前途的选项，然后回到更重要的工作，如写奥斯卡颁奖典礼的获奖感言。

由于现实世界的数据包含超过两个特征，决策树很快就会变得比这复杂得多，包含更多的节点、分支和叶子。在接下来的部分中，你将学习一种流行的算法，它可以自动构建决策树模型。

C5.0 决策树算法

决策树有很多实现方式，但最著名的一种实现是C5.0 算法。该算法由计算机科学家 J. Ross Quinlan 开发，是他之前算法C4.5的改进版，而C4.5本身则是他迭代二分法 3（ID3）算法的改进。尽管 Quinlan 将 C5.0 推向商业客户（详情见www.rulequest.com/），但该算法的单线程版本的源代码已经公开，因此被像 R 这样的程序所采用。

注意

更令人困惑的是，一个流行的基于 Java 的开源替代方案J48，它是 C4.5 的替代品，已包含在 R 的RWeka包中。由于 C5.0、C4.5 和 J48 之间的差异很小，本章中的原理将适用于这三种方法，且这些算法应被视为同义。

C5.0 算法已成为生成决策树的行业标准，因为它对大多数类型的问题在开箱即用时表现出色。与其他先进的机器学习模型（如第七章中描述的黑箱方法 – 神经网络和支持向量机）相比，C5.0 生成的决策树通常表现几乎相同，但更易于理解和部署。此外，如下表所示，该算法的弱点相对较小，并且大多可以避免：

优势	弱点

|

• 一种多用途的分类器，能够在大多数问题上表现良好

• 高度自动化的学习过程，能够处理数值型或名义型特征，以及缺失数据

• 排除不重要的特征

• 可用于小型和大型数据集

• 生成的模型可以在没有数学背景的情况下进行解释（对于相对较小的树）

• 比其他复杂模型更高效

|

• 决策树模型通常偏向于对具有大量层级的特征进行划分

• 容易出现过拟合或欠拟合问题

• 由于依赖于轴对齐的划分，可能在建模某些关系时遇到困难

• 训练数据的微小变化可能会导致决策逻辑的巨大变化

• 大型树可能难以解释，且它们做出的决策可能显得不合直觉

|

为了简化问题，我们之前的决策树示例忽略了机器如何运用分治策略的数学原理。让我们更详细地探讨这个问题，研究这种启发式方法在实践中的工作原理。

选择最佳的划分

决策树面临的第一个挑战是识别应该在哪个特征上进行划分。在之前的示例中，我们寻找了一种划分数据的方法，使得划分后的数据主要包含单一类别的示例。一个示例子集仅包含单一类别的程度被称为纯度，任何仅由单一类别组成的子集都被称为纯。

有多种纯度度量方法可以用来识别最佳的决策树划分候选。C5.0 使用熵，这是一个借自信息论的概念，用于量化一个类值集合中的随机性或无序性。熵高的集合非常多样化，几乎不能提供关于其他可能属于该集合的项目的信息，因为没有明显的共同点。决策树希望找到减少熵的划分，从而最终增加组内的同质性。

通常情况下，熵以比特为单位进行度量。如果只有两个可能的类别，熵值的范围为 0 到 1。对于n个类别，熵的范围是从 0 到log2。在每种情况下，最小值表示样本完全同质，而最大值表示数据尽可能多样化，且没有任何一个群体占据主导地位。

在数学符号中，熵被定义如下：

在这个公式中，对于给定的数据片段*(S)，术语c表示类别的数量，而p[i]表示属于类别级别i*的值的比例。例如，假设我们有一个数据分区，其中包含两个类别：红色（60%）和白色（40%）。我们可以按如下方式计算熵：

> -0.60 * log2(0.60) - 0.40 * log2(0.40)
[1] 0.9709506

我们可以检查所有可能的二类排列的熵。如果我们知道一个类别中示例的比例为x，那么另一个类别的比例就是*(1 – x)。通过使用curve()函数，我们可以绘制出所有可能的x*值下的熵：

> curve(-x * log2(x) - (1 - x) * log2(1 - x),
 col = "red", xlab = "x", ylab = "Entropy", lwd = 4)

这会产生如下图所示的结果：

正如在x = 0.50时熵的峰值所示，50-50 的划分会导致最大熵。当一个类别逐渐主导另一个类别时，熵会减少到零。

为了使用熵来确定最佳划分特征，算法会计算在每个可能的特征上进行划分后同质性变化的量，这个度量被称为信息增益。特征F的信息增益是通过计算划分前的片段熵*(S[1])与划分后分区的熵(S[2])*之间的差值来得到的：

一个复杂的地方是，划分后数据被分成了多个分区。因此，计算*Entropy(S[2])*的函数需要考虑所有分区的总熵。它通过根据每个分区中记录所占比例来加权每个分区的熵。可以用以下公式来表示：

简单来说，划分后的总熵是每个n个分区的熵之和，加权每个分区的示例比例（w[i]）。

信息增益越高，特征在此特征上进行划分后，生成同质群体的效果越好。如果信息增益为零，则说明在该特征上进行划分不会减少熵。另一方面，最大信息增益等于划分前的熵。这意味着划分后的熵为零，表示该划分结果形成了完全同质的群体。

之前的公式假设了名义特征，但决策树也使用信息增益对数值特征进行分裂。为此，一个常见的做法是测试不同的分裂方法，将值划分为大于或小于某个数值阈值的组。这将数值特征转换为一个二级类别特征，从而可以像往常一样计算信息增益。选择具有最大信息增益的数值切分点进行分裂。

注意

尽管 C5.0 使用了信息增益，但信息增益并不是构建决策树时可以使用的唯一分裂准则。其他常用的准则包括基尼指数、卡方统计量和增益比。有关这些（以及更多）准则的回顾，请参考 Mingers J. 决策树归纳的选择度量的实证比较。机器学习。1989; 3:319-342。

剪枝决策树

决策树可以无限地生长，选择分裂特征并将数据分割成越来越小的部分，直到每个示例都被完美分类，或者算法无法再找到特征进行分裂。然而，如果树过度生长，许多决策将变得过于具体，模型将过拟合训练数据。剪枝决策树的过程涉及缩小树的大小，以便它能够更好地泛化到未见过的数据。

解决此问题的一个方法是，当树达到一定的决策数量或决策节点仅包含少量示例时，停止树的生长。这被称为早期停止或预剪枝决策树。由于树避免了不必要的工作，这是一个有吸引力的策略。然而，这种方法的一个缺点是，无法知道树是否会错过那些微妙但重要的模式，这些模式如果树生长到更大规模时可能会学习到。

另一种方法，称为后剪枝，包括先生长一棵故意过大的树，并通过剪枝叶子节点将树的大小减少到一个更合适的水平。这通常比预剪枝更有效，因为在没有先生长树的情况下很难确定决策树的最优深度。稍后对树进行剪枝可以确保算法发现了所有重要的数据结构。

注意

剪枝操作的实现细节非常技术性，超出了本书的范围。如需了解一些可用方法的比较，请参阅 Esposito F, Malerba D, Semeraro G. 决策树剪枝方法的比较分析。IEEE 模式分析与机器智能学报。1997;19: 476-491。

C5.0 算法的一个优点是它在修剪过程中有明确的方向——它会使用相当合理的默认设置自动做出许多决策。其整体策略是后期修剪树形结构。它首先生成一个过拟合训练数据的大树，然后删除那些对分类错误影响较小的节点和分支。在某些情况下，整个分支会被移动到树的更高位置，或被更简单的决策所替代。这些移植分支的过程分别被称为子树提升和子树替换。

平衡决策树的过拟合与欠拟合是一项艺术，但如果模型准确性至关重要，花时间调整不同的修剪选项，看看是否能提高测试数据的表现，是值得投入的。正如你将很快看到的，C5.0 算法的一个优点是它非常容易调整训练选项。

示例 – 使用 C5.0 决策树识别高风险银行贷款

2007-2008 年的全球金融危机突显了银行业务中透明度和严格性的重要性。由于信贷供应受到限制，银行收紧了贷款系统，并转向机器学习，以更准确地识别高风险贷款。

由于决策树具有高准确性和用通俗语言制定统计模型的能力，因此在银行业得到了广泛应用。由于许多国家的政府组织严格监控贷款实践，银行高层必须能够解释为什么一个申请人被拒绝贷款，而其他申请人却被批准。这些信息对于希望了解为什么自己的信用评级不合格的客户也非常有用。

自动化信用评分模型可能被用来在电话和网络上即时批准信用申请。在本节中，我们将使用 C5.0 决策树开发一个简单的信用批准模型。我们还将看到如何调整模型结果，以最小化可能导致机构经济损失的错误。

步骤 1 – 收集数据

我们的信用模型背后的理念是识别那些能够预测较高违约风险的因素。因此，我们需要获取大量过去银行贷款的数据，了解这些贷款是否发生了违约，以及有关申请人的信息。

具有这些特征的数据可以在由汉斯·霍夫曼（Hans Hofmann）捐赠给 UCI 机器学习数据仓库的一个数据集中找到 (archive.ics.uci.edu/ml)。该数据集包含来自德国一家信用机构的贷款信息。

提示

本章中展示的数据集与原始数据集略有修改，目的是消除一些预处理步骤。为了跟随示例操作，请从 Packt Publishing 的网站下载 credit.csv 文件并将其保存到你的 R 工作目录中。

信用数据集包括 1,000 个贷款实例，另外还有一组数值型和名义型特征，表示贷款和贷款申请者的特点。一个类别变量表示贷款是否违约。让我们看看是否能发现一些预测这一结果的模式。

步骤 2 – 探索和准备数据

正如我们之前所做的那样，我们将使用read.csv()函数导入数据。我们将忽略stringsAsFactors选项，因此使用默认值TRUE，因为数据中的大多数特征都是名义型的：

> credit <- read.csv("credit.csv")

str()函数的前几行输出如下：

> str(credit)
'data.frame':1000 obs. of  17 variables:
 $ checking_balance : Factor w/ 4 levels "< 0 DM","> 200 DM",..
 $ months_loan_duration: int  6 48 12 ...
 $ credit_history      : Factor w/ 5 levels "critical","good",..
 $ purpose             : Factor w/ 6 levels "business","car",..
 $ amount              : int  1169 5951 2096 ...

我们看到期望的 1,000 条观察数据和 17 个特征，这些特征是因子和整数数据类型的组合。

让我们看一下table()函数输出的几个可能预测违约的贷款特征。申请者的支票和储蓄账户余额被记录为分类变量：

> table(credit$checking_balance)
 < 0 DM   > 200 DM 1 - 200 DM    unknown 
 274         63        269        394
> table(credit$savings_balance)
 < 100 DM > 1000 DM  100 - 500 DM 500 - 1000 DM   unknown 
 603        48           103            63       183

支票和储蓄账户余额可能是预测贷款违约状态的重要指标。请注意，由于贷款数据来自德国，因此货币记录为德国马克（DM）。

贷款的某些特征是数值型的，例如贷款的期限和请求的信用金额：

> summary(credit$months_loan_duration)
 Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
 4.0    12.0    18.0    20.9    24.0    72.0 
> summary(credit$amount)
 Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
 250    1366    2320    3271    3972   18420

贷款金额从 250 DM 到 18,420 DM 不等，期限从 4 个月到 72 个月不等，贷款的中位数期限为 18 个月，金额为 2,320 DM。

default向量表示贷款申请者是否未能按约定的付款条款履约并进入违约。该数据集中有 30%的贷款进入了违约状态：

> table(credit$default)
 no yes 
700 300

高违约率对银行来说是不利的，因为这意味着银行不太可能完全收回其投资。如果我们成功了，我们的模型将能够识别出高违约风险的申请者，从而允许银行拒绝这些申请的信用请求。

数据准备 – 创建随机的训练和测试数据集

正如我们在前面的章节中所做的那样，我们将把数据分成两部分：一个用于构建决策树的训练数据集和一个用于评估模型在新数据上表现的测试数据集。我们将使用 90%的数据用于训练，10%的数据用于测试，这将为我们提供 100 条记录来模拟新申请者。

与前几章使用的数据是按随机顺序排序的不同，我们简单地将数据集分为两部分，取前 90%的记录用于训练，剩余的 10%用于测试。相比之下，信用数据集没有经过随机排序，因此采用之前的方法是不明智的。假设银行按贷款金额排序数据，最大贷款位于文件末尾。如果我们使用前 90%的数据用于训练，剩余 10%用于测试，我们将仅在小额贷款上训练模型，而在大额贷款上测试模型。显然，这可能会造成问题。

我们将通过使用信用数据的随机样本来解决这个问题。随机样本仅仅是一个随机选择记录子集的过程。在 R 中，sample()函数用于执行随机抽样。然而，在执行之前，一个常见的做法是设置种子值，这样可以确保随机化过程遵循一个可以在以后复制的序列。看起来这似乎违背了生成随机数的目的，但这样做是有原因的。通过set.seed()函数提供种子值可以确保如果将来重复分析，可以得到相同的结果。

提示

你可能会想，怎么一个所谓的随机过程可以设置种子来产生相同的结果呢？这是因为计算机使用一种名为伪随机数生成器的数学函数来创建看似非常随机的随机数序列，但实际上只要知道序列中前一个值，它们是可以预测的。实际上，现代伪随机数序列与真正的随机序列几乎无法区分，但它们的优势在于计算机可以快速、轻松地生成这些序列。

以下命令使用sample()函数从 1 到 1000 的整数序列中随机选择 900 个值。请注意，set.seed()函数使用了一个任意值123。如果省略这个种子值，你的训练和测试集划分将与本章其余部分所示的结果不同：

> set.seed(123)
> train_sample <- sample(1000, 900)

正如预期的那样，生成的train_sample对象是一个包含 900 个随机整数的向量：

> str(train_sample)
 int [1:900] 288 788 409 881 937 46 525 887 548 453 ... 

通过使用这个向量从信用数据中选择行，我们可以将其分为我们所需的 90%的训练数据集和 10%的测试数据集。请记住，在选择测试记录时使用的破折号运算符告诉 R 选择那些不在指定行中的记录；换句话说，测试数据仅包含那些不在训练样本中的行。

> credit_train <- credit[train_sample, ]
> credit_test  <- credit[-train_sample, ]

如果一切顺利，我们应该在每个数据集中都有大约 30%的违约贷款：

> prop.table(table(credit_train$default))
 no       yes 
0.7033333 0.2966667 

> prop.table(table(credit_test$default))
 no  yes 
0.67 0.33

这似乎是一个相当均匀的划分，所以我们现在可以构建我们的决策树。

提示

如果你的结果不完全匹配，请确保在创建train_sample向量之前立即运行了命令set.seed(123)。

第 3 步 – 在数据上训练模型

我们将使用C50包中的 C5.0 算法来训练我们的决策树模型。如果你还没有安装该包，可以通过install.packages("C50")来安装，并使用library(C50)将其加载到 R 会话中。

以下语法框列出了构建决策树时最常用的一些命令。与我们之前使用的机器学习方法相比，C5.0 算法提供了更多定制模型以适应特定学习问题的方式，但也提供了更多选项。一旦加载了C50包，?C5.0Control命令将显示帮助页面，以获得有关如何精细调整算法的更多细节。

对于我们信用审批模型的第一次迭代，我们将使用默认的 C5.0 配置，如下所示的代码。credit_train的第 17 列是default类别变量，因此我们需要将其从训练数据框中排除，但作为分类的目标因子向量提供：

> credit_model <- C5.0(credit_train[-17], credit_train$default)

credit_model对象现在包含一个 C5.0 决策树。我们可以通过输入它的名称来查看树的一些基本数据：

> credit_model

Call:
C5.0.default(x = credit_train[-17], y = credit_train$default)

Classification Tree
Number of samples: 900 
Number of predictors: 16 

Tree size: 57 

Non-standard options: attempt to group attributes

上述文本显示了有关决策树的一些基本事实，包括生成该树的函数调用、特征数量（标记为predictors）和用于生成树的示例（标记为samples）。还列出了树的大小为 57，表示该树有 57 个决策层级——比我们之前考虑的示例树要大得多！

要查看树的决策，我们可以在模型上调用summary()函数：

> summary(credit_model)

这将产生以下输出：

上述输出显示了决策树中的一些初步分支。前三行可以用简单的语言表示为：

1. 如果支票账户余额未知或大于 200 德国马克，则分类为“ unlikely to default”（不太可能违约）。

2. 否则，如果支票账户余额小于零德国马克或在 1 至 200 德国马克之间。

3. 如果信用历史完美或非常好，则分类为“ likely to default”（可能违约）。

括号中的数字表示满足该决策标准的示例数量，以及被该决策错误分类的数量。例如，在第一行，412/50表示在达到该决策的 412 个示例中，50 个被错误分类为不太可能违约。换句话说，尽管模型预测相反，实际上有 50 个申请人违约了。

提示

有时，决策树会产生一些逻辑上没有意义的决策。例如，为什么信用历史非常好的申请人可能会违约，而支票账户余额未知的申请人不太可能违约？像这样的矛盾规则有时会出现。它们可能反映了数据中的真实模式，或者可能是统计异常。无论是哪种情况，调查这些奇怪的决策，看看树的逻辑是否适用于业务使用，都是很重要的。

在决策树之后，summary(credit_model)的输出显示了一个混淆矩阵，这是一个交叉表，表示模型在训练数据中错误分类的记录：

Evaluation on training data (900 cases):

 Decision Tree 
 ---------------- 
 Size      Errors 
 56  133(14.8%)   <<

 (a)   (b)    <-classified as
 ----  ----
 598    35    (a): class no
 98   169    (b): class yes

错误输出指出，模型正确分类了 900 个训练实例中除了 133 个实例之外的所有实例，错误率为 14.8%。总共有 35 个实际的"no"被错误分类为"yes"（假阳性），而 98 个"yes"被错误分类为"no"（假阴性）。

决策树因其容易将模型过度拟合训练数据而闻名。因此，报告的训练数据错误率可能过于乐观，特别重要的是要在测试数据集上评估决策树。

第四步 – 评估模型性能

为了将我们的决策树应用于测试数据集，我们使用predict()函数，如下面的代码行所示：

> credit_pred <- predict(credit_model, credit_test)

这创建了一个预测类别值的向量，我们可以使用gmodels包中的CrossTable()函数将其与实际类别值进行比较。将prop.c和prop.r参数设置为FALSE可以从表格中移除列和行的百分比。剩余的百分比（prop.t）表示单元格中记录占总记录数的比例：

> library(gmodels)
> CrossTable(credit_test$default, credit_pred,
 prop.chisq = FALSE, prop.c = FALSE, prop.r = FALSE,
 dnn = c('actual default', 'predicted default'))

这导致了以下表格：

在 100 个测试贷款申请记录中，我们的模型正确预测了 59 个未违约和 14 个违约，准确率为 73%，错误率为 27%。这比其在训练数据上的表现稍差，但考虑到模型在未见数据上的表现通常较差，这并不意外。还要注意，模型只正确预测了测试数据中 33 个实际贷款违约中的 14 个，正确率为 42%。不幸的是，这种类型的错误可能是一个非常昂贵的错误，因为银行在每次违约时都会损失资金。让我们看看是否可以通过更多的努力来改进结果。

第五步 – 改进模型性能

我们模型的错误率可能太高，无法在实时信用评分应用中部署。事实上，如果模型对每个测试案例都预测为"无违约"，它的正确率将是 67%，这与我们模型的结果差不多，但所需的努力要小得多！从 900 个样本中预测贷款违约似乎是一个具有挑战性的问题。

更糟糕的是，我们的模型在识别那些确实违约的申请人时表现尤其糟糕。幸运的是，有几种简单的方法可以调整 C5.0 算法，这可能有助于提高模型的整体性能，并改善那些更为昂贵的错误类型。

提高决策树的准确率

C5.0 算法通过增加自适应提升（adaptive boosting）改进了 C4.5 算法。这是一个过程，其中构建了多个决策树，并且这些树对每个实例的最佳类别进行投票。

注意

提升的理念主要基于 Rob Schapire 和 Yoav Freund 的研究。欲了解更多信息，请尝试在网上搜索他们的出版物或他们的近期教材《Boosting: Foundations and Algorithms》。MIT 出版社（2012 年）。

由于提升可以更广泛地应用于任何机器学习算法，它将在本书后续章节中详细介绍，第十一章，提高模型性能。现在，我们只需说提升基于这样一个概念：通过将多个表现较弱的学习者结合起来，可以创建一个比任何单独的学习者都强大的团队。每个模型都有独特的优缺点，它们在解决某些问题时可能表现得更好或更差。因此，使用多个具有互补优缺点的学习者的组合，可以显著提高分类器的准确性。

C5.0()函数使得在我们的 C5.0 决策树中添加提升变得非常简单。我们只需要添加一个额外的trials参数，指示要在提升团队中使用的单独决策树的数量。trials参数设置了上限；如果算法识别到额外的试验似乎并没有改善准确性，它将停止添加树。我们将从 10 次试验开始，这是一个事实上的标准，研究表明这样可以将测试数据的错误率降低大约 25%：

> credit_boost10 <- C5.0(credit_train[-17], credit_train$default,
 trials = 10)

在检查结果模型时，我们可以看到一些额外的线条被添加进来，表明了变化：

> credit_boost10
Number of boosting iterations: 10 
Average tree size: 47.5

在这 10 次迭代中，我们的树的大小缩小了。如果您愿意，可以通过在命令提示符下键入summary(credit_boost10)来查看所有 10 棵树。它还列出了模型在训练数据上的表现：

> summary(credit_boost10)

 (a)   (b)    <-classified as
 ----  ----
 629     4    (a): class no
 30   237    (b): class yes

分类器在 900 个训练样本上犯了 34 个错误，错误率为 3.8%。这相比我们在添加提升前注意到的 13.9%的训练误差率有了相当大的改进！然而，是否能在测试数据上看到类似的改进仍然有待观察。让我们来看一下：

> credit_boost_pred10 <- predict(credit_boost10, credit_test)
> CrossTable(credit_test$default, credit_boost_pred10,
 prop.chisq = FALSE, prop.c = FALSE, prop.r = FALSE,
 dnn = c('actual default', 'predicted default'))

结果表格如下：

在这里，我们将总错误率从提升前的 27%降低到了提升后模型的 18%。这看起来并不是一个很大的增益，但实际上它比我们预期的 25%的降低要大。另一方面，模型在预测违约方面仍然表现不佳，只有*20/33 = 61%*的预测是正确的。没有看到更大改进的原因可能是我们相对较小的训练数据集，或者这可能只是一个非常难以解决的问题。

话虽如此，如果提升（boosting）可以如此轻松地添加，为什么不默认将其应用于每个决策树呢？原因有二。首先，如果建立一个决策树需要大量的计算时间，那么构建多个树可能在计算上是不可行的。其次，如果训练数据非常嘈杂，那么提升可能根本不会带来改进。不过，如果需要更高的准确度，尝试一下还是值得的。

使错误的成本比其他错误更高

向可能违约的申请人发放贷款可能是一个昂贵的错误。减少错误负样本数量的一种解决方案可能是拒绝更多边缘申请人，假设银行从高风险贷款中获得的利息远远不能弥补若贷款完全无法偿还时所遭受的巨额损失。

C5.0 算法允许我们为不同类型的错误分配惩罚，以避免决策树犯更多代价更高的错误。惩罚值被指定在成本矩阵中，矩阵定义了每个错误相对于其他预测的成本。

要开始构建成本矩阵，我们需要首先指定维度。由于预测值和实际值都可以取“是”或“否”两种值，我们需要描述一个 2 x 2 矩阵，使用两个向量的列表，每个向量包含两个值。同时，我们还将为矩阵的维度命名，以避免日后混淆：

> matrix_dimensions <- list(c("no", "yes"), c("no", "yes"))
> names(matrix_dimensions) <- c("predicted", "actual")

检查新对象时，表明我们的维度已正确设置：

> matrix_dimensions
$predicted
[1] "no"  "yes"

$actual
[1] "no"  "yes"

接下来，我们需要通过提供四个值来为各种类型的错误分配惩罚，以填充矩阵。由于 R 通过从上到下依次填充列来填充矩阵，我们需要按特定顺序提供这些值：

• 预测为否，实际为否

• 预测为是，实际为否

• 预测为否，实际为是

• 预测为是，实际为是

假设我们认为贷款违约对银行的成本是错失机会的四倍。那么我们的惩罚值可以定义为：

> error_cost <- matrix(c(0, 1, 4, 0), nrow = 2,
 dimnames = matrix_dimensions)

这将创建以下矩阵：

> error_cost
 actual
predicted no yes
 no   0   4
 yes  1   0

根据这个矩阵的定义，当算法正确地分类为“否”或“是”时没有任何成本，但错误负样本的成本为 4，而错误正样本的成本为 1。为了了解这如何影响分类，我们可以将其应用到决策树中，使用C5.0()函数的costs参数。其他步骤与我们之前做的一样：

> credit_cost <- C5.0(credit_train[-17], credit_train$default,
 costs = error_cost)
> credit_cost_pred <- predict(credit_cost, credit_test)
> CrossTable(credit_test$default, credit_cost_pred,
 prop.chisq = FALSE, prop.c = FALSE, prop.r = FALSE,
 dnn = c('actual default', 'predicted default'))

这将生成以下混淆矩阵：

与我们提升后的模型相比，这个版本的错误更多：这里的错误率为 37%，而提升模型的错误率为 18%。然而，错误的类型却大不相同。之前的模型仅有 42%和 61%的违约案例被正确分类，而在这个模型中，79%的实际违约被预测为非违约。这种以增加假阳性为代价减少假阴性的权衡，若我们的成本估算准确，可能是可以接受的。

理解分类规则

分类规则以逻辑的 if-else 语句形式表示知识，赋予未标记示例一个类别。它们通过前件和后件来指定；这些形成一个假设，声明“如果发生这个，那么就会发生那个”。一个简单的规则可能是，“如果硬盘发出点击声，那么它即将发生故障。”前件包含某些特征值的组合，而后件则指定当规则条件满足时，应该分配的类别值。

规则学习者通常以类似于决策树学习者的方式使用。像决策树一样，它们可以用于生成未来行动的知识应用，如：

• 确定导致机械设备硬件故障的条件

• 描述群体的关键特征用于客户细分

• 寻找股票市场上股价大幅下跌或上涨之前的条件

另一方面，规则学习者在某些任务中相比树形结构提供了一些独特的优势。与必须通过一系列决策从上到下应用的树不同，规则是一种命题，可以像陈述事实一样被读取。此外，出于稍后会讨论的原因，规则学习者的结果可能比基于相同数据构建的决策树更简单、直接且更容易理解。

提示

正如你将在本章后面看到的，规则可以通过决策树生成。那么，为什么还要使用一个单独的规则学习算法呢？原因在于，决策树给任务带来了一些特定的偏差，而规则学习者通过直接识别规则来避免这些偏差。

规则学习者通常应用于特征主要或完全为名义型的数据问题。即使罕见事件仅在特征值之间的某种特定交互作用下发生，它们也能很好地识别这些罕见事件。

分而治之

分类规则学习算法利用一种叫做分而治之的启发式方法。这个过程包括识别一个覆盖训练数据子集的规则，然后将这个子集与剩余数据分离。当规则被添加时，数据的其他子集也会被分离，直到整个数据集被覆盖，且没有更多的例子留下。

想象规则学习过程的一种方式是通过创建越来越具体的规则来逐步深入数据，以识别类别值。假设你的任务是创建规则来判断一个动物是否是哺乳动物。你可以将所有动物的集合描绘成一个大空间，如下图所示：

规则学习者首先利用可用的特征来找到同质的群体。例如，使用一个表示物种是通过陆地、海洋还是空中移动的特征，第一个规则可能建议所有陆地动物都是哺乳动物：

你注意到这个规则有什么问题吗？如果你是一个动物爱好者，你可能已经意识到青蛙是两栖动物，而不是哺乳动物。因此，我们的规则需要更具体一点。我们可以进一步深入，假设哺乳动物必须在陆地上行走并且有尾巴：

可以定义一个额外的规则来分离蝙蝠，它是唯一剩下的哺乳动物。因此，这个子集可以与其他数据分开。

可以定义一个额外的规则来区分蝙蝠，它是唯一剩下的哺乳动物。一个可能的特征是蝙蝠与其他剩余动物的区别在于它们有毛发。通过使用基于这个特征的规则，我们就正确地识别了所有的动物：

到这一点时，由于所有训练实例都已经被分类，规则学习过程将停止。我们总共学到了三个规则：

• 会走在陆地上并且有尾巴的动物是哺乳动物。

• 如果动物没有毛发，它就不是哺乳动物。

• 否则，动物就是哺乳动物。

上面的例子说明了规则如何逐渐消耗更大、更大的数据片段，最终对所有实例进行分类。

由于规则似乎覆盖了数据的部分内容，分离并征服算法也被称为覆盖算法，而由此产生的规则则被称为覆盖规则。在下一节中，我们将通过研究一个简单的规则学习算法来了解覆盖规则是如何在实践中应用的。然后我们将研究一个更复杂的规则学习者，并将这两者应用于实际问题。

1R 算法

假设有一个电视游戏节目，节目中有一个轮盘，轮盘上有十个大小均等的彩色分区。三个分区是红色的，三个是蓝色的，四个是白色的。在转动轮盘之前，你需要选择其中一种颜色。当轮盘停止时，如果显示的颜色与你的预测匹配，你将赢得一大笔现金奖励。你应该选择哪种颜色？

如果你选择白色，当然更有可能赢得奖品——这是轮盘上最常见的颜色。显然，这个游戏节目有点荒谬，但它展示了最简单的分类器ZeroR，一个实际上不学习任何规则的规则学习器（因此得名）。对于每一个未标记的样本，无论其特征值如何，都会预测最常见的类别。

1R 算法（单规则或OneR）通过选择一个规则改进了 ZeroR。尽管这看起来可能过于简化，但它往往比你想象的表现更好。正如实证研究所示，对于许多实际任务，这个算法的准确率可以接近更复杂算法的表现。

注意

对 1R 算法出乎意料的表现进行深入了解，请参见 Holte RC. 非常简单的分类规则在大多数常用数据集上表现良好。《机器学习》1993 年；11:63-91。

1R 算法的优缺点如下表所示：

优点	缺点

|

• 生成一个简单易懂、可读性强的经验法则

• 经常表现得出奇的好

• 可以作为更复杂算法的基准

|

• 仅使用一个特征

• 可能过于简化

|

这个算法的工作原理很简单。对于每个特征，1R 将数据根据特征的相似值划分成组。然后，对于每个子集，算法预测多数类。计算基于每个特征的规则的错误率，选择错误最少的规则作为最终的单一规则。

以下表格展示了这个算法如何作用于我们在本节中早些时候看到的动物数据：

对于通过何种方式移动特征，数据集被划分为三组：空中、陆地和海洋。空中和海洋组的动物被预测为非哺乳动物，而陆地组的动物则被预测为哺乳动物。这导致了两个错误：蝙蝠和青蛙。是否有毛发特征将动物分为两组。有毛发的被预测为哺乳动物，而没有毛发的则被预测为非哺乳动物。统计了三个错误：猪、大象和犀牛。由于通过何种方式移动特征导致的错误较少，1R 算法将基于通过何种方式移动返回以下“单一规则”：

• 如果动物通过空气移动，它不是哺乳动物

• 如果动物通过陆地移动，它是哺乳动物

• 如果动物通过海洋移动，它不是哺乳动物

算法在这里停止，已找到最重要的单一规则。

显然，这个规则学习算法对于某些任务可能过于基础。你希望医疗诊断系统只考虑一个症状，还是希望自动驾驶系统仅根据一个因素来决定停车或加速？对于这些类型的任务，可能需要更复杂的规则学习器。我们将在接下来的章节中了解一个。

RIPPER 算法

早期的规则学习算法存在一些问题。首先，它们以速度慢而闻名，这使得它们在处理日益增多的大型数据集时效果不佳。其次，它们在噪声数据上通常容易出现不准确的情况。

解决这些问题的第一步由 Johannes Furnkranz 和 Gerhard Widmer 于 1994 年提出。他们的增量减少误差修剪（IREP）算法结合了预修剪和后修剪方法，这些方法使得规则变得非常复杂，然后在从完整数据集中分离实例之前进行修剪。尽管这一策略提高了规则学习者的性能，但决策树通常仍然表现得更好。

注意

关于 IREP 的更多信息，参见 Furnkranz J, Widmer G. 增量减少误差修剪。1994 年《第 11 届国际机器学习大会论文集》：70-77。

规则学习算法在 1995 年迈出了新的一步，当时 William W. Cohen 提出了重复增量修剪以产生误差减少（RIPPER）算法，该算法在 IREP 的基础上进行了改进，生成的规则能够匹配或超越决策树的性能。

注意

关于 RIPPER 的更多细节，参见 Cohen WW. 快速有效的规则归纳。1995 年《第 12 届国际机器学习大会论文集》：115-123。

如下表所示，RIPPER 的优缺点与决策树基本相当。主要的优点是，它们可能生成一个稍微更简洁的模型：

优势	劣势

|

• 生成易于理解的人类可读规则

• 在大型和噪声数据集上效率高

• 通常比可比的决策树生成更简单的模型

|

• 可能会生成似乎违背常识或专家知识的规则

• 不适合处理数值数据

• 可能不如更复杂的模型表现得好

|

RIPPER 算法是从多个规则学习算法的迭代中发展而来的，它是规则学习的高效启发式算法的拼凑。由于其复杂性，技术实现的详细讨论超出了本书的范围。然而，可以将其大致理解为一个三步过程：

1. 增长

2. 修剪

3. 优化

增长阶段使用分离与征服技术，贪婪地向规则中添加条件，直到它完美分类数据子集或没有更多的属性可用于拆分。与决策树类似，信息增益标准用于确定下一个拆分属性。当增加规则的特异性不再减少熵时，规则会立即被修剪。步骤一和步骤二会反复进行，直到达到停止标准，此时使用各种启发式方法优化整个规则集。

RIPPER 算法能够创建比 1R 算法更复杂的规则，因为它可以考虑多个特征。这意味着它可以创建具有多个前提的规则，例如“如果动物会飞并且有毛发，那么它是哺乳动物”。这提高了算法处理复杂数据的能力，但就像决策树一样，这也意味着规则可能会迅速变得更难理解。

注意

分类规则学习者的发展并未止步于 RIPPER。新的规则学习算法正在快速提出。文献调查显示了诸如 IREP++、SLIPPER、TRIPPER 等多种算法。

来自决策树的规则

分类规则也可以直接从决策树中获得。从一个叶节点开始，沿着分支回到根节点，你将获得一系列决策。这些可以组合成一条规则。下图展示了如何从决策树构建规则来预测电影的成功：

从根节点到每个叶子节点的路径，规则将是：

1. 如果明星数量较少，那么电影将会是票房惨败。

2. 如果明星数量较多且预算较高，那么电影将会是主流热片。

3. 如果明星数量较多且预算较低，那么电影将会是口碑成功。

由于接下来的部分将会阐明的原因，使用决策树生成规则的主要缺点是，结果的规则通常比规则学习算法学到的规则更复杂。决策树采用的分治策略与规则学习者的偏差不同。另一方面，有时从树中生成规则在计算上更加高效。

提示

C5.0()函数在C50包中会生成一个使用分类规则的模型，前提是你在训练模型时指定rules = TRUE。

是什么让树和规则显得贪心？

决策树和规则学习者被称为贪心学习者，因为它们按先到先得的方式使用数据。决策树使用的分治策略和规则学习者使用的分离策略都试图一次做出一个分割，首先找到最同质的分割，然后是下一个最好的分割，依此类推，直到所有实例都被分类。

贪婪方法的缺点是，贪婪算法并不保证为特定数据集生成最佳、最准确或最少数量的规则。通过提前采摘低悬果实，贪婪学习者可能会迅速找到一个对某个数据子集准确的单一规则；然而，采取这种做法，学习者可能会错失开发一个更加细致的规则集的机会，而这个规则集在整个数据集上具有更好的整体准确性。然而，如果不使用贪婪方法进行规则学习，那么对于除了最小的数据集以外的所有情况，规则学习将变得计算上不可行。

尽管树和规则都采用贪婪学习启发式，但它们在构建规则的方式上存在微妙的差异。也许区分它们最好的方式是注意到，一旦分治法在某个特征上进行分裂，分裂产生的分区就不能被重新征服，只能进一步细分。通过这种方式，树会永久地受到过去决策历史的限制。相反，一旦分离与征服法找到一个规则，任何不被所有规则条件覆盖的示例都可以被重新征服。

为了说明这一对比，考虑前面的例子，其中我们构建了一个规则学习器来判断一种动物是否是哺乳动物。规则学习器识别了三个规则，能够完美地分类示例动物：

• 上陆行走且有尾巴的动物是哺乳动物（熊、猫、狗、大象、猪、兔子、老鼠、犀牛）

• 如果动物没有毛发，那么它就不是哺乳动物（鸟类、电鳗、鱼类、青蛙、昆虫、鲨鱼）

• 否则，该动物是哺乳动物（蝙蝠）

相反，基于相同数据构建的决策树可能会提出四个规则来实现相同的完美分类：

• 如果动物上陆行走且有毛发，那么它是哺乳动物（熊、猫、狗、大象、猪、兔子、老鼠、犀牛）

• 如果一种动物上陆行走且没有毛发，那么它就不是哺乳动物（青蛙）

• 如果动物不上陆行走且有毛发，那么它是哺乳动物（蝙蝠）

• 如果动物不上陆行走且没有毛发，那么它就不是哺乳动物（鸟类、昆虫、鲨鱼、鱼类、电鳗）

这两种方法产生不同结果的原因在于青蛙在“上陆行走”决策后发生的变化。规则学习器允许青蛙被“没有毛发”决策重新分类，而决策树无法修改现有的分区，因此必须将青蛙归入自己的规则中。

一方面，由于规则学习器可以重新审视那些曾被考虑但最终未被先前规则覆盖的案例，规则学习器通常会找到比决策树生成的规则集更加简洁的规则集。另一方面，这种数据的重复使用意味着规则学习器的计算成本可能比决策树略高。

示例——使用规则学习器识别有毒蘑菇

每年，许多人因食用有毒野生蘑菇而生病，有时甚至死亡。由于许多蘑菇在外观上非常相似，偶尔甚至经验丰富的蘑菇采摘者也会中毒。

与识别有毒植物（如毒橡树或常春藤）不同，识别野生蘑菇是否有毒或可食用并没有像“叶子三片，留它们”这样的明确规则。更复杂的是，许多传统规则，如“有毒蘑菇是鲜艳的颜色”，提供了危险或误导性的信息。如果有简单、明确且一致的规则来识别有毒蘑菇，它们可能会拯救采蘑菇者的生命。

由于规则学习算法的一个优势是它们生成易于理解的规则，因此它们似乎非常适合这项分类任务。然而，这些规则的有效性将取决于它们的准确性。

第 1 步——收集数据

为了识别区分有毒蘑菇的规则，我们将使用卡内基梅隆大学 Jeff Schlimmer 的蘑菇数据集。原始数据集可以在 UCI 机器学习库中免费获得（archive.ics.uci.edu/ml）。

数据集包含来自《北美蘑菇观鸟指南》（1981 年版）中 23 种带柄蘑菇的 8,124 个蘑菇样本的信息。在该指南中，每种蘑菇的种类都被标定为“肯定可食用”，“肯定有毒”或“可能有毒，不建议食用”。对于该数据集，后者与“肯定有毒”类别合并，形成了两个类别：有毒和无毒。UCI 网站上的数据字典描述了蘑菇样本的 22 个特征，包括盖形、盖色、气味、鳃的大小和颜色、柄形状以及栖息地等特征。

提示

本章使用的是经过略微修改的蘑菇数据。如果你打算跟着例子一起操作，请从 Packt Publishing 网站下载mushrooms.csv文件并将其保存到你的 R 工作目录中。

第 2 步——探索并准备数据

我们首先使用read.csv()导入数据进行分析。由于所有 22 个特征和目标类别都是名义型的，在本例中，我们将设置stringsAsFactors = TRUE并利用自动因子转换功能：

> mushrooms <- read.csv("mushrooms.csv", stringsAsFactors = TRUE)

str(mushrooms)命令的输出显示数据包含 8,124 个观察值，涵盖 23 个变量，正如数据字典所描述的那样。虽然str()的输出大部分没有特别之处，但有一个特征值得一提。你是否注意到以下行中的veil_type变量有什么不同之处？

$ veil_type : Factor w/ 1 level "partial": 1 1 1 1 1 1 ...

如果你觉得一个因子只有一个水平很奇怪，那么你是对的。数据字典列出了该特征的两个水平：部分和普遍。然而，我们数据中的所有示例都被归类为部分。很可能这个数据元素被错误地编码了。无论如何，由于面纱类型在样本间没有变化，它不会为预测提供任何有用信息。我们将通过以下命令从分析中删除该变量：

> mushrooms$veil_type <- NULL

通过将NULL赋值给面纱类型向量，R 会从mushrooms数据框中删除该特征。

在深入分析之前，我们应该快速查看一下数据集中蘑菇type类别变量的分布：

> table(mushrooms$type)
 edible poisonous 
 4208      3916

大约 52%的蘑菇样本（N = 4,208）是可食用的，而 48%（N = 3,916）是有毒的。

为了本实验的目的，我们将蘑菇数据集中的 8,214 个样本视为所有可能野生蘑菇的完整集合。这是一个重要的假设，因为这意味着我们不需要将一些样本从训练数据中排除用于测试。我们不是试图开发涵盖未知蘑菇类型的规则；我们只是试图找到准确描述已知蘑菇类型完整集合的规则。因此，我们可以在相同的数据上构建并测试模型。

步骤 3 – 在数据上训练模型

如果我们在这些数据上训练一个假设的 ZeroR 分类器，它会预测什么？由于 ZeroR 忽略所有特征，仅仅预测目标的众数，用简单的话来说，它的规则会说所有蘑菇都是可食用的。显然，这不是一个很有用的分类器，因为它会让采蘑菇的人生病或死亡，因为几乎一半的蘑菇样本有可能是有毒的。我们的规则必须比这个做得好得多，才能提供可以发布的安全建议。同时，我们需要简单易记的规则。

由于简单规则通常具有极强的预测性，让我们看看一个非常简单的规则学习器在蘑菇数据上的表现。最后，我们将应用 1R 分类器，它将识别目标类别中最具预测性的单一特征，并使用它来构建一套规则。

我们将使用RWeka包中的 1R 实现，名为OneR()。你可能记得我们在第一章，机器学习介绍中，作为安装和加载包的教程一部分，安装了RWeka。如果你没有按照这些说明安装该包，你将需要使用install.packages("RWeka")命令，并确保你的系统上安装了 Java（参阅安装说明获取更多详情）。完成这些步骤后，输入library(RWeka)加载该包：

OneR()实现使用 R 公式语法来指定要训练的模型。公式语法使用~运算符（称为波浪号）来表达目标变量与其预测变量之间的关系。要学习的类别变量放在波浪号的左侧，预测特征写在右侧，用+运算符分隔。如果你想建模y类与预测变量x1和x2之间的关系，你可以写成y ~ x1 + x2。如果你想在模型中包含所有变量，可以使用特殊术语.。例如，y ~ .指定了y与数据集中的所有其他特征之间的关系。

提示

R 公式语法在许多 R 函数中被广泛使用，并提供了一些强大的功能来描述预测变量之间的关系。我们将在后面的章节中探索其中的一些特性。如果你急于了解，可以使用?formula命令查看文档。

使用type ~ .公式，我们将允许我们的第一个OneR()规则学习器在构建预测类型的规则时，考虑蘑菇数据中的所有可能特征：

> mushroom_1R <- OneR(type ~ ., data = mushrooms)

要检查它创建的规则，我们可以输入分类器对象的名称，在这个例子中是mushroom_1R：

> mushroom_1R

odor:
 almond  -> edible
 anise  -> edible
 creosote  -> poisonous
 fishy  -> poisonous
 foul  -> poisonous
 musty  -> poisonous
 none  -> edible
 pungent  -> poisonous
 spicy  -> poisonous
(8004/8124 instances correct)

在输出的第一行，我们看到选中了“气味”特征用于规则生成。气味的类别，如杏仁、茴香等，指定了蘑菇是否可能是可食用或有毒的规则。例如，如果蘑菇闻起来有腥味、腐臭味、霉味、刺激味、辣味或柏油味，那么蘑菇很可能是有毒的。另一方面，闻起来更宜人的气味，如杏仁和茴香，或者没有气味的蘑菇，通常被预测为可食用。对于蘑菇采集的野外指南来说，这些规则可以总结为一个简单的经验法则：“如果蘑菇闻起来不舒服，那么它很可能是有毒的。”

第 4 步 – 评估模型性能

输出的最后一行指出，规则正确预测了 8,124 个蘑菇样本中 8,004 个的可食性，约占 99％。我们可以使用summary()函数获得分类器的更多细节，如下例所示：

> summary(mushroom_1R)

=== Summary ===
Correctly Classified Instances        8004  98.5229 %
Incorrectly Classified Instances       120  1.4771 %
Kappa statistic                          0.9704
Mean absolute error                      0.0148
Root mean squared error                  0.1215
Relative absolute error                  2.958  %
Root relative squared error             24.323  %
Coverage of cases (0.95 level)          98.5229 %
Mean rel. region size (0.95 level)      50      %
Total Number of Instances             8124 

=== Confusion Matrix ===
 a    b   <-- classified as
 4208    0 |    a = edible
 120 3796 |    b = poisonous

标记为Summary的部分列出了多种不同的方式来衡量我们 1R 分类器的性能。我们将在第十章中详细介绍其中许多统计量，评估模型性能，因此现在暂时不讨论这些内容。

标有混淆矩阵的部分与之前使用的相似。在这里，我们可以看到规则出现错误的地方。右侧显示了关键，其中a = 可食用，b = 有毒。表格的列表示预测的蘑菇类别，行则将 4,208 个可食用蘑菇与 3,916 个有毒蘑菇分开。通过检查表格，我们可以看到，尽管 1R 分类器没有将任何可食用蘑菇错误分类为有毒蘑菇，但它却将 120 个有毒蘑菇错误分类为可食用的——这可是一个非常危险的错误！

考虑到学习者仅使用了一个特征，模型表现得相当不错；如果在采摘蘑菇时避免不愉快的气味，他们几乎可以避免去医院的风险。尽管如此，当涉及到生命时，“接近”并不足够，更不用说当读者生病时，野外指南的出版商可能会因为面临诉讼而不高兴。让我们看看能否增加一些规则，开发出一个更好的分类器。

步骤 5 – 改进模型性能

对于更复杂的规则学习器，我们将使用JRip()，它是一个基于 Java 的 RIPPER 规则学习算法实现。与之前使用的 1R 实现一样，JRip()包含在RWeka包中。如果尚未安装，确保使用library(RWeka)命令加载该包：

如语法框所示，训练JRip()模型的过程与我们之前训练OneR()模型的方式非常相似。这是RWeka包的一个优点；不同算法的语法一致，使得比较多个不同模型的过程变得非常简单。

让我们像训练OneR()模型一样训练JRip()规则学习器，让它从所有可用的特征中选择规则：

> mushroom_JRip <- JRip(type ~ ., data = mushrooms)

要检查规则，输入分类器的名称：

> mushroom_JRip

JRIP rules:
===========
(odor = foul) => type=poisonous (2160.0/0.0)
(gill_size = narrow) and (gill_color = buff) => type=poisonous (1152.0/0.0)
(gill_size = narrow) and (odor = pungent) => type=poisonous (256.0/0.0)
(odor = creosote) => type=poisonous (192.0/0.0)
(spore_print_color = green) => type=poisonous (72.0/0.0)
(stalk_surface_below_ring = scaly) and (stalk_surface_above_ring = silky) => type=poisonous (68.0/0.0)
(habitat = leaves) and (cap_color = white) => type=poisonous (8.0/0.0)
(stalk_color_above_ring = yellow) => type=poisonous (8.0/0.0)
 => type=edible (4208.0/0.0)
Number of Rules : 9

JRip()分类器从蘑菇数据中学习了共九条规则。可以将这些规则看作是一系列的 if-else 语句，类似于编程逻辑。前面三条规则可以表示为：

• 如果气味刺鼻，则蘑菇类型为有毒

• 如果鳃的大小狭窄且鳃的颜色为浅棕色，则蘑菇类型为有毒

• 如果鳃的大小狭窄且气味刺鼻，则蘑菇类型为有毒

最后，第九条规则意味着任何未被前面八条规则涵盖的蘑菇样本都是可食用的。按照我们编程逻辑的示例，可以这样理解：

• 否则，蘑菇是可食用的

每条规则旁边的数字表示该规则涵盖的实例数量以及误分类实例的数量。值得注意的是，使用这九条规则时，没有任何蘑菇样本被误分类。因此，最后一条规则涵盖的实例数量正好等于数据中可食用蘑菇的数量（N = 4,208）。

以下图表大致说明了规则如何应用于蘑菇数据。如果你把椭圆内的所有内容想象为所有种类的蘑菇，规则学习器识别了特征或特征集，将同质的群体从更大的群体中分离出来。首先，算法找到了一大群具有独特恶臭味的有毒蘑菇。接着，它找到了更小、更具体的有毒蘑菇群体。通过识别每种有毒蘑菇的覆盖规则，剩余的所有蘑菇都被确定为可食用。感谢大自然，每种蘑菇的特性足够独特，以至于分类器能够实现 100%的准确率。

总结

本章介绍了两种分类方法，它们使用所谓的“贪婪”算法，根据特征值将数据进行划分。决策树采用分而治之的策略，创建类似流程图的结构，而规则学习器则通过分离和征服数据，识别逻辑的 if-else 规则。两种方法生成的模型都可以在没有统计背景的情况下进行解释。

一种流行且高度可配置的决策树算法是 C5.0。我们使用 C5.0 算法创建了一个决策树来预测贷款申请人是否会违约。通过使用提升和成本敏感错误的选项，我们能够提高准确率，避免那些会使银行蒙受更大损失的高风险贷款。

我们还使用了两种规则学习器，1R 和 RIPPER，来开发规则以识别有毒蘑菇。1R 算法使用一个特征，在识别潜在致命的蘑菇样本时达到了 99%的准确率。另一方面，使用更复杂的 RIPPER 算法生成的九条规则正确地识别了每个蘑菇的可食性。

本章仅仅触及了树和规则应用的表面。在第六章，数值数据预测—回归方法中，我们将学习回归树和模型树等技术，它们使用决策树进行数值预测而非分类。在第十一章，提升模型性能中，我们将探讨如何通过将决策树组合成一个叫做随机森林的模型来提升其性能。在第八章，寻找模式—使用关联规则的市场篮子分析中，我们将看到如何利用关联规则（分类规则的一种相关形式）来识别交易数据中的项目组。

第六章：数值数据预测——回归方法

数学关系帮助我们理解日常生活中的许多方面。例如，体重是摄入卡路里的函数，收入通常与教育和工作经验有关，民意调查数据帮助我们估计总统候选人连任的几率。

当这些关系用精确的数字表示时，我们获得了更多的清晰度。例如，每天额外摄入 250 千卡的热量可能导致每月增加近 1 公斤体重；每增加一年的工作经验，年薪可能增加$1,000；而在经济强劲时，总统连任的可能性更大。显然，这些方程式并不完美适用于每种情况，但我们可以预期它们在平均情况下是相对正确的。

本章通过超越前面所述的分类方法，引入了估计数值数据之间关系的技术，从而扩展了我们的机器学习工具包。在研究几个现实世界的数值预测任务时，您将学习到：

• 回归中使用的基本统计原理，一种模拟数值关系的大小和强度的技术

• 如何准备回归分析的数据，并估算和解读回归模型

• 一对混合技术，称为回归树和模型树，它们将决策树分类器适应于数值预测任务

基于统计学领域的大量研究，本章所用的方法在数学方面略重于之前所讲的内容，但别担心！即使你的代数技能有些生疏，R 语言会帮你完成繁重的计算。

理解回归

回归分析关注的是指定单一数值因变量（即要预测的值）与一个或多个数值自变量（即预测变量）之间的关系。顾名思义，因变量取决于自变量的值。回归的最简单形式假设自变量与因变量之间的关系呈直线。

注意

“回归”一词用于描述将线条拟合到数据的过程，源自 19 世纪末弗朗西斯·高尔顿爵士在遗传学研究中的发现。他发现，极端矮小或极端高大的父亲，往往有身高更接近平均值的儿子。他将这一现象称为“回归到均值”。

你可能还记得从基础代数中，直线可以用斜率-截距形式来定义，类似于 y = a + bx。在这种形式中，字母y表示因变量，x表示自变量。斜率项b指定了直线在x每增加一个单位时上升的量。正值定义了向上倾斜的直线，而负值则定义了向下倾斜的直线。项a被称为截距，因为它指定了直线交叉或截取垂直y轴的点。它表示当x = 0时y的值。

回归方程使用类似于斜率-截距形式的数据模型。机器的任务是确定* a 和 b 的值，使得指定的直线能最好地将提供的x值与y*值关联起来。可能并不总是存在一个完美关联这些值的单一函数，因此机器还必须有某种方式来量化误差范围。我们稍后会深入讨论这一点。

回归分析通常用于建模数据元素之间的复杂关系，估计处理对结果的影响，并进行未来的外推。尽管它可以应用于几乎任何任务，但一些特定的应用案例包括：

• 检验群体和个体在其测量特征上的变异性，广泛应用于经济学、社会学、心理学、物理学和生态学等多个学科的科学研究。

• 定量分析事件与响应之间的因果关系，例如临床药物试验、工程安全测试或市场调研中的因果关系。

• 识别可以用来根据已知标准预测未来行为的模式，如预测保险理赔、自然灾害损失、选举结果和犯罪率等。

回归方法也用于统计假设检验，该方法通过观察到的数据来判断一个前提是否可能为真或为假。回归模型对关系的强度和一致性的估计提供了可以用来评估观察结果是否仅仅由偶然因素造成的信息。

注意

假设检验是极其复杂的，超出了机器学习的范畴。如果你对这个话题感兴趣，入门级统计学教材是一个很好的起点。

回归分析并不等同于单一算法。相反，它是一个涵盖大量方法的总称，这些方法可以适应几乎任何机器学习任务。如果你只能选择一种方法，回归方法是一个不错的选择。有人可以将一生献给这一领域，或许仍然有许多要学习的内容。

本章将仅关注最基本的线性回归模型——那些使用直线的模型。当只有一个自变量时，称为简单线性回归。当有两个或更多自变量时，称为多元线性回归，或简称为“多元回归”。这两种技术都假设因变量是连续量度的。

回归还可以用于其他类型的因变量，甚至一些分类任务。例如，逻辑回归用于建模二元分类结果，而泊松回归——以法国数学家西门·泊松命名——用于建模整数计数数据。被称为多项式逻辑回归的方法则用于建模分类结果，因此它可以用于分类。所有回归方法遵循相同的基本原理，因此在理解了线性回归后，学习其他方法相对简单。

提示

许多专业的回归方法属于广义线性模型（GLM）类。使用 GLM 时，线性模型可以通过使用链接函数来推广到其他模式，这些函数为x和y之间的关系指定了更复杂的形式。这样，回归就可以应用于几乎任何类型的数据。

我们将从简单线性回归的基本情况开始。尽管名字中有“简单”二字，这种方法并不简单，仍能解决复杂问题。在下一节中，我们将看到如何通过使用简单线性回归模型来避免一次悲剧性的工程灾难。

简单线性回归

1986 年 1 月 28 日，美国航天飞机挑战者号的七名机组成员在火箭助推器发生故障后丧生，导致灾难性的解体。在事后，专家们将发射温度视为潜在的罪魁祸首。负责密封火箭接头的橡胶 O 型环从未在 40ºF（4ºC）以下的温度下进行过测试，而发射当天的天气异常寒冷，气温低于冰点。

从后见之明来看，这起事故成为了数据分析和可视化重要性的案例研究。尽管目前不清楚火箭工程师和决策者在发射前掌握了哪些信息，但不可否认的是，如果有更好的数据并加以谨慎使用，很可能能够避免这场灾难。

注意

本节分析基于 Dalal SR, Fowlkes EB, Hoadley B. 航天飞机的风险分析：挑战者号发射前的故障预测。美国统计学会学报，1989 年；84：945-957。有关数据如何改变结果的一个视角，参见 Tufte ER. 视觉解释：图像与数量、证据与叙述。Graphics Press，1997 年。一个反观点请参见 Robison W, Boisioly R, Hoeker D, Young, S. 代表与误代表：Tufte 与 Morton Thiokol 工程师对挑战者号的看法。科学与工程伦理，2002 年；8：59-81。

火箭工程师几乎肯定知道，低温可能使部件变脆，无法正确密封，从而导致更高的危险燃料泄漏的可能性。然而，考虑到继续发射的政治压力，他们需要数据来支持这一假设。一个能够展示温度与 O 型环故障之间关联的回归模型，并能够根据预期的发射温度预测故障的可能性，可能会非常有帮助。

为了构建回归模型，科学家们可能使用了来自 23 次成功航天飞机发射的数据，包括发射温度和部件故障数据。部件故障表示两种问题之一。第一个问题叫做侵蚀，当过热烧坏 O 型环时就会发生这种情况。第二个问题叫做泄漏，当热气体通过或“冲过”密封不良的 O 型环时，就会发生泄漏。由于航天飞机总共有六个主要 O 型环，因此每次飞行最多可能发生六个故障事件。尽管火箭可以在发生一个或多个故障事件的情况下生还，或者在发生一个故障事件的情况下就失败，但每增加一个故障事件，灾难性失败的概率就会增加。

以下散点图展示了前 23 次发射中检测到的主要 O 型环故障与发射温度的关系：

观察图表，可以发现一个明显的趋势。发生在较高温度下的发射，O 型环故障事件较少。

此外，最冷的发射温度（53º F）发生了两次故障事件，这是在另外一次发射中才达到的水平。掌握了这些信息后，可以看出挑战者号计划在比这低 20 多度的温度下发射，这似乎令人担忧。那么我们到底应该多担心呢？为了回答这个问题，我们可以借助简单线性回归。

简单线性回归模型通过使用一个由方程定义的直线，定义了一个因变量与单一自变量之间的关系，方程形式如下：

不要被希腊字符吓到，这个方程仍然可以通过之前描述的斜截式形式来理解。截距，α（阿尔法），描述了直线与y轴的交点，而斜率，β（贝塔），描述了当x增加时，y的变化。对于航天飞机发射数据，斜率告诉我们每升高 1 度发射温度，O 型环故障的预期减少量。

提示

希腊字符常常在统计学领域中用来表示统计函数的参数变量。因此，进行回归分析时，需要找到参数估计值，即α和β的估计值。阿尔法和贝塔的参数估计值通常用a和b表示，尽管你可能会发现这些术语和符号有时可以互换使用。

假设我们知道航天飞机发射数据的回归方程中估计的回归参数为：a = 3.70和b = -0.048。

因此，完整的线性方程是y = 3.70 – 0.048x。暂时忽略这些数字是如何得出的，我们可以像这样将直线绘制在散点图上：

如图所示，在华氏 60 度时，我们预测 O 型环故障略低于 1 个；在华氏 70 度时，我们预计约有 0.3 个故障。如果我们将模型外推到 31 度——挑战者号发射时的预测温度——我们预计将有大约3.70 - 0.048 * 31 = 2.21个 O 型环故障事件。假设每个 O 型环故障都有相同的概率导致灾难性的燃料泄漏，那么在 31 度时，挑战者号发射的风险几乎是 60 度典型发射的三倍，是 70 度发射的八倍以上。

注意，直线并不是精确通过每一个数据点的。相反，它大致穿过数据，部分预测值高于直线，部分低于直线。在下一部分中，我们将学习为什么选择了这条特定的直线。

普通最小二乘估计

为了确定α和β的最优估计值，使用了一种叫做普通最小二乘法（OLS）的估计方法。在 OLS 回归中，选择的斜率和截距是为了最小化平方误差的和，即预测的y值与实际的y值之间的垂直距离。这些误差被称为残差，并在以下图示中为若干点进行了说明：

用数学术语表达，OLS 回归的目标可以表示为最小化以下方程：

用通俗的话来说，这个方程将e（误差）定义为实际的y值与预测的y值之间的差异。误差值被平方并在所有数据点上求和。

提示

位于* y 项上方的插入符号（^）是统计符号中常用的功能。它表示该项是对真实 y 值的估计。这被称为 y *帽，发音就像你头上戴的帽子。

• a 的解依赖于 b *的值。可以使用以下公式来获得：

提示

为了理解这些方程式，您需要了解另一个统计符号。出现在* x 和 y 项上方的横杠表示 x 或 y 的平均值。这被称为 x 横杠或 y *横杠，发音就像你去喝酒的地方。

尽管证明超出了本书的范围，但可以使用微积分来证明，导致最小平方误差的* b *值为：

如果我们将这个方程拆解成其组成部分，我们可以稍微简化它。* b 的分母应该看起来很熟悉；它非常类似于 x 的方差，表示为 Var(x)。正如我们在第二章《管理与理解数据》中学到的，方差涉及找到 x *均值的平均平方偏差。可以表示为：

分子涉及将每个数据点从均值* x 值的偏差与该点从均值 y 值的偏差相乘。这类似于 x 和 y 的协方差函数，表示为 Cov(x, y)*。协方差公式为：

如果我们将协方差函数除以方差函数，* n 项将被取消，我们可以将 b *的公式重写为：

给定这个重述，通过使用内置的 R 函数，计算* b *的值变得很容易。让我们将其应用于火箭发射数据，来估算回归线。

提示

如果您想跟随这些示例，请从 Packt Publishing 网站下载challenger.csv文件，并使用launch <- read.csv("challenger.csv")命令加载到数据框中。

假设我们的航天飞机发射数据存储在一个名为launch的数据框中，独立变量* x 是温度，依赖变量 y 是distress_ct。然后我们可以使用 R 的cov()和var()函数来估算 b *：

> b <- cov(launch$temperature, launch$distress_ct) /
 var(launch$temperature)
> b
[1] -0.04753968

从这里我们可以使用mean()函数估算* a *：

> a <- mean(launch$distress_ct) - b * mean(launch$temperature)
> a
[1] 3.698413

手动估计回归方程并不理想，因此 R 提供了自动执行此计算的函数。我们将很快使用这些方法。首先，我们将通过学习一种测量线性关系强度的方法来扩展我们对回归的理解，然后我们将看到如何将线性回归应用于包含多个自变量的数据。

相关性

相关性是两个变量之间的一个数字，表示它们的关系与一条直线的贴合程度。通常，相关性指的是 Pearson 相关系数，它由 20 世纪的数学家 Karl Pearson 开发。相关性范围在 -1 到 +1 之间。极端值表示完美的线性关系，而接近零的相关性则表示缺乏线性关系。

以下公式定义了 Pearson 相关性：

提示

这里引入了更多的希腊符号。第一个符号（看起来像小写字母 p）是 rho，用于表示 Pearson 相关性统计量。那些看起来像横过的 q 字母是希腊字母 sigma，它们表示 x 或 y 的标准差。

使用这个公式，我们可以计算发射温度与 O 环密封件故障事件数量之间的相关性。回顾一下，协方差函数是cov()，标准差函数是sd()。我们将把结果存储在r中，这是一个常用于表示估计相关性的字母：

> r <- cov(launch$temperature, launch$distress_ct) /
 (sd(launch$temperature) * sd(launch$distress_ct))
> r
[1] -0.5111264

或者，我们可以使用 R 的相关性函数cor()：

> cor(launch$temperature, launch$distress_ct)
[1] -0.5111264

温度与受损 O 环数量之间的相关性为 -0.51。负相关意味着温度的升高与 O 环受损数量的减少相关联。对于研究 O 环数据的 NASA 工程师来说，这将是一个非常明显的指标，表明低温发射可能会存在问题。相关性还告诉我们温度与 O 环受损之间关系的相对强度。由于 -0.51 距离最大负相关 -1 还差一半，这意味着它们之间有着适度强的负线性关联。

有各种经验法则用于解释相关性的强度。一种方法是将 0.1 到 0.3 之间的值标为“弱”，0.3 到 0.5 之间的值标为“中等”，大于 0.5 的值标为“强”（这些也适用于负相关的类似范围）。然而，这些阈值对于某些目的来说可能太宽松了。通常，相关性必须结合上下文进行解释。对于涉及人类的数据，0.5 的相关性可能被认为是极高的，而对于由机械过程生成的数据，0.5 的相关性可能较弱。

提示

你可能听过“相关性不代表因果关系”这一说法。其根源在于，相关性仅描述了两个变量之间的关联，但可能存在其他未测量的解释。例如，死亡率与每天观看电影的时间之间可能有很强的关联，但在医生开始建议我们都多看电影之前，我们需要排除另一种解释——年轻人看更多的电影，并且死亡的可能性较小。

衡量两个变量之间的相关性为我们提供了一种快速评估独立变量和依赖变量之间关系的方法。随着我们开始使用更多预测变量定义回归模型，这一点将变得愈加重要。

多元线性回归

大多数实际分析都有多个独立变量。因此，你很可能会使用多元线性回归来进行大多数数值预测任务。多元线性回归的优缺点见下表：

优点	缺点

|

• 到目前为止，建模数值数据最常见的方法

• 可以适应几乎任何建模任务

• 提供特征与结果之间关系的强度和大小的估计

|

• 对数据做出了强假设

• 模型的形式必须由用户提前指定

• 不处理缺失数据

• 仅适用于数值特征，因此分类数据需要额外的处理

• 需要一定的统计学知识才能理解该模型

|

我们可以将多元回归理解为简单线性回归的扩展。两者的目标相似——寻找能够最小化线性方程预测误差的β系数值。关键的不同点在于，针对额外的独立变量，有了更多的项。

多元回归方程通常遵循以下方程的形式。依赖变量 y 被指定为截距项 α 与估计的 β 值与每个 i 特征的 x 值之积的和。这里加入了一个误差项（用希腊字母 epsilon 表示），以提醒我们预测并不完美。这代表了前面提到的残差项：

让我们暂时考虑一下估计的回归参数的解释。你会注意到，在前面的方程中，为每个特征提供了一个系数。这使得每个特征对 y 值的影响可以单独估计。换句话说，y 会因每个 x[i] 单位增加而改变 β[i] 的数量。截距 α 则是当所有独立变量为零时 y 的期望值。

由于截距项α与其他任何回归参数并无不同，它有时也被表示为β[0]（读作 beta 零），如以下公式所示：

就像之前一样，截距与任何自变量x并无直接关系。然而，出于很快会变得清晰的原因，帮助理解的方式是将*β[0]想象成与一个常数项x[0]相乘，其中x[0]*的值为 1：

为了估计回归参数的值，每个观测值的因变量y必须通过前述形式的回归方程与观测值的自变量x相关联。以下图展示了这种结构：

前面图示中展示的大量行列数据，可以用粗体字体的矩阵符号来简洁地表示，表明每个术语代表多个值：

现在，因变量是一个向量Y，每一行对应一个示例。自变量已经合并成一个矩阵X，其中每列对应一个特征，另外还有一列全为'1'的值，代表截距项。每列有一个对应示例的行。回归系数β和残差ε现在也都是向量。

目标是解出β，即回归系数向量，它最小化预测值与实际Y值之间的平方误差和。寻找最优解需要使用矩阵代数；因此，推导过程需要比本书所能提供的更为详细的关注。不过，如果你愿意相信他人的研究，β向量的最佳估计可以通过以下方式计算：

这个解决方案使用了一对矩阵运算——T表示矩阵X的转置，而负指数表示矩阵的逆。使用 R 的内置矩阵运算，我们可以实现一个简单的多元回归学习器。接下来，我们将这个公式应用于挑战者发射数据。

提示

如果你对前面的矩阵运算不熟悉，维基百科上关于转置和矩阵逆的页面提供了详细的介绍，即使没有强大的数学背景，也能很好地理解。

使用以下代码，我们可以创建一个名为reg()的基本回归函数，它接受一个参数y和一个参数x，并返回一个估计的 beta 系数向量：

reg <- function(y, x) {
 x <- as.matrix(x)
 x <- cbind(Intercept = 1, x)
 b <- solve(t(x) %*% x) %*% t(x) %*% y
 colnames(b) <- "estimate"
 print(b)
}

这里创建的 reg() 函数使用了我们之前没有用过的几个 R 命令。首先，由于我们将使用该函数处理数据框中的列集，as.matrix() 函数用于将数据框转换为矩阵形式。接下来，cbind() 函数用于将一个额外的列绑定到 x 矩阵上；命令 Intercept = 1 指示 R 将新列命名为 Intercept 并用重复的 1 填充该列。然后，对 x 和 y 对象执行一系列矩阵操作：

• solve() 求解矩阵的逆

• t() 用于转置矩阵

• %*% 用于两个矩阵的乘法

通过将这些组合在一起，我们的函数将返回一个向量 b，其中包含线性模型中与 x 相关的估计参数。函数中的最后两行给 b 向量命名，并将结果打印到屏幕上。

让我们将该函数应用于航天飞机发射数据。如以下代码所示，该数据集包括三个特征和故障计数（distress_ct），这是我们关心的结果：

> str(launch)
'data.frame':      23 obs. of  4 variables:
 $ distress_ct         : int  0 1 0 0 0 0 0 0 1 1 ...
 $ temperature         : int  66 70 69 68 67 72 73 70 57 63 ...
 $ field_check_pressure: int  50 50 50 50 50 50 100 100 200 ...
 $ flight_num          : int  1 2 3 4 5 6 7 8 9 10 ...

我们可以通过将结果与温度与 O 型环故障之间的简单线性回归模型进行比较来确认我们的函数是否正确。我们之前发现该模型的参数为 a = 3.70 和 b = -0.048。由于温度在发射数据的第三列中，我们可以按如下方式运行 reg() 函数：

> reg(y = launch$distress_ct, x = launch[2])
 estimate
Intercept    3.69841270
temperature -0.04753968

这些值与我们之前的结果完全一致，因此我们可以使用该函数来构建一个多元回归模型。我们将像之前一样应用它，不过这次我们指定三列数据，而不仅仅是其中的一列：

> reg(y = launch$distress_ct, x = launch[2:4])
 estimate
Intercept             3.527093383
temperature          -0.051385940
field_check_pressure  0.001757009
flight_num            0.014292843

该模型预测了温度、现场检查压力和发射 ID 号与 O 型环故障事件的关系。与简单线性回归模型类似，温度变量的系数为负，表明随着温度升高，预期的 O 型环故障事件数减少。现场检查压力是指发射前对 O 型环施加的压力。虽然检查压力最初是 50 psi，但在某些发射中，它被提高到 100 psi 和 200 psi，这使得一些人认为它可能是 O 型环腐蚀的原因。该系数为正，但较小。飞行编号被包含在内，以考虑航天飞机的使用年限。随着使用年限的增加，其部件可能变得更脆弱或更容易发生故障。飞行编号与故障次数之间的微弱正相关可能反映了这一事实。

到目前为止，我们只触及了线性回归建模的表面。虽然我们的工作有助于我们准确理解回归模型是如何构建的，但 R 的函数还包括一些额外的功能，必要时可用于更复杂的建模任务和诊断输出，这些对于帮助模型解释和评估拟合度是必需的。让我们将回归知识应用于一个更具挑战性的学习任务。

示例 – 使用线性回归预测医疗费用

为了使健康保险公司盈利，它需要收取的年度保费高于用于支付受益人医疗费用的支出。因此，保险公司在开发能够准确预测被保险人群体医疗费用的模型上投入了大量的时间和金钱。

医疗费用很难估算，因为最昂贵的疾病较为罕见且看似随机。然而，某些疾病在特定人群中更为常见。例如，吸烟者比非吸烟者更容易患上肺癌，而肥胖者可能更容易患上心脏病。

本次分析的目标是使用患者数据来估算特定人群的平均医疗费用。这些估算结果可以用于创建精算表，根据预期的治疗费用来调整年度保费价格。

第 1 步 – 收集数据

本次分析将使用一个模拟数据集，其中包含美国患者的假设医疗费用数据。这些数据是使用美国人口普查局的统计数据为本书创建的，因此大致反映了现实世界的情况。

提示

如果你希望互动地跟随本书的内容，可以从 Packt Publishing 网站下载insurance.csv文件，并将其保存在 R 的工作文件夹中。

insurance.csv文件包括 1,338 个目前已参加保险计划的受益人示例，文件中包含患者的特征以及该年度总医疗费用。特征包括：

• age：一个整数，表示主要受益人的年龄（不包括 64 岁以上的人群，因为他们通常由政府提供保险）。

• sex：保险持有人的性别，可以是男性或女性。

• bmi：身体质量指数（BMI），它表示一个人的体重与身高的关系。BMI 等于体重（以千克为单位）除以身高（以米为单位）的平方。理想的 BMI 范围是 18.5 到 24.9。

• children：一个整数，表示保险计划所涵盖的子女/受抚养人数量。

• smoker：一个是或否的分类变量，表示被保险人是否定期吸烟。

• region：受益人在美国的居住地，分为四个地理区域：东北、东南、西南或西北。

考虑这些变量如何与计费医疗费用相关是非常重要的。例如，我们可能预期，老年人和吸烟者面临较大医疗费用的风险。与许多其他机器学习方法不同，在回归分析中，特征之间的关系通常由用户指定，而不是自动检测。我们将在下一节探讨这些潜在的关系。

第 2 步 – 探索和准备数据

如我们之前所做的，我们将使用read.csv()函数加载数据进行分析。我们可以安全地使用stringsAsFactors = TRUE，因为将这三个名义变量转换为因子是合适的：

> insurance <- read.csv("insurance.csv", stringsAsFactors = TRUE)

str()函数确认数据格式符合我们的预期：

> str(insurance)
'data.frame':     1338 obs. of  7 variables:
 $ age     : int  19 18 28 33 32 31 46 37 37 60 ...
 $ sex     : Factor w/ 2 levels "female","male": 1 2 2 2 2 1 ...
 $ bmi     : num  27.9 33.8 33 22.7 28.9 25.7 33.4 27.7 ...
 $ children: int  0 1 3 0 0 0 1 3 2 0 ...
 $ smoker  : Factor w/ 2 levels "no","yes": 2 1 1 1 1 1 1 1 ...
 $ region  : Factor w/ 4 levels "northeast","northwest",..: ...
 $ expenses: num  16885 1726 4449 21984 3867 ...

我们模型的因变量是expenses，它表示每个人在一年内为保险计划支付的医疗费用。在构建回归模型之前，检查数据是否符合正态分布通常是有帮助的。尽管线性回归并不严格要求因变量正态分布，但当这一假设成立时，模型往往拟合得更好。让我们看一下汇总统计信息：

> summary(insurance$expenses)
 Min. 1st Qu.  Median    Mean 3rd Qu.    Max.
 1122    4740    9382   13270   16640   63770

因为均值大于中位数，这意味着保险费用的分布是右偏的。我们可以通过直方图来直观确认这一点：

> hist(insurance$expenses)

输出结果如下所示：

正如预期的那样，图表显示了右偏的分布。它还表明，尽管分布的尾部远超这些峰值，但我们数据中的大多数人年医疗费用在零到$15,000 之间。虽然这个分布并不适合线性回归，但提前了解这一弱点有助于我们设计出更适合的模型。

在我们解决这个问题之前，还有另一个问题需要处理。回归模型要求每个特征都是数值型的，但我们在数据框中有三个因子型特征。例如，性别变量分为男性和女性两个级别，而吸烟者分为是和否。从summary()输出中，我们知道region变量有四个级别，但我们需要进一步查看它们的分布：

> table(insurance$region)
northeast northwest southeast southwest
 324       325       364       325

在这里，我们看到数据几乎均匀地划分到四个地理区域中。稍后我们将看到 R 的线性回归函数如何处理这些因子变量。

探索特征之间的关系——相关矩阵

在拟合回归模型之前，确定自变量与因变量及相互之间的关系是很有帮助的。相关矩阵提供了这些关系的快速概览。给定一组变量，它为每一对关系提供相关性。

要为保险数据框中的四个数值变量创建相关矩阵，可以使用cor()命令：

> cor(insurance[c("age", "bmi", "children", "expenses")])
 age        bmi   children   expenses
age      1.0000000 0.10934101 0.04246900 0.29900819
bmi      0.1093410 1.00000000 0.01264471 0.19857626
children 0.0424690 0.01264471 1.00000000 0.06799823
expenses 0.2990082 0.19857626 0.06799823 1.00000000

在每一行和列的交点处，会列出由该行和列所表示的变量的相关性。对角线上的值始终是1.0000000，因为一个变量与它自身的相关性总是完美的。对角线上的值上下是对称的，因此相关性是对称的。换句话说，cor(x, y)等于cor(y, x)。

矩阵中的相关性都不算强，但有一些显著的关联。例如，age和bmi之间似乎存在弱正相关，意味着随着年龄的增长，体重通常会增加。age与expenses、bmi与expenses、以及children与expenses之间也有中等程度的正相关。这些关联意味着随着年龄、体重和孩子数量的增加，预期的保险费用也会上升。我们将在构建最终回归模型时更清晰地揭示这些关系。

可视化特征间的关系 – 散点图矩阵

使用散点图来可视化数值特征之间的关系也很有帮助。尽管我们可以为每一种可能的关系创建散点图，但对于大量特征来说，这样做可能会变得很繁琐。

另一种选择是创建一个散点图矩阵（有时简写为SPLOM），它仅仅是一个按网格排列的散点图集合。它用于检测三个或更多变量之间的模式。散点图矩阵并不是真正的多维可视化，因为每次仅检查两个特征。但它提供了一个大致的感知，显示数据可能的相互关系。

我们可以使用 R 的图形功能为四个数值特征（age、bmi、children和expenses）创建一个散点图矩阵。pairs()函数在默认的 R 安装中提供，能够生成基本的散点图矩阵功能。要调用该函数，只需提供要展示的数据框。这里，我们将insurance数据框限制为四个感兴趣的数值变量：

> pairs(insurance[c("age", "bmi", "children", "expenses")])

这会生成如下图表：

在散点图矩阵中，每一行和每一列的交点处展示的是由行和列所表示的变量的散点图。对角线以上和以下的图表是转置的，因为x轴和y轴已被交换。

你在这些图表中注意到任何模式吗？尽管有些看起来像是随机的点云，但有些似乎展示了一些趋势。age与expenses之间的关系显示了几条相对直的线，而bmi与expenses的图表则有两组明显不同的点。在其他图表中很难发现趋势。

如果我们在图表中添加更多信息，它将变得更加有用。可以使用psych包中的pairs.panels()函数创建一个增强版的散点图矩阵。如果你没有安装该包，请输入install.packages("psych")进行安装，并通过library(psych)命令加载它。然后，我们可以像之前那样创建一个散点图矩阵：

> pairs.panels(insurance[c("age", "bmi", "children", "expenses")])

这会生成一个稍微更具信息量的散点图矩阵，如下所示：

对角线以上的散点图已被相关矩阵替代。对角线上的直方图显示了每个特征值的分布情况。最后，对角线下方的散点图现在呈现了更多的视觉信息。

每个散点图上椭圆形的物体是一个相关椭圆。它提供了相关性强度的可视化。椭圆中心的点表示 x 和 y 轴变量的均值位置。变量之间的相关性通过椭圆的形状来表示；它越被拉伸，相关性越强。几乎完全圆形的椭圆，如 bmi 和 children 之间的关系，表明相关性非常弱（在这种情况下为 0.01）。

散点图上绘制的曲线称为局部加权回归曲线（loess curve）。它表示 x 和 y 轴变量之间的一般关系。通过示例最容易理解。age 和 children 的曲线呈倒 U 形，峰值出现在中年左右。这意味着样本中最年长和最年轻的人群在保险计划中有的孩子比中年人少。由于这一趋势是非线性的，仅凭相关性是无法得出这一结论的。另一方面，age 和 bmi 的 loess 曲线是逐渐向上的直线，表明随着年龄的增长，体重指数也在增加，但我们已经通过相关矩阵推断出这一点。

步骤 3 – 在数据上训练模型

要在 R 中拟合线性回归模型，可以使用 lm() 函数。该函数包含在 stats 包中，默认情况下应该随着 R 安装一起加载。lm() 的语法如下：

以下命令拟合一个线性回归模型，将六个自变量与总医疗费用相关联。R 的公式语法使用波浪号字符 ~ 来描述模型；因变量 expenses 位于波浪号的左侧，而自变量位于右侧，并由 + 符号分隔。由于回归模型的截距项默认假定存在，因此无需明确指定：

> ins_model <- lm(expenses ~ age + children + bmi + sex +
 smoker + region, data = insurance)

由于 . 字符可以用来指定所有特征（排除公式中已指定的特征），因此以下命令等同于之前的命令：

> ins_model <- lm(expenses ~ ., data = insurance)

构建模型后，只需输入模型对象的名称即可查看估计的贝塔系数：

> ins_model

Call:
lm(formula = expenses ~ ., data = insurance)

Coefficients:
 (Intercept)              age          sexmale 
 -11941.6            256.8           -131.4 
 bmi         children        smokeryes 
 339.3            475.7          23847.5 
regionnorthwest  regionsoutheast  regionsouthwest 
 -352.8          -1035.6           -959.3

理解回归系数是相当简单的。截距是当所有自变量为零时，expenses的预测值。正如这里的情况，截距通常单独没有太大意义，因为所有特征不可能都有零值。例如，由于不存在年龄为零且 BMI 为零的人，因此截距没有实际意义。出于这个原因，在实践中，截距通常会被忽略。

贝塔系数表示每增加一个特征值，假设其他所有值保持不变，医疗费用的预估增加。例如，每增加一年年龄，平均预计医疗费用增加$256.80，假设其他条件不变。类似地，每增加一个孩子，平均每年医疗费用增加$475.70，而每增加一个单位的 BMI，平均每年医疗费用增加$339.30，其他条件不变。

你可能会注意到，尽管我们在模型公式中只指定了六个特征，但报告的系数除了截距之外还有八个。这是因为lm()函数自动对我们在模型中包含的每个因子类型变量应用了一种称为虚拟编码的技术。

虚拟编码通过为特征的每个类别创建一个二元变量，通常称为虚拟 变量，将一个名义特征处理为数值型特征。如果观察值属于指定的类别，则虚拟变量为1，否则为0。例如，sex特征有两个类别：male和female。这将被拆分为两个二元变量，R 分别命名为sexmale和sexfemale。对于sex = male的观察值，sexmale = 1且sexfemale = 0；反之，如果sex = female，则sexmale = 0且sexfemale = 1。对于具有三个或更多类别的变量也适用相同的编码方式。例如，R 将四类别特征region拆分为四个虚拟变量：regionnorthwest、regionsoutheast、regionsouthwest和regionnortheast。

在向回归模型中添加虚拟变量时，始终会有一个类别被排除作为参考类别。然后，系数会相对于该参考类别进行解释。在我们的模型中，R 自动排除了sexfemale、smokerno和regionnortheast变量，将女性非吸烟者所在的东北地区作为参考组。因此，相对于女性，每年男性的医疗费用少$131.40，而吸烟者每年比非吸烟者多花费$23,847.50。模型中每个三大地区的系数都是负数，这意味着参考组——东北地区的平均费用最高。

提示

默认情况下，R 使用因子变量的第一个级别作为参考。如果你希望使用其他级别，可以使用relevel()函数手动指定参考组。欲了解更多信息，请在 R 中使用?relevel命令。

线性回归模型的结果是合乎逻辑的：老年、吸烟和肥胖倾向于与额外的健康问题相关，而额外的家庭成员抚养负担可能导致医生就诊次数和预防性护理（如疫苗接种和年度体检）的增加。然而，我们目前并不清楚模型拟合数据的效果如何。我们将在下一节回答这个问题。

步骤 4 – 评估模型性能

我们通过输入ins_model获得的参数估计值告诉我们自变量与因变量之间的关系，但它们并没有告诉我们模型与数据的拟合度。为了评估模型性能，我们可以对存储的模型使用summary()命令：

> summary(ins_model)

这会产生以下输出。请注意，输出已经标注了说明性标签：

summary()输出一开始可能让人感到困惑，但基本概念很容易掌握。如前输出中的编号标签所示，输出提供了评估模型性能或拟合度的三种关键方法：

1. 残差部分提供了我们预测中错误的汇总统计信息，其中一些显然相当大。由于残差等于真实值减去预测值，因此最大误差为 29981.7，表明模型在至少一个观测值中低估了近$30,000 的支出。另一方面，50$2,850.90 和真实值之下$1,383.90 之间。

2. 对于每个估计的回归系数，p 值，用Pr(>|t|)表示，提供了一个估算值，表示在给定估计值的情况下，真实系数为零的概率。小的 p 值表明真实系数非常不可能为零，这意味着该特征极不可能与因变量没有关系。请注意，某些 p 值有星号（***），这些星号对应脚注，表示估计值所达到的显著性水平。这个水平是一个阈值，在建立模型之前选择，用来指示“真实”的发现，而不是仅仅由于偶然因素；小于显著性水平的 p 值被视为统计显著。如果模型中有很少这样的项，可能需要引起关注，因为这表明所使用的特征对结果的预测能力较弱。这里，我们的模型有多个高度显著的变量，而且它们似乎与结果在逻辑上相关。

3. 多重 R 方值（也叫决定系数）提供了衡量我们模型整体解释因变量值的能力。它类似于相关系数，数值越接近 1.0，模型对数据的解释越完美。由于 R 方值为 0.7494，我们知道模型解释了因变量近 75% 的变化。由于具有更多特征的模型总是能解释更多变化，调整 R 方值通过对具有大量独立变量的模型进行惩罚来修正 R 方值。它对于比较具有不同数量解释变量的模型表现非常有用。

根据前述的三个性能指标，我们的模型表现得相当不错。回归模型在实际数据中往往会有较低的 R 方值，0.75 的值实际上已经相当不错。一些误差的大小有些令人担忧，但考虑到医疗费用数据的特性，这并不令人惊讶。然而，正如下一节所示，我们可能通过稍微不同的方式指定模型来提高模型的性能。

第 5 步 – 提高模型性能

如前所述，回归建模与其他机器学习方法的一个关键区别在于，回归通常将特征选择和模型指定留给用户。因此，如果我们了解某个特征与结果之间的关系，我们可以利用这些信息来指导模型的指定，从而可能提高模型的性能。

模型指定 – 添加非线性关系

在线性回归中，自变量与因变量之间的关系假定为线性，但这不一定为真。例如，年龄对医疗支出的影响可能在所有年龄值中并非恒定；对于最年长的群体，治疗可能变得不成比例地昂贵。

如果你记得的话，典型的回归方程通常类似于此：

为了考虑非线性关系，我们可以向回归模型中添加更高阶的项，将模型视为多项式。实际上，我们将模拟这样的关系：

这两个模型之间的区别在于将估计一个额外的贝塔系数，旨在捕捉 x 的平方项的影响。这使得年龄的影响可以作为年龄平方的函数来衡量。

要将非线性年龄加入模型，我们只需创建一个新的变量：

> insurance$age2 <- insurance$age²

然后，当我们生成我们改进的模型时，我们将使用expenses ~ age + age2形式将age和age2都添加到lm()公式中。这将使模型能够分离年龄对医疗支出的线性和非线性影响。

转换 - 将数值变量转换为二元指示器

假设我们有一种直觉，即某个特征的影响不是累积的，而是只有在达到特定阈值后才会产生影响。例如，BMI 对于正常体重范围内的个体可能没有医疗支出的影响，但对于肥胖者（即 BMI 达到或超过 30）可能与更高的费用密切相关。

我们可以通过创建一个二元肥胖指示变量来建模这种关系，如果 BMI 至少为 30，则为 1，否则为 0。然后，这个二元特征的估计β值将指示 BMI 达到或超过 30 的个体对医疗支出的平均净影响，相对于 BMI 低于 30 的个体。

要创建该特征，我们可以使用ifelse()函数，该函数对向量中的每个元素测试指定的条件，并根据条件是真还是假返回一个值。对于 BMI 大于或等于 30，我们将返回1，否则返回0：

> insurance$bmi30 <- ifelse(insurance$bmi >= 30, 1, 0)

我们接下来可以在我们改进的模型中包含bmi30变量，可以替换原来的bmi变量或者同时加入，这取决于我们是否认为肥胖的影响会在单独的线性 BMI 效应之外发生。如果没有充分的理由做出改变，我们将在最终模型中都包含它们。

提示

如果您在决定是否包含一个变量时遇到困难，一个常见的做法是将其包含在内并检查 P 值。如果变量在统计上不显著，您有理由在将来排除它。

模型规范化 - 添加交互效应

到目前为止，我们只考虑了每个特征对结果的个体贡献。如果某些特征共同影响因变量呢？例如，吸烟和肥胖可能分别具有有害效应，但合并效应可能比单独每个因素的总和更糟糕是合理的假设。

当两个特征具有联合效应时，这被称为交互作用。如果我们怀疑两个变量之间存在交互作用，我们可以通过将它们的交互项添加到模型中来测试这一假设。交互效应使用 R 公式语法指定。为了让肥胖指示器(bmi30)和吸烟指示器(smoker)交互，我们会编写如下形式的公式：expenses ~ bmi30*smoker。

*运算符是一个快捷方式，指示 R 模型为expenses ~ bmi30 + smokeryes + bmi30:smokeryes。在扩展形式中，冒号:运算符表示bmi30:smokeryes是两个变量之间的交互作用。请注意，扩展形式还自动包括了bmi30和smoker变量以及交互作用。

提示

交互项在模型中绝不能单独包含，而不添加每个交互变量。如果你总是使用*运算符来创建交互作用，那么这个问题就不会出现，因为 R 会自动添加所需的组件。

综合起来——一个改进的回归模型

基于对医学费用可能与患者特征相关的一些主题知识，我们开发了我们认为更准确的回归公式。总结改进的内容，我们：

• 为年龄添加了一个非线性项

• 创建了一个用于肥胖的指标

• 指定了肥胖与吸烟之间的交互作用

我们将像之前一样使用lm()函数来训练模型，但这次我们将添加新构建的变量和交互项：

> ins_model2 <- lm(expenses ~ age + age2 + children + bmi + sex +
 bmi30*smoker + region, data = insurance)

接下来，我们总结结果：

> summary(ins_model2)

输出如下所示：

模型拟合统计量有助于确定我们的修改是否提高了回归模型的性能。与我们的第一个模型相比，R 平方值从 0.75 提高到了大约 0.87。同样，考虑到模型复杂度增加的调整后 R 平方值也从 0.75 提高到了 0.87。我们的模型现在解释了 87%的医疗费用变异。此外，我们关于模型功能形式的理论似乎得到了验证。高阶age2项是统计显著的，肥胖指标bmi30也是如此。肥胖与吸烟的交互作用表明有巨大的影响；除了吸烟本身导致的超过$13,404 的费用外，肥胖吸烟者每年还要额外花费$19,810。这可能表明吸烟加剧了与肥胖相关的疾病。

注意

严格来说，回归建模对数据做出了一些强假设。这些假设对于数值预测来说并不那么重要，因为模型的价值并不依赖于它是否真正捕捉到了潜在的过程——我们关心的是其预测的准确性。然而，如果你希望从回归模型系数中得出明确的推论，就有必要进行诊断测试，确保回归假设没有被违反。关于这一主题的优秀介绍，参见 Allison PD*《多元回归：基础》*，Pine Forge Press，1998 年。

理解回归树和模型树

如果你记得第五章，分而治之 – 使用决策树和规则进行分类，决策树构建的模型就像一个流程图，其中决策节点、叶节点和分支定义了一系列决策，用于分类示例。这样的树也可以通过对树生长算法进行小的调整，用于数值预测。在本节中，我们只会讨论用于数值预测的树与用于分类的树的不同之处。

用于数值预测的树分为两类。第一类称为回归树，它们在 1980 年代作为开创性的分类与回归树（CART）算法的一部分被引入。尽管名称中有回归二字，回归树并没有使用本章前面描述的线性回归方法，而是基于达到叶节点的示例的平均值来进行预测。

注意

CART 算法在 Breiman L, Friedman JH, Stone CJ, Olshen RA 的《分类与回归树》一书中有详细描述。Belmont, CA: Chapman and Hall; 1984 年。

第二类用于数值预测的树被称为模型树。它们比回归树晚几年引入，虽然不太为人所知，但可能更强大。模型树的生长方式与回归树非常相似，但在每个叶节点上，会基于到达该节点的示例构建一个多元线性回归模型。根据叶节点的数量，模型树可能会构建几十个甚至几百个这样的模型。这可能使得模型树比等效的回归树更难理解，但好处是，它们可能会产生更精确的模型。

注意

最早的模型树算法，M5，在 Quinlan JR 的《使用连续类别进行学习》中有描述。1992 年，澳大利亚第五届联合人工智能会议论文集：343-348。

向树中添加回归

可以执行数值预测的树提供了一种引人注目但常被忽视的回归建模替代方案。回归树和模型树相对于更常见的回归方法的优缺点在下表中列出：

优点	缺点

|

• 结合了决策树的优点和建模数值数据的能力

• 无需用户提前指定模型

• 使用自动特征选择，允许该方法应用于非常大量的特征

• 可能比线性回归更适合某些类型的数据

• 不需要统计学知识来解释模型

|

• 不像线性回归那样知名

• 需要大量的训练数据

• 难以确定单个特征对结果的整体净影响

• 大型树可能比回归模型更难以解释

|

虽然传统的回归方法通常是数值预测任务的首选，但在某些情况下，数值决策树提供了明显的优势。例如，决策树可能更适合具有许多特征或特征与结果之间存在许多复杂的非线性关系的任务。这些情况对回归造成了挑战。回归建模还对数值数据的分布方式做出了假设，而这些假设在现实世界的数据中经常被违反。而这种情况在树中并非如此。

用于数值预测的树的构建方式与分类中的方式基本相同。从根节点开始，根据能够在拆分后最大增加结果均匀性的特征使用分治策略对数据进行分区。在分类树中，您会记得均匀性是通过熵来衡量的，但对于数值数据，熵是未定义的。相反，对于数值决策树，均匀性是通过方差、标准偏差或与均值的绝对偏差等统计量来衡量的。

一种常见的拆分标准称为标准偏差减少（SDR）。它由以下公式定义：

在这个公式中，sd(T)函数指的是集合T中值的标准偏差，而T[1]、T[2]、...、T[n]是由特征拆分产生的值集合。|T|项指的是集合T中的观测数量。本质上，该公式通过比较拆分前的标准偏差与加权后的标准偏差来衡量标准偏差的减少。

例如，考虑以下情况，树正在决定是否在二进制特征 A 或 B 上执行拆分：

使用提议的拆分结果得到的组，我们可以计算如下的 A 和 B 的 SDR。这里使用的length()函数返回向量中的元素数量。请注意，整体组 T 被命名为tee，以避免覆盖 R 的内置T()和t()函数：

> tee <- c(1, 1, 1, 2, 2, 3, 4, 5, 5, 6, 6, 7, 7, 7, 7)
> at1 <- c(1, 1, 1, 2, 2, 3, 4, 5, 5)
> at2 <- c(6, 6, 7, 7, 7, 7)
> bt1 <- c(1, 1, 1, 2, 2, 3, 4)
> bt2 <- c(5, 5, 6, 6, 7, 7, 7, 7)
> sdr_a <- sd(tee) - (length(at1) / length(tee) * sd(at1) +
 length(at2) / length(tee) * sd(at2))
> sdr_b <- sd(tee) - (length(bt1) / length(tee) * sd(bt1) +
 length(bt2) / length(tee) * sd(bt2))

让我们比较 A 的 SDR 和 B 的 SDR：

> sdr_a
[1] 1.202815
> sdr_b
[1] 1.392751

特征 A 的拆分的 SDR 约为 1.2，而特征 B 的拆分的 SDR 约为 1.4。由于 B 的标准偏差减少更多，决策树将首先使用 B。这导致比 A 更均匀的集合。

假设树在此处停止生长，并仅使用此一次拆分。回归树的工作完成了。它可以根据特征 B 上的示例值将示例预测为属于组T[1]或T[2]。如果示例位于T[1]中，则模型将预测mean(bt1) = 2，否则将预测mean(bt2) = 6.25。

与此不同，模型树会更进一步。利用落在组* T[1]中的七个训练示例和落在 T[2]中的八个训练示例，模型树可以建立一个线性回归模型，将结果与特征 A 进行比较。请注意，特征 B 在建立回归模型时没有帮助，因为所有叶节点中的示例 B 值相同——它们是根据 B 值被分配到 T[1]或 T[2]*中的。然后，模型树可以使用这两个线性模型中的任何一个为新示例做出预测。

为了进一步说明这两种方法之间的差异，让我们通过一个现实世界的例子来说明。

示例 – 使用回归树和模型树估计葡萄酒的质量

酿酒业是一个具有挑战性且竞争激烈的行业，具有巨大的盈利潜力。然而，葡萄酒厂的盈利能力受到多种因素的影响。作为一种农产品，天气和生长环境等变量会影响葡萄酒的质量。瓶装和生产过程也会对味道产生影响，可能是更好也可能是更差。即便是产品的营销方式，从瓶子设计到定价，也会影响消费者对味道的认知。

因此，酿酒行业在数据收集和机器学习方法上的投资非常巨大，这些方法有助于酿酒决策科学。例如，机器学习已经被用于发现不同地区葡萄酒化学成分的关键差异，或者识别导致葡萄酒味道更甜的化学因素。

最近，机器学习已被用来协助葡萄酒质量的评分——这是一项著名的困难任务。著名酒评家的评论往往决定了产品最终是上架还是下架，尽管即使是专家在盲测中评分时也常常表现出不一致性。

在这个案例研究中，我们将使用回归树和模型树创建一个系统，能够模仿专家对葡萄酒的评分。由于树结构生成的模型易于理解，这可以帮助酿酒师识别出对高评分葡萄酒至关重要的因素。或许更重要的是，该系统不受品酒时人类因素的影响，比如评分者的情绪或味觉疲劳。因此，计算机辅助的葡萄酒测试可能会带来更好的产品，并且评分更加客观、一致和公平。

步骤 1 – 收集数据

为了开发葡萄酒评级模型，我们将使用 P. Cortez, A. Cerdeira, F. Almeida, T. Matos 和 J. Reis 捐赠给 UCI 机器学习数据库的数据（archive.ics.uci.edu/ml）。这些数据包括来自葡萄牙的红葡萄酒和白葡萄酒样本—葡萄牙是世界领先的葡萄酒生产国之一。由于影响葡萄酒高评分的因素可能在红葡萄酒和白葡萄酒之间有所不同，因此在本次分析中，我们只会考察更受欢迎的白葡萄酒。

提示

要跟随本示例，请从 Packt Publishing 网站下载whitewines.csv文件，并将其保存到您的 R 工作目录中。如果您希望自行探索数据，redwines.csv文件也可以下载。

白葡萄酒数据包括 4,898 个葡萄酒样本的 11 个化学属性信息。每个样本通过实验室分析测量了酸度、糖分、氯化物、硫磺、酒精、pH 值和密度等特性。然后，这些样本通过不低于三名评委的盲品进行质量评分，评分范围从 0（非常差）到 10（优秀）。如果评委们对评分有异议，则使用中位数值。

Cortez 的研究评估了三种机器学习方法在建模葡萄酒数据方面的能力：多元回归、人工神经网络和支持向量机。我们在本章早些时候已经介绍了多元回归，而我们将在第七章中学习神经网络和支持向量机，黑箱方法 – 神经网络与支持向量机。研究发现，支持向量机比线性回归模型提供了显著更好的结果。然而，与回归模型不同，支持向量机模型较难解释。通过使用回归树和模型树，我们可能能够在改进回归结果的同时，仍然保持一个易于理解的模型。

注意

要了解更多关于这里描述的葡萄酒研究的信息，请参阅 Cortez P, Cerdeira A, Almeida F, Matos T, Reis J. 通过物理化学性质的数据挖掘建模葡萄酒偏好。决策支持系统*。2009; 47:547-553。

第 2 步 – 探索和准备数据

如常，我们将使用read.csv()函数将数据加载到 R 中。由于所有特征都是数值型的，因此我们可以安全地忽略stringsAsFactors参数：

> wine <- read.csv("whitewines.csv")

葡萄酒数据包括 11 个特征和质量结果，如下所示：

> str(wine)
'data.frame':  4898 obs. of  12 variables:
 $ fixed.acidity       : num  6.7 5.7 5.9 5.3 6.4 7 7.9 ...
 $ volatile.acidity    : num  0.62 0.22 0.19 0.47 0.29 0.12 ...
 $ citric.acid         : num  0.24 0.2 0.26 0.1 0.21 0.41 ...
 $ residual.sugar      : num  1.1 16 7.4 1.3 9.65 0.9 ...
 $ chlorides           : num  0.039 0.044 0.034 0.036 0.041 ...
 $ free.sulfur.dioxide : num  6 41 33 11 36 22 33 17 34 40 ...
 $ total.sulfur.dioxide: num  62 113 123 74 119 95 152 ...
 $ density             : num  0.993 0.999 0.995 0.991 0.993 ...
 $ pH                  : num  3.41 3.22 3.49 3.48 2.99 3.25 ...
 $ sulphates           : num  0.32 0.46 0.42 0.54 0.34 0.43 ...
 $ alcohol             : num  10.4 8.9 10.1 11.2 10.9 ...
 $ quality             : int  5 6 6 4 6 6 6 6 6 7 ...

与其他类型的机器学习模型相比，树模型的一个优势是能够处理多种类型的数据而无需预处理。这意味着我们不需要对特征进行归一化或标准化。

然而，需要稍微努力检查结果变量的分布，以便更好地评估模型的表现。例如，假设酒的质量变化很小，或者酒的质量呈双峰分布：要么非常好，要么非常差。为了检查这些极端情况，我们可以使用直方图检查质量分布：

> hist(wine$quality)

这将生成以下图形：

酒的质量值似乎遵循一个相当正常的钟形分布，集中在六的值附近。这是直观合理的，因为大多数酒的质量是中等的；很少有酒特别差或者特别好。尽管这里没有展示结果，但检查 summary(wine) 输出中的异常值或其他潜在的数据问题也是很有用的。尽管树模型对杂乱数据相当稳健，但检查严重问题总是明智的。现在，我们假设数据是可靠的。

我们的最后一步是将数据划分为训练集和测试集。由于 wine 数据集已经是随机排序的，我们可以按照以下方式将其分成两组连续的行：

> wine_train <- wine[1:3750, ]
> wine_test <- wine[3751:4898, ]

为了与 Cortez 使用的条件相符，我们分别使用了 75% 和 25% 的数据集进行训练和测试。我们将评估基于树的模型在测试数据上的表现，以查看我们是否能获得与先前研究相当的结果。

第三步 – 在数据上训练模型

我们将首先训练一棵回归树模型。尽管几乎任何决策树的实现都可以用于回归树建模，但 rpart（递归分区）包提供了最忠实的回归树实现，正如 CART 团队所描述的那样。作为 CART 的经典 R 实现，rpart 包也有很好的文档支持，并提供了可视化和评估 rpart 模型的函数。

使用 install.packages("rpart") 命令安装 rpart 包。然后可以使用 library(rpart) 命令将其加载到 R 会话中。以下语法将使用默认设置训练一棵树，这些设置通常工作得相当好。如果需要更精细的设置，请使用 ?rpart.control 命令查看控制参数的文档。

使用 R 公式接口，我们可以指定 quality 作为结果变量，并使用点符号使得 wine_train 数据框中的所有其他列作为预测变量。生成的回归树模型对象命名为 m.rpart，以便与我们稍后训练的模型树区分开来：

> m.rpart <- rpart(quality ~ ., data = wine_train)

要获取树的基本信息，只需输入模型对象的名称：

> m.rpart
n= 3750

node), split, n, deviance, yval
 * denotes terminal node

 1) root 3750 2945.53200 5.870933 
 2) alcohol< 10.85 2372 1418.86100 5.604975 
 4) volatile.acidity>=0.2275 1611  821.30730 5.432030 
 8) volatile.acidity>=0.3025 688  278.97670 5.255814 *
 9) volatile.acidity< 0.3025 923  505.04230 5.563380 *
 5) volatile.acidity< 0.2275 761  447.36400 5.971091 *
 3) alcohol>=10.85 1378 1070.08200 6.328737 
 6) free.sulfur.dioxide< 10.5 84   95.55952 5.369048 *
 7) free.sulfur.dioxide>=10.5 1294  892.13600 6.391036 
 14) alcohol< 11.76667 629  430.11130 6.173291 
 28) volatile.acidity>=0.465 11   10.72727 4.545455 *
 29) volatile.acidity< 0.465 618  389.71680 6.202265 *
 15) alcohol>=11.76667 665  403.99400 6.596992 *

对于树中的每个节点，都会列出到达该决策点的示例数量。例如，所有 3,750 个示例从根节点开始，其中 2,372 个示例满足alcohol < 10.85，1,378 个示例满足alcohol >= 10.85。由于酒精在树中首先被使用，它是预测葡萄酒质量的最重要因素。

由*指示的节点是终端节点或叶节点，这意味着它们会产生预测结果（在此列为yval）。例如，节点 5 的yval为 5.971091。当使用该树进行预测时，任何酒精含量alcohol < 10.85且挥发性酸度volatile.acidity < 0.2275的葡萄酒样本，预测的质量值为 5.97。

可以使用summary(m.rpart)命令获取树的拟合结果的更详细摘要，包括每个节点的均方误差和特征重要性的总体度量。

可视化决策树

尽管可以仅通过前面的输出理解树的结构，但通常使用可视化会更容易理解。Stephen Milborrow 的rpart.plot包提供了一个易于使用的函数，可以生成高质量的决策树图。

注意

有关rpart.plot的更多信息，包括该函数可以生成的不同类型决策树图的其他示例，请访问作者的网站：www.milbo.org/rpart-plot/。

在使用install.packages("rpart.plot")命令安装该包后，rpart.plot()函数可以从任何rpart模型对象生成树形图。以下命令绘制了我们之前构建的回归树：

> library(rpart.plot)
> rpart.plot(m.rpart, digits = 3)

生成的树形图如下所示：

除了控制图中数字位数的digits参数外，许多可视化的其他方面也可以进行调整。以下命令展示了几种有用的选项：fallen.leaves参数强制叶节点对齐到图的底部，而type和extra参数则影响决策和节点的标签方式：

> rpart.plot(m.rpart, digits = 4, fallen.leaves = TRUE,
 type = 3, extra = 101)

这些更改的结果是一个外观截然不同的树形图：

像这样的可视化图可能有助于回归树结果的传播，因为即使没有数学背景的人也能轻松理解。在这两种情况下，叶节点中显示的数字是到达该节点的示例的预测值。因此，将图表展示给葡萄酒生产者，可能有助于确定预测高评分葡萄酒的关键因素。

第 4 步 – 评估模型性能

要使用回归树模型对测试数据进行预测，我们使用predict()函数。默认情况下，它返回结果变量的估计数值，我们将其保存在名为p.rpart的向量中：

> p.rpart <- predict(m.rpart, wine_test)

快速查看我们预测的摘要统计数据表明了一个潜在的问题；预测值的范围明显比真实值要窄：

> summary(p.rpart)
 Min. 1st Qu.  Median    Mean 3rd Qu.    Max.
 4.545   5.563   5.971   5.893   6.202   6.597
> summary(wine_test$quality)
 Min. 1st Qu.  Median    Mean 3rd Qu.    Max.
 3.000   5.000   6.000   5.901   6.000   9.000

这一发现表明模型未能正确识别极端情况，尤其是最好的和最差的葡萄酒。另一方面，在第一四分位数和第三四分位数之间，我们可能做得很好。

预测值与实际质量值之间的相关性提供了一种简单的方法来衡量模型的表现。请记住，cor()函数可用于衡量两个等长向量之间的关系。我们将使用它来比较预测值与真实值的对应情况：

> cor(p.rpart, wine_test$quality)
[1] 0.5369525

0.54 的相关性是完全可以接受的。然而，相关性仅仅衡量预测值与真实值之间的关联程度；它并不能衡量预测值与真实值之间的差距。

使用平均绝对误差衡量性能

另一种衡量模型表现的方法是考虑其预测与真实值之间的平均偏差。这一度量被称为平均绝对误差（MAE）。MAE 的计算公式如下，其中n表示预测次数，*e[i]*表示第 i 个预测的误差：

正如名字所示，这个方程计算误差的绝对值的平均值。由于误差仅仅是预测值和实际值之间的差异，我们可以创建一个简单的MAE()函数，如下所示：

> MAE <- function(actual, predicted) {
 mean(abs(actual - predicted))
}

我们的预测的 MAE 为：

> MAE(p.rpart, wine_test$quality)
[1] 0.5872652

这意味着，平均而言，我们模型的预测值与真实质量评分之间的差异大约为 0.59。在从零到 10 的质量评分尺度上，这似乎表明我们的模型表现相当不错。

另一方面，请记住，大多数葡萄酒既不特别好也不特别差；典型的质量评分大约在五到六之间。因此，仅仅预测平均值的分类器根据这个标准仍然可能表现得相当好。

训练数据中的平均质量评分如下：

> mean(wine_train$quality)
[1] 5.870933

如果我们对每个葡萄酒样本都预测值为 5.87，我们的平均绝对误差将仅约为 0.67：

> MAE(5.87, wine_test$quality)
[1] 0.6722474

我们的回归树（MAE = 0.59）的平均预测结果比填充平均值（MAE = 0.67）更接近真实质量评分，但差距不大。相比之下，Cortez 报告了神经网络模型的 MAE 为 0.58，支持向量机的 MAE 为 0.45。这表明仍有改进的空间。

步骤 5 – 改进模型性能

为了提高学习器的性能，我们可以尝试构建一个模型树。回想一下，模型树通过将叶节点替换为回归模型，从而改进了回归树。这通常会比回归树产生更准确的结果，因为回归树在叶节点使用的只是一个单一的预测值。

当前模型树的最新技术是M5'算法（M5-prime），由 Y. Wang 和 I.H. Witten 提出，它是 J.R. Quinlan 于 1992 年提出的原始 M5 模型树算法的变种。

注意

欲了解更多有关 M5'算法的信息，请参阅 Wang Y, Witten IH. 用于预测连续类别的模型树诱导。欧洲机器学习会议海报论文集，1997 年。

M5 算法可以通过 R 中的RWeka包和M5P()函数来实现。该函数的语法如下表所示。如果你还没有安装RWeka包，请确保先安装。由于它依赖于 Java，安装说明已包含在第一章中，引入机器学习。

我们将使用与回归树相同的语法来拟合模型树：

> library(RWeka)
> m.m5p <- M5P(quality ~ ., data = wine_train)

可以通过输入树的名称来检查树本身。在这种情况下，树非常大，只有前几行输出会显示出来：

> m.m5p
M5 pruned model tree:
(using smoothed linear models)

alcohol <= 10.85 :
|   volatile.acidity <= 0.238 :
|   |   fixed.acidity <= 6.85 : LM1 (406/66.024%)
|   |   fixed.acidity >  6.85 :
|   |   |   free.sulfur.dioxide <= 24.5 : LM2 (113/87.697%)

你会注意到，这些分裂与我们之前构建的回归树非常相似。酒精是最重要的变量，其次是挥发性酸度和游离二氧化硫。然而，一个关键的区别是，节点的终止不是数值预测，而是一个线性模型（在此显示为LM1和LM2）。

线性模型本身会在输出中显示。例如，LM1的模型将在接下来的输出中显示。这些值的解释方式与我们在本章前面构建的多重回归模型完全相同。每个数字表示相关特征对预测酒质的净效应。固定酸度的系数为0.266，意味着酸度增加 1 单位时，预计酒质将增加0.266：

LM num: 1
quality =
 0.266 * fixed.acidity
 - 2.3082 * volatile.acidity
 - 0.012 * citric.acid
 + 0.0421 * residual.sugar
 + 0.1126 * chlorides
 + 0 * free.sulfur.dioxide
 - 0.0015 * total.sulfur.dioxide
 - 109.8813 * density
 + 0.035 * pH
 + 1.4122 * sulphates
 - 0.0046 * alcohol
 + 113.1021

需要注意的是，由LM1估算的效果仅适用于达到该节点的酒样；在这个模型树中总共构建了 36 个线性模型，每个模型都对固定酸度和其他 10 个特征的影响进行了不同的估算。

若要查看模型拟合训练数据的统计信息，可以对 M5P 模型应用summary()函数。然而，请注意，由于这些统计信息是基于训练数据的，因此它们仅应作为粗略的诊断依据：

> summary(m.m5p)

=== Summary ===

Correlation coefficient                  0.6666
Mean absolute error                      0.5151
Root mean squared error                  0.6614
Relative absolute error                 76.4921 %
Root relative squared error             74.6259 %
Total Number of Instances             3750 

相反，我们将查看模型在未见过的测试数据上的表现。predict()函数为我们提供了一组预测值：

> p.m5p <- predict(m.m5p, wine_test)

模型树似乎预测了比回归树更广泛的值范围：

> summary(p.m5p)
 Min. 1st Qu.  Median    Mean 3rd Qu.    Max.
 4.389   5.430   5.863   5.874   6.305   7.437

相关性似乎也显著更高：

> cor(p.m5p, wine_test$quality)
[1] 0.6272973

此外，该模型略微减少了平均绝对误差：

> MAE(wine_test$quality, p.m5p)
[1] 0.5463023

尽管我们在回归树的基础上并没有取得太大进展，但我们超越了 Cortez 发布的神经网络模型的表现，并且我们逐渐接近支持向量机模型发布的平均绝对误差值 0.45，所有这些都是通过使用一种更简单的学习方法实现的。

提示

毫不意外，我们确认预测葡萄酒质量是一个困难的问题；毕竟，葡萄酒品鉴本质上是主观的。如果你想要额外练习，可以在阅读第十一章，提高模型性能后再次尝试解决这个问题，里面涵盖了一些可能带来更好结果的技术。

总结

在本章中，我们研究了两种建模数值数据的方法。第一种方法是线性回归，它涉及将直线拟合到数据上。第二种方法使用决策树进行数值预测。后者有两种形式：回归树，它通过叶节点处示例的平均值来进行数值预测；以及模型树，它在每个叶节点构建一个回归模型，这是一种混合方法，在某些方面是两者的最佳结合。

我们使用线性回归模型来计算不同人口群体的预期医疗费用。由于特征与目标变量之间的关系通过估算的回归模型得到了很好的描述，我们能够识别出某些人群，如吸烟者和肥胖者，可能需要被收取更高的保险费率，以覆盖高于平均水平的医疗费用。

回归树和模型树被用来根据可衡量的特征对葡萄酒的主观质量进行建模。在此过程中，我们了解到回归树提供了一种简单的方法来解释特征与数值结果之间的关系，但更复杂的模型树可能更为准确。在这个过程中，我们学习了几种评估数值模型表现的方法。

与本章涉及的机器学习方法（这些方法能够清晰地理解输入与输出之间的关系）形成鲜明对比的是，下章将介绍一些导致近乎难以理解的模型的方法。其优点是，这些方法是极为强大的技术——是最强大的股票分类器之一——可以应用于分类和数值预测问题。

第七章：黑箱方法——神经网络与支持向量机

已故科幻作家阿瑟·C·克拉克曾写道：“任何足够先进的技术都无法与魔法区分开。”本章涵盖了一对初看起来像魔法的机器学习方法。尽管它们极其强大，但其内部原理可能难以理解。

在工程学中，这些被称为黑箱过程，因为将输入转化为输出的机制被一个虚拟的盒子所遮掩。例如，闭源软件的黑箱故意隐藏专有算法，政治立法的黑箱根植于官僚流程，而香肠制造的黑箱则涉及一些故意的（但美味的）无知。在机器学习的情况下，黑箱则源于其运作所依赖的复杂数学。

虽然它们可能不易理解，但盲目应用黑箱模型是危险的。因此，在本章中，我们将一窥黑箱内部，并调查拟合此类模型所涉及的统计香肠制作过程。你将发现：

• 神经网络模仿动物大脑的结构来模拟任意函数

• 支持向量机使用多维表面来定义特征与结果之间的关系

• 尽管它们的复杂性，依然可以轻松应用于现实世界问题

如果幸运的话，你会意识到，解决黑箱机器学习方法并不需要统计学的黑带——完全没有必要感到畏惧！

了解神经网络

人工神经网络（ANN）使用一种基于我们对生物大脑如何响应来自感官输入的刺激的理解所衍生的模型，来模拟输入信号与输出信号之间的关系。就像大脑通过一组叫做神经元的相互连接的细胞构建一个巨大的并行处理器，ANN 则利用一组人工神经元或节点来解决学习问题。

人脑由约 850 亿个神经元组成，形成了一个能够表示大量知识的网络。正如你所预料的，这一数量远远超过其他生物的脑量。例如，一只猫大约有 10 亿个神经元，一只老鼠大约有 7500 万个神经元，而一只蟑螂只有大约 100 万个神经元。相比之下，许多 ANN 包含的神经元要少得多，通常只有几百个，因此我们目前离制造人工大脑还远——即便是一个拥有 10 万个神经元的果蝇大脑，也远远超出了现有的最先进 ANN 的能力。

虽然完全模拟一只蟑螂的大脑可能不可行，但神经网络仍然可以提供一个足够的启发式模型来模拟其行为。假设我们开发了一个算法，能够模仿蟑螂在被发现时逃跑的反应。如果机器人蟑螂的行为让人信服，那么它的大脑是否与活体生物一样复杂重要吗？这个问题正是有争议的图灵测试的基础。图灵测试由开创性计算机科学家艾伦·图灵于 1950 年提出，旨在通过判断一个人类是否无法将机器的行为与活体生物区分开，来评定机器是否具备智能。

基本的人工神经网络（ANNs）已经被使用超过 50 年，用来模拟大脑解决问题的方法。最初，这包括学习简单的函数，比如逻辑与（AND）函数或逻辑或（OR）函数。这些早期的练习主要是为了帮助科学家理解生物大脑如何运作。然而，随着近年来计算机性能的不断增强，人工神经网络的复杂性也大幅增加，以至于现在它们常常被应用于更实际的问题，包括：

• 语音和手写识别程序，例如语音邮件转录服务和邮政邮件分拣机所使用的技术

• 智能设备的自动化，如办公楼的环境控制、自动驾驶汽车和自驾无人机

• 精密的天气和气候模式、抗拉强度、流体动力学以及许多其他科学、社会或经济现象的模型

广义来说，人工神经网络是多才多艺的学习者，几乎可以应用于任何学习任务：分类、数值预测，甚至是无监督的模式识别。

提示

无论是否值得，人工神经网络学习器常常在媒体中被大肆宣传。例如，谷歌最近开发的“人工大脑”因其能够识别 YouTube 上的猫视频而受到推崇。这种炒作可能与人工神经网络的独特性关系不大，而更多地与人工神经网络因其与活体思维的相似性而引人入胜有关。

人工神经网络最适用于那些输入数据和输出数据定义明确或至少相对简单，但输入与输出之间的过程极其复杂的问题。作为一种黑箱方法，它们在这类黑箱问题中表现优异。

从生物神经元到人工神经元

因为人工神经网络（ANNs）是故意设计为人类大脑活动的概念模型，因此首先理解生物神经元的功能是很有帮助的。如以下图所示，细胞的树突通过生化过程接收传入的信号。这个过程使得冲动信号可以根据其相对重要性或频率进行加权。当细胞体开始积累传入的信号时，会达到一个阈值，在此阈值下细胞会激发，并通过电化学过程将输出信号沿着轴突传递。在轴突的末端，电信号再次被处理为化学信号，并通过一个被称为突触的小间隙传递给相邻的神经元。

单个人工神经元的模型可以用非常类似于生物模型的方式理解。如以下图所示，一个有向网络图定义了输入信号（由树突接收的 x 变量）和输出信号（y 变量）之间的关系。就像生物神经元一样，每个树突的信号根据其重要性被加权（w 值）——暂时忽略这些权重是如何确定的。输入信号被细胞体汇总，并根据激活函数 f 传递信号：

一个典型的人工神经元有n个输入树突，可以通过以下公式表示。w 权重使得每个n个输入（用*x[i]*表示）能够对输入信号的总和贡献更多或更少。总和会被激活函数 f(x) 使用，结果信号 y(x) 就是输出轴突的信号：

神经网络使用这样定义的神经元作为构建块来构建复杂的数据模型。尽管神经网络有许多变种，但每个变种都可以通过以下特征来定义：

• 激活函数，它将神经元的综合输入信号转换为一个单一的输出信号，然后将其进一步广播到网络中

• 网络拓扑（或架构），它描述了模型中神经元的数量、层数以及它们如何连接的方式

• 训练算法，它指定了如何设置连接权重，以便根据输入信号来抑制或激发神经元

让我们看一下每个类别内的一些变化，看看它们如何被用来构建典型的神经网络模型。

激活函数

激活函数是人工神经元处理传入信息并将其传递到整个网络的机制。就像人工神经元是基于生物神经元模型的，激活函数也是基于自然界的设计模型的。

在生物学的情况下，激活函数可以被想象为一个过程，涉及将所有输入信号加总并确定是否达到触发阈值。如果是，神经元将传递信号；否则，它什么也不做。在人工神经网络（ANN）术语中，这被称为阈值激活函数，因为它只有在达到指定的输入阈值后才会产生输出信号。

以下图显示了一个典型的阈值函数；在这种情况下，当输入信号的总和至少为零时，神经元触发。因为它的形状类似楼梯，有时被称为单位阶跃激活函数。

尽管阈值激活函数因其与生物学的相似性而具有趣味性，但在人工神经网络中很少使用。摆脱了生物化学的限制，ANN 的激活函数可以根据它们展示期望的数学特性和准确建模数据之间关系的能力来选择。

或许最常用的替代函数是sigmoid 激活函数（更具体来说，是logistic sigmoid），如以下图所示。请注意，在公式中，e是自然对数的底数（约为 2.72）。尽管它与阈值激活函数共享类似的阶梯或“S”形状，但输出信号不再是二进制的；输出值可以落在 0 到 1 的范围内。此外，sigmoid 是可微分的，这意味着可以计算整个输入范围内的导数。如你稍后所学，这一特性对于创建高效的 ANN 优化算法至关重要。

尽管 sigmoid 可能是最常用的激活函数，且通常作为默认选项使用，但一些神经网络算法允许选择其他替代函数。以下图展示了这种激活函数的选择：

区分这些激活函数的主要细节是输出信号的范围。通常，这个范围是（0，1）、（-1，+1）或（-∞，+∞）。激活函数的选择会影响神经网络的偏置，使其能够更适合某些类型的数据，从而构建专门的神经网络。例如，线性激活函数会使神经网络非常类似于线性回归模型，而高斯激活函数则会导致一个称为径向基函数（RBF）网络的模型。每个模型都有其更适合某些学习任务的优点。

重要的是要认识到，对于许多激活函数来说，影响输出信号的输入值范围相对较窄。例如，在 Sigmoid 函数的情况下，当输入信号低于*-5或高于+5*时，输出信号始终接近 0 或 1。这种信号压缩会导致在高低端出现饱和信号，就像将吉他放大器音量调得过高导致声音失真一样，因为音波的峰值被削波。由于这会将输入值压缩到更小的输出范围，像 Sigmoid 这样的激活函数有时被称为压缩函数。

解决压缩问题的方法是对所有神经网络输入进行变换，使得特征的值落在接近 0 的一个小范围内。通常，这涉及对特征进行标准化或归一化。通过限制输入值的范围，激活函数将在整个范围内起作用，从而防止像家庭收入这样的高值特征主导像家庭中孩子数量这样的低值特征。一个额外的好处是，模型的训练速度可能更快，因为算法可以更快速地在有效的输入值范围内进行迭代。

提示

尽管理论上神经网络可以通过多次迭代调整其权重来适应非常动态的特征，但在极端情况下，许多算法会在此之前停止迭代。如果你的模型预测结果不合理，请仔细检查是否正确标准化了输入数据。

网络拓扑

神经网络的学习能力根植于其拓扑结构，即互联神经元的模式和结构。尽管网络架构形式多种多样，但它们可以通过三个关键特征来区分：

• 层数

• 网络中的信息是否允许反向传播

• 网络中每一层的节点数量

拓扑结构决定了网络能够学习的任务复杂性。通常，更大且更复杂的网络能够识别更微妙的模式和复杂的决策边界。然而，网络的能力不仅取决于网络的大小，还与单元的排列方式有关。

层数

为了定义拓扑结构，我们需要一种术语来区分基于其在网络中位置的人工神经元。接下来的图示例了一个非常简单网络的拓扑结构。一组被称为输入节点的神经元直接接收来自输入数据的未处理信号。每个输入节点负责处理数据集中一个单独的特征；该特征的值将通过对应节点的激活函数进行转换。输入节点发送的信号被输出节点接收，输出节点利用自己的激活函数生成最终预测（此处表示为p）。

输入节点和输出节点按层的形式分组排列。由于输入节点以接收到的数据的原始形式进行处理，因此网络只有一组连接权重（此处标记为w[1]，w[2]，和w[3]）。因此，它被称为单层网络。单层网络可以用于基本的模式分类，尤其是对于线性可分的模式，但大多数学习任务需要更复杂的网络。

正如你可能预期的那样，创建更复杂网络的一种明显方法是通过增加额外的层。如图所示，多层网络增加了一个或多个隐藏层，它们在信号到达输出节点之前处理来自输入节点的信号。大多数多层网络是全连接的，这意味着一个层中的每个节点都与下一层中的每个节点相连，但这并不是必须的。

信息流动的方向

你可能已经注意到，在之前的例子中，箭头用来表示信号仅朝一个方向流动。那些输入信号从一个连接到另一个连接，持续朝一个方向传输直到到达输出层的网络被称为前馈网络。

尽管存在信息流动的限制，前馈网络仍然提供了令人惊讶的灵活性。例如，可以改变每一层的节点数和层数，可以同时建模多个结果，或可以应用多个隐藏层。具有多个隐藏层的神经网络被称为深度神经网络（DNN），而训练这种网络的实践有时被称为深度学习。

相比之下，递归网络（或反馈网络）允许信号通过循环在两个方向上传播。这一特性，更加接近生物神经网络的工作方式，使得学习极为复杂的模式成为可能。加入短期记忆或延迟，极大增强了递归网络的能力。特别是，这使得网络能够理解一段时间内的事件序列。这可以用于股市预测、语音理解或天气预报等应用。一个简单的递归网络如下所示：

尽管具有潜力，递归网络仍然主要是理论性的，实际中很少使用。另一方面，前馈网络已广泛应用于现实问题中。实际上，多层前馈网络，有时被称为多层感知器（MLP），是事实上的标准人工神经网络拓扑结构。如果有人提到他们正在拟合神经网络，那么他们很可能是在指 MLP。

每层的节点数量

除了层数和信息传递方向的变化外，神经网络的复杂性还可以通过每层的节点数量来变化。输入节点的数量由输入数据中的特征数量预先确定。类似地，输出节点的数量由要建模的结果数量或结果中的类别层次数预先确定。然而，隐藏层节点的数量则由用户在训练模型之前决定。

不幸的是，没有可靠的规则来确定隐藏层中的神经元数量。合适的数量取决于输入节点的数量、训练数据的数量、噪声数据的数量以及学习任务的复杂性等许多因素。

一般来说，具有更多网络连接的复杂网络拓扑结构可以学习更复杂的问题。更多的神经元将导致一个更贴近训练数据的模型，但这也带来了过拟合的风险；它可能在未来数据上表现不佳。大型神经网络还可能在计算上非常昂贵，并且训练速度较慢。

最佳做法是使用最少的节点，以在验证数据集上获得足够的性能。在大多数情况下，即使只有少量的隐藏节点——通常只有几个——神经网络也能提供极大的学习能力。

提示

已经证明，至少具有一个隐藏层且神经元足够的神经网络是通用函数逼近器。这意味着神经网络可以用于在有限区间内，以任意精度逼近任何连续函数。

使用反向传播训练神经网络

网络拓扑结构本身是空白的，尚未学习到任何内容。就像一个新生儿，它必须通过经验来训练。当神经网络处理输入数据时，神经元之间的连接会被加强或削弱，类似于婴儿大脑在经历环境时的发育过程。网络的连接权重会根据观察到的模式随时间调整。

通过调整连接权重训练神经网络非常耗费计算资源。因此，尽管人工神经网络（ANNs）在此之前已被研究了几十年，但直到 1980 年代中后期，当一种高效的 ANN 训练方法被发现时，ANN 才被应用于现实世界的学习任务。这种算法通过反向传播误差的策略训练网络，现在通常被称为反向传播（backpropagation）。

注意

巧合的是，多个研究团队在大致相同的时间独立地发现并发布了反向传播算法。其中，最常被引用的工作之一是：Rumelhart DE, Hinton GE, Williams RJ. 通过反向传播误差学习表示。自然杂志，1986 年；323：533-566。

尽管相对于许多其他机器学习算法仍然非常缓慢，反向传播方法却促使了对人工神经网络（ANNs）的重新关注。因此，使用反向传播算法的多层前馈网络现在在数据挖掘领域中非常常见。这些模型具有以下优点和缺点：

优势	劣势

|

• 可以适用于分类或数值预测问题

• 能够建模比几乎任何其他算法更复杂的模式

• 对数据的潜在关系假设较少

|

• 计算量极大且训练速度慢，特别是当网络拓扑结构复杂时

• 非常容易导致过拟合训练数据

• 结果是一个复杂的黑箱模型，难以理解，甚至不可能理解。

|

在其最一般的形式中，反向传播算法通过两个过程的多个循环进行迭代。每个循环被称为一个周期（epoch）。因为网络不包含任何先验（已有的）知识，所以起始权重通常是随机设置的。然后，算法在两个过程中迭代，直到达到停止标准。反向传播算法中的每个周期包括：

• 一个前向阶段，在该阶段，神经元从输入层到输出层依次激活，并在过程中应用每个神经元的权重和激活函数。当达到最后一层时，会生成一个输出信号。

• 一个后向阶段，在该阶段，网络在前向阶段产生的输出信号与训练数据中的真实目标值进行比较。网络输出信号与真实值之间的差异导致一个误差，该误差向后传播至网络，以修改神经元之间的连接权重，减少未来的误差。

随着时间的推移，网络会利用向后传播的信息来减少网络的总误差。然而，仍然有一个问题：由于每个神经元的输入和输出之间的关系是复杂的，算法如何确定权重应变化多少？这个问题的答案涉及到一个名为梯度下降的技术。从概念上讲，它的工作原理类似于被困在丛林中的探险者如何找到通往水源的路径。通过检查地形并不断朝着最大下坡的方向前进，探险者最终会到达最低的山谷，这里很可能就是河床。

在类似的过程中，反向传播算法使用每个神经元激活函数的导数来识别各输入权重方向上的梯度——因此，具有可微分激活函数非常重要。梯度指示了权重变化时，误差会如何陡峭地增加或减少。该算法将尝试通过一个称为学习率的量来改变权重，以实现误差的最大减少。学习率越大，算法越快地沿着梯度下降，这可以减少训练时间，但也可能有超越最低点的风险。

尽管这个过程看起来复杂，但在实际应用中却很容易操作。我们将运用对多层前馈网络的理解来解决一个现实问题。

示例 – 使用 ANNs 建模混凝土强度

在工程领域，准确估算建筑材料的性能至关重要。这些估算对于制定建筑物、桥梁和道路施工材料的安全指南是必要的。

估算混凝土的强度是一个特别具有挑战性的任务。尽管混凝土几乎在每个建筑项目中都有使用，但其性能因多种不同的成分以及复杂的相互作用而存在很大差异。因此，很难准确预测最终产品的强度。一个能够根据输入材料的组成来可靠预测混凝土强度的模型，可能会导致更安全的建筑实践。

第一步 – 收集数据

对于本分析，我们将使用由叶怡诚提供给 UCI 机器学习数据集库的混凝土抗压强度数据（archive.ics.uci.edu/ml）。由于他成功地利用神经网络建模这些数据，我们也将尝试使用简单的神经网络模型在 R 中复制他的工作。

注意

欲了解叶氏方法在此学习任务中的应用，请参阅：叶怡诚. 使用人工神经网络建模高性能混凝土的强度. 水泥与混凝土研究. 1998; 28:1797-1808.

根据网站信息，混凝土数据集包含 1,030 个混凝土样本，这些样本有八个特征，描述了混合物中所用成分的情况。这些特征被认为与最终的抗压强度有关，包含了水泥、矿渣、粉煤灰、水、减水剂、粗骨料和细骨料的数量（单位为千克每立方米），以及养护时间（以天为单位）。

提示

要跟随这个示例，请从 Packt Publishing 网站下载concrete.csv文件，并将其保存到您的 R 工作目录中。

第 2 步 – 探索和准备数据

和往常一样，我们将通过使用read.csv()函数将数据加载到 R 对象中，并确认它符合预期的结构：

> concrete <- read.csv("concrete.csv")
> str(concrete)
'data.frame':   1030 obs. of  9 variables:
 $ cement      : num  141 169 250 266 155 ...
 $ slag        : num  212 42.2 0 114 183.4 ...
 $ ash         : num  0 124.3 95.7 0 0 ...
 $ water       : num  204 158 187 228 193 ...
 $ superplastic: num  0 10.8 5.5 0 9.1 0 0 6.4 0 9 ...
 $ coarseagg   : num  972 1081 957 932 1047 ...
 $ fineagg     : num  748 796 861 670 697 ...
 $ age         : int  28 14 28 28 28 90 7 56 28 28 ...
 $ strength    : num  29.9 23.5 29.2 45.9 18.3 ...

数据框中的九个变量对应着我们预期的八个特征和一个结果，尽管一个问题已经显现出来。神经网络在输入数据缩放到接近零的狭窄范围时表现最佳，而在这里，我们看到的值从零到一千多不等。

通常，解决这个问题的方法是使用归一化或标准化函数对数据进行重新缩放。如果数据遵循钟形曲线（如第二章，管理和理解数据中所描述的正态分布），那么使用 R 的内置scale()函数进行标准化可能是有意义的。另一方面，如果数据呈均匀分布或严重偏离正态分布，则可能更适合将数据归一化到 0-1 范围。在这种情况下，我们将使用后者。

在第三章，懒学习——使用最近邻进行分类中，我们定义了自己的normalize()函数，如下所示：

> normalize <- function(x) {
 return((x - min(x)) / (max(x) - min(x)))
 }

执行此代码后，我们的normalize()函数可以使用lapply()函数应用于混凝土数据框的每一列，如下所示：

> concrete_norm <- as.data.frame(lapply(concrete, normalize))

为了确认归一化操作已生效，我们可以看到，最小强度和最大强度现在分别为 0 和 1：

> summary(concrete_norm$strength)
 Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
 0.0000  0.2664  0.4001  0.4172  0.5457  1.0000

相比之下，原始的最小值和最大值分别为 2.33 和 82.60：

> summary(concrete$strength)
 Min. 1st Qu.  Median    Mean 3rd Qu.    Max. 
 2.33   23.71   34.44   35.82   46.14   82.60

提示

在训练模型之前对数据应用的任何转换，之后都必须逆向应用，以便转换回原始的计量单位。为了方便重新缩放，最好保存原始数据，或者至少保存原始数据的汇总统计信息。

根据 Yeh 在原始出版物中的做法，我们将数据分为训练集和测试集，训练集占 75%的样本，测试集占 25%。我们使用的 CSV 文件已经是随机排序的，因此我们只需要将其分成两部分：

> concrete_train <- concrete_norm[1:773, ]
> concrete_test <- concrete_norm[774:1030, ]

我们将使用训练数据集来构建神经网络，并使用测试数据集评估模型如何将结果推广到未来。由于神经网络容易发生过拟合，因此这个步骤非常重要。

第三步 – 在数据上训练模型

为了建模混凝土中所用原料与成品强度之间的关系，我们将使用一个多层前馈神经网络。Stefan Fritsch 和 Frauke Guenther 的neuralnet包提供了这样网络的标准且易于使用的实现。它还提供了一个函数用于绘制网络拓扑。因此，neuralnet的实现是学习神经网络的一个强有力选择，尽管这并不意味着它不能用于完成实际工作——它是一个非常强大的工具，正如你很快会看到的那样。

提示

还有几个其他常用的包可以用来训练 R 中的 ANN 模型，每个包都有独特的优缺点。由于它作为 R 标准安装的一部分，nnet包可能是引用最多的 ANN 实现。它使用的算法比标准的反向传播算法稍微复杂一些。另一个强大的选择是RSNNS包，它提供了完整的神经网络功能，但缺点是它更难学习。

由于neuralnet不包含在基础 R 中，你需要通过输入install.packages("neuralnet")来安装它，并使用library(neuralnet)命令加载它。所包含的neuralnet()函数可以用于训练用于数值预测的神经网络，语法如下：

我们将从训练最简单的多层前馈网络开始，只有一个隐藏节点：

> concrete_model <- neuralnet(strength ~ cement + slag + ash + water + superplastic + coarseagg + fineagg + age, data = concrete_train)

然后我们可以使用结果模型对象上的plot()函数来可视化网络拓扑：

> plot(concrete_model)

在这个简单模型中，每个八个特征都有一个输入节点，接着是一个隐藏节点和一个输出节点，输出节点预测混凝土强度。每个连接的权重也被描绘出来，偏置项（由标记为1的节点表示）也被显示。偏置项是数值常数，它允许指示节点的值上下移动，类似于线性方程中的截距。

提示

拥有一个隐藏节点的神经网络可以被看作是我们在第六章中学习的线性回归模型的远亲，预测数值数据 – 回归方法。每个输入节点和隐藏节点之间的权重类似于回归系数，而偏置项的权重类似于截距。

在图的底部，R 报告了训练步骤的数量和一个名为平方误差和（SSE）的误差度量，正如你所预期的那样，它是预测值减去实际值的平方和。较低的 SSE 意味着更好的预测性能。这有助于估计模型在训练数据上的表现，但不能告诉我们它在未见数据上的表现。

第 4 步 – 评估模型性能

网络拓扑图为我们提供了一窥人工神经网络黑匣子的视角，但并未提供关于模型如何拟合未来数据的详细信息。要在测试数据集上生成预测，我们可以使用如下的compute()函数：

> model_results <- compute(concrete_model, concrete_test[1:8])

compute()函数的工作方式与我们迄今使用的predict()函数有所不同。它返回一个包含两个组件的列表：$neurons，用于存储网络中每层的神经元，以及$net.result，用于存储预测值。我们将需要后者：

> predicted_strength <- model_results$net.result

因为这是一个数值预测问题，而不是分类问题，我们不能使用混淆矩阵来检查模型的准确性。相反，我们必须衡量我们预测的混凝土强度与真实值之间的相关性。这提供了两个变量之间线性关联强度的洞察。

请记住，cor()函数用于获取两个数值向量之间的相关性：

> cor(predicted_strength, concrete_test$strength)
 [,1]
[1,] 0.8064655576

提示

如果你的结果不同，不要惊慌。由于神经网络以随机权重开始，预测结果可能因模型而异。如果你想完全匹配这些结果，请尝试在构建神经网络之前使用set.seed(12345)。

接近 1 的相关性表明两个变量之间有很强的线性关系。因此，这里约为 0.806 的相关性表明了相当强的关系。这意味着我们的模型做得相当不错，即使只有一个隐藏节点。

鉴于我们只使用了一个隐藏节点，我们很可能可以提高模型的性能。让我们试着再做得更好一些。

第 5 步 – 改进模型性能

鉴于更复杂的拓扑结构的网络能够学习更复杂的概念，让我们看看当我们将隐藏节点数增加到五时会发生什么。我们像之前一样使用neuralnet()函数，但添加了hidden = 5参数：

> concrete_model2 <- neuralnet(strength ~ cement + slag +
 ash + water + superplastic +
 coarseagg + fineagg + age,
 data = concrete_train, hidden = 5)

再次绘制网络，我们看到连接数量大幅增加。我们可以看到这如何影响了性能，如下所示：

> plot(concrete_model2)

注意，报告的误差（再次由 SSE 测量）已从前一个模型的 5.08 降至 1.63。此外，训练步骤的数量从 4,882 增加到 86,849，这并不奇怪，因为模型变得更加复杂。更复杂的网络需要更多的迭代来找到最优权重。

应用相同的步骤将预测值与真实值进行比较，我们现在得到的相关性约为 0.92，相比之前单个隐藏节点得到的 0.80，这是一个相当可观的改进：

> model_results2 <- compute(concrete_model2, concrete_test[1:8])
> predicted_strength2 <- model_results2$net.result
> cor(predicted_strength2, concrete_test$strength)
 [,1]
[1,] 0.9244533426

有趣的是，在最初的出版物中，Yeh 报告使用一个非常相似的神经网络获得了 0.885 的平均相关性。这意味着，凭借相对较少的努力，我们就能够匹配专家的表现。如果你想更深入地练习神经网络，可能可以尝试应用本章早些时候学到的原理，看看它如何影响模型性能。也许可以尝试使用不同数量的隐藏节点，应用不同的激活函数等等。?neuralnet帮助页面提供了更多关于可调参数的信息。

理解支持向量机

支持向量机（SVM）可以被想象为一个表面，它在多维空间中为数据点创建边界，这些数据点代表示例及其特征值。SVM 的目标是创建一个平坦的边界，称为超平面，该超平面将空间划分为相对均匀的两侧。通过这种方式，SVM 学习结合了第三章中介绍的基于实例的最近邻学习，懒惰学习——使用最近邻进行分类，以及第六章中描述的线性回归建模，预测数值数据——回归方法。这种结合极为强大，使 SVM 能够建模高度复杂的关系。

尽管驱动 SVM 的基础数学已经存在几十年，但它们最近爆发式地流行起来。当然，这源于它们的先进性能，但也可能与获奖的 SVM 算法在多个流行且得到良好支持的库中实现有关，涵盖了包括 R 在内的多种编程语言。因此，SVM 被更广泛的用户群体采用，原本这些用户可能无法应用实现 SVM 所需的复杂数学。好消息是，尽管数学可能很困难，但基本概念是可以理解的。

SVM 可以适应几乎任何类型的学习任务，包括分类和数值预测。该算法的许多关键成功案例出现在模式识别中。显著的应用包括：

• 在生物信息学领域，对微阵列基因表达数据进行分类，以识别癌症或其他遗传疾病

• 文本分类，例如识别文档使用的语言或根据主题对文档进行分类

• 检测稀有但重要的事件，如燃烧发动机故障、安全漏洞或地震

SVM 最容易理解当它用于二元分类时，这也是该方法传统应用的方式。因此，在剩余的部分中，我们将仅专注于 SVM 分类器。但不要担心，因为你在这里学到的原则同样适用于将 SVM 适应到其他学习任务，比如数值预测。

超平面分类

正如前面所述，SVM 使用一种称为超平面的边界来将数据分成类似类值的群体。例如，下图展示了在二维和三维中分隔圆圈和方块的超平面。因为圆圈和方块可以通过直线或平面完美分开，所以它们被称为线性可分。首先，我们将仅考虑这种情况，在其中这是真实的，但 SVM 也可以扩展到数据点不是线性可分的问题。

小贴士

为了方便起见，在 2D 空间中，超平面通常被描绘为一条线，但这只是因为在超过两个维度的空间中描绘空间是困难的。实际上，超平面是高维空间中的一个平面——这是一个可能很难理解的概念。

在二维空间中，SVM 算法的任务是识别一个分隔两个类别的直线。如下图所示，在圆圈和方块的群体之间有多种分隔线的选择。这些可能性标记为a、b和c。算法是如何选择的？

那个问题的答案涉及寻找最大间隔超平面（MMH），它在两个类别之间创建了最大的分离。虽然分隔圆圈和方块的三条线都可以正确分类所有数据点，但最大间隔的线可能对未来数据的泛化效果最好。最大间隔将提高即使在随机噪声的情况下，点仍然保持在边界的正确一侧的机会。

支持向量（在接下来的图中用箭头标出）是每个类别中距离 MMH 最近的点；每个类别必须至少有一个支持向量，但可以有多个。仅使用支持向量，就可以定义 MMH。这是 SVM 的一个关键特征；支持向量提供了一种非常紧凑的方式来存储分类模型，即使特征数量极大。

识别支持向量的算法依赖于向量几何，并涉及一些超出本书范围的相当棘手的数学。然而，该过程的基本原理是相当简单的。

注意

有关支持向量机（SVM）数学原理的更多信息，可以参考经典论文：Cortes C, Vapnik V. Support-vector network. Machine Learning. 1995; 20:273-297。初学者级别的讨论可以参考：Bennett KP, Campbell C. Support vector machines: hype or hallelujah. SIGKDD Explorations. 2003; 2:1-13。更深入的探讨可以参阅：Steinwart I, Christmann A. Support Vector Machines. New York: Springer; 2008。

线性可分数据的案例

在假设类别是线性可分的情况下，最容易理解如何找到最大间隔。在这种情况下，最大间隔超平面（MMH）距离两组数据点的外边界尽可能远。这些外边界被称为 凸包。最大间隔超平面是两条凸包之间最短线段的垂直平分线。使用一种叫做 二次优化 的技术的复杂计算机算法能够通过这种方式找到最大间隔。

另一种（但等效的）方法涉及通过搜索每个可能的超平面空间，以找到一对将数据点分成同质组的平行平面，并且这两个平面之间的距离尽可能远。用一个比喻来说，可以把这个过程想象成在试图找到一个能顺利通过楼梯到卧室的最厚的床垫。

为了理解这个搜索过程，我们需要准确地定义什么是超平面。在 n 维空间中，使用以下方程：

如果你不熟悉这种符号，上面字母上的箭头表示它们是向量，而不是单个数字。特别地，w 是一个 n 维的权重向量，即 {w[1], w[2], ..., w[n]}，而 b 是一个单一数字，称为 偏置。偏置在概念上等同于 第六章 中讨论的斜截式方程中的截距项，数值数据预测 - 回归方法。

提示

如果你很难想象这个平面，不要担心细节。只需将这个方程视为指定一个表面的方式，就像在二维空间中使用斜截式方程 (y = mx + b) 来指定直线一样。

使用这个公式，过程的目标是找到一组权重，指定两个超平面，如下所示：

我们还需要要求这些超平面是按以下方式指定的：使得一个类别的所有点都位于第一个超平面之上，而另一个类别的所有点都位于第二个超平面之下。只要数据是线性可分的，这是可能的。

向量几何定义了这两个平面之间的距离如下：

在这里，||w|| 表示欧几里得范数（即原点到向量 w 的距离）。由于 ||w|| 位于分母，为了最大化距离，我们需要最小化 ||w||。这个任务通常会被重新表达为一组约束条件，如下所示：

虽然这看起来有些混乱，但实际上从概念上理解并不太复杂。基本上，第一行表示我们需要最小化欧几里得范数（平方后除以 2 以简化计算）。第二行表示这是受约束（s.t.）的，即每个 y[i] 数据点必须被正确分类。请注意，y 表示类别值（转化为 +1 或 -1），而倒立的 "A" 是 "对于所有" 的简写。

与寻找最大间隔的另一种方法一样，解决这个问题的任务最好交给二次优化软件。尽管这可能会占用大量处理器资源，但专门的算法能够迅速解决这些问题，即使是在相当大的数据集上。

非线性可分数据的情况

在我们分析了 SVM 背后的理论后，你可能会想知道一个问题：如果数据不是线性可分的，会发生什么？解决这个问题的方法是使用松弛变量，它创建了一个软间隔，允许一些数据点落在间隔的错误一侧。接下来的图示说明了两个数据点落在直线的错误一侧，并显示了相应的松弛项（用希腊字母 Xi 表示）：

对所有违反约束的点应用一个成本值（记为 C），与其寻找最大间隔，算法会尝试最小化总成本。因此，我们可以将优化问题修改为：

如果你仍然感到困惑，不用担心，你不是唯一一个。幸运的是，SVM 软件包会很高兴地为你优化这一过程，而无需你理解技术细节。需要理解的重要部分是成本参数 C的引入。修改该值将调整惩罚，例如，数据点落在超平面错误的一侧。成本参数越大，优化过程越会努力实现 100% 的分类准确率。另一方面，较低的成本参数则会将重点放在更宽的整体间隔上。为了创建一个能很好地泛化到未来数据的模型，平衡这两者是非常重要的。

使用核函数处理非线性空间

在许多实际应用中，变量之间的关系是非线性的。正如我们刚刚发现的，SVM 仍然可以通过添加松弛变量来训练这类数据，从而允许某些示例被误分类。然而，这并不是处理非线性问题的唯一方式。SVM 的一个关键特点是它们能够使用一种称为 核技巧 的过程将问题映射到一个更高维度的空间。在这样做的过程中，非线性关系可能突然变得非常线性。

尽管这看起来像是废话，但实际上通过示例很容易说明。在下图中，左侧的散点图描绘了天气类别（晴天或雪天）与两个特征：纬度和经度之间的非线性关系。图中心的点属于雪天类别，而图边缘的点都是晴天。这些数据可能来自一组天气报告，其中一些来自山顶附近的气象站，另一些则来自山脚下的气象站。

在图的右侧，应用了核技巧后，我们通过一个新维度——海拔高度来观察数据。随着这个特征的加入，类别现在可以完美地线性分开。这之所以可能，是因为我们从一个新的视角来看待数据。在左图中，我们是从鸟瞰视角观察这座山，而在右图中，我们是从地面远处观察这座山。在这里，趋势非常明显：雪天出现在更高的海拔。

使用非线性核的 SVM 通过增加数据的额外维度来创造分离。实际上，核技巧包括构建新特征的过程，这些特征表达了度量特征之间的数学关系。例如，海拔特征可以数学地表达为纬度和经度之间的交互作用——该点离这两个尺度的中心越近，海拔就越高。这使得 SVM 能够学习到原始数据中没有明确度量的概念。

使用非线性核的 SVM 是非常强大的分类器，尽管它们确实有一些缺点，如下表所示：

优势	弱点

|

• 可用于分类或数值预测问题

• 不容易受到噪声数据的影响，也不容易发生过拟合

• 可能比神经网络更容易使用，尤其是因为存在多个得到良好支持的 SVM 算法

• 由于其高准确率和在数据挖掘竞赛中的高调获胜，SVM 正在获得越来越多的关注

|

• 寻找最佳模型需要测试不同核函数和模型参数的组合

• 训练可能较慢，特别是当输入数据集具有大量特征或示例时

• 结果是一个复杂的黑箱模型，难以解释，甚至可能无法解释

|

核函数通常具有以下形式。由希腊字母 phi 表示的函数，即 ϕ(x)，是将数据映射到另一个空间。因此，一般的核函数对特征向量 x[i] 和 x[j] 进行某种变换，并使用 点积将它们结合，点积操作将两个向量转换为一个单一的数字。

使用这种形式，已经为许多不同领域的数据开发了核函数。以下列出了几种最常用的核函数。几乎所有的 SVM 软件包都会包含这些核函数及其他许多核函数。

线性核完全不对数据进行转换。因此，它可以简单地表示为特征的点积：

多项式核的度数为 d，对数据进行简单的非线性变换：

Sigmoid 核导致的 SVM 模型与使用 sigmoid 激活函数的神经网络有些类似。希腊字母 kappa 和 delta 被用作核参数：

高斯 RBF 核类似于 RBF 神经网络。RBF 核在许多类型的数据上表现良好，且被认为是许多学习任务的合理起点：

没有可靠的规则可以将特定的核函数与某个学习任务匹配。匹配程度很大程度上依赖于要学习的概念、训练数据的数量以及特征之间的关系。通常，需要通过在验证数据集上训练和评估多个 SVM 来进行一些试错法。话虽如此，在许多情况下，核函数的选择是任意的，因为性能可能会有所波动。为了实践这种方法，让我们将 SVM 分类应用到一个现实问题中。

示例 – 使用 SVM 进行 OCR（光学字符识别）

图像处理是许多类型的机器学习算法面临的困难任务。像素模式与更高层次概念之间的关系极为复杂且难以定义。例如，人类很容易识别出一张脸、一只猫或字母 "A"，但将这些模式定义为严格的规则却十分困难。此外，图像数据往往存在噪声。图像的捕获方式可能因光线、方向和主体位置的不同而有所变化。

支持向量机（SVM）非常适合应对图像数据的挑战。它们能够学习复杂的模式而不对噪声过于敏感，从而能够高精度地识别视觉模式。此外，SVM 的一个关键弱点——黑箱模型表示——对于图像处理的影响较小。如果一个 SVM 能区分猫和狗，那么它是如何做到的就不那么重要了。

在本节中，我们将开发一个类似于光学字符识别（OCR）软件中核心使用的模型，这类软件通常与桌面文档扫描仪捆绑在一起。此类软件的目的是通过将打印或手写文本转换为电子形式来处理纸质文档，以便保存到数据库中。当然，由于手写风格和打印字体的多样性，这是一个困难的问题。尽管如此，软件用户仍然期望完美，因为错误或错别字可能在商业环境中导致尴尬或代价高昂的错误。让我们看看我们的 SVM 是否能够完成这项任务。

步骤 1 – 收集数据

当 OCR 软件首次处理文档时，它会将纸张分割成一个矩阵，使得网格中的每个单元格包含一个单独的字形，字形只是指一个字母、符号或数字的术语。接下来，软件将尝试将每个单元格的字形与它识别的所有字符集合进行匹配。最后，单个字符将被重新组合成单词，并且可以选择通过字典对文档语言中的单词进行拼写检查。

在本次练习中，我们假设我们已经开发了算法，将文档划分为每个包含单个字符的矩形区域。我们还假设文档仅包含英文字母字符。因此，我们将模拟一个过程，涉及将字形与 26 个字母（A 到 Z）中的一个进行匹配。

为此，我们将使用 W. Frey 和 D. J. Slate 捐赠给 UCI 机器学习数据集库（archive.ics.uci.edu/ml）的数据集。该数据集包含 20,000 个示例，展示了使用 20 种不同随机形变和失真的黑白字体打印的 26 个英文字母的大写字母。

注意

有关此数据集的更多信息，请参考 Slate DJ, Frey W.《使用霍兰德风格适应分类器的字母识别*》。机器学习.* 1991; 6:161-182。

以下图像由 Frey 和 Slate 发布，提供了一些打印字形的示例。这些字形被扭曲后，计算机很难识别，但人类却能轻松识别：

步骤 2 – 探索和准备数据

根据 Frey 和 Slate 提供的文档，当字形被扫描到计算机中时，它们会被转换为像素，并记录 16 个统计属性。

这些属性测量了字形的水平和垂直维度、黑色（与白色）像素的比例，以及像素的平均水平和垂直位置。可以推测，字形中黑色像素在不同区域的浓度差异应能提供一种区分 26 个字母的方法。

提示

为了跟随这个例子，请从 Packt Publishing 网站下载letterdata.csv文件，并将其保存到你的 R 工作目录中。

将数据读入 R 后，我们确认已经接收到包含定义每个字母类别的 16 个特征的数据。正如预期的那样，字母类有 26 个水平：

> letters <- read.csv("letterdata.csv")
> str(letters)
'data.frame':  20000 obs. of 17 variables:
 $ letter: Factor w/ 26 levels "A","B","C","D",..
 $ xbox  : int  2 5 4 7 2 4 4 1 2 11 ...
 $ ybox  : int  8 12 11 11 1 11 2 1 2 15 ...
 $ width : int  3 3 6 6 3 5 5 3 4 13 ...
 $ height: int  5 7 8 6 1 8 4 2 4 9 ...
 $ onpix : int  1 2 6 3 1 3 4 1 2 7 ...
 $ xbar  : int  8 10 10 5 8 8 8 8 10 13 ...
 $ ybar  : int  13 5 6 9 6 8 7 2 6 2 ...
 $ x2bar : int  0 5 2 4 6 6 6 2 2 6 ...
 $ y2bar : int  6 4 6 6 6 9 6 2 6 2 ...
 $ xybar : int  6 13 10 4 6 5 7 8 12 12 ...
 $ x2ybar: int  10 3 3 4 5 6 6 2 4 1 ...
 $ xy2bar: int  8 9 7 10 9 6 6 8 8 9 ...
 $ xedge : int  0 2 3 6 1 0 2 1 1 8 ...
 $ xedgey: int  8 8 7 10 7 8 8 6 6 1 ...
 $ yedge : int  0 4 3 2 5 9 7 2 1 1 ...
 $ yedgex: int  8 10 9 8 10 7 10 7 7 8 ...

记住，SVM 学习器要求所有特征都是数值型的，而且每个特征都要缩放到一个相对较小的区间。在这种情况下，每个特征都是整数，因此我们不需要将任何因子转换为数字。另一方面，这些整数变量的某些范围似乎相当宽泛，这表明我们需要对数据进行归一化或标准化。然而，我们现在可以跳过这一步，因为我们将使用的 R 包在拟合 SVM 模型时会自动执行重新缩放。

由于数据准备工作已经基本完成，我们可以直接进入机器学习过程的训练和测试阶段。在之前的分析中，我们随机地将数据分为训练集和测试集。尽管我们可以在这里这样做，但 Frey 和 Slate 已经对数据进行了随机化，因此建议使用前 16,000 条记录（80%）来构建模型，使用接下来的 4,000 条记录（20%）进行测试。根据他们的建议，我们可以按以下方式创建训练和测试数据框：

> letters_train <- letters[1:16000, ]
> letters_test  <- letters[16001:20000, ]

数据准备好后，让我们开始构建分类器。

步骤 3 – 在数据上训练模型

在 R 中拟合 SVM 模型时，有几个出色的包可供选择。来自维也纳大学（TU Wien）统计学系的e1071包提供了一个 R 接口，连接到获奖的 LIBSVM 库，这是一个广泛使用的用 C++编写的开源 SVM 程序。如果你已经熟悉 LIBSVM，你可能希望从这里开始。

注意

有关 LIBSVM 的更多信息，请访问作者的官方网站：www.csie.ntu.edu.tw/~cjlin/libsvm/。

同样，如果你已经使用 SVMlight 算法，来自多特蒙德大学（TU Dortmund）统计学系的klaR包提供了直接在 R 中使用该 SVM 实现的功能。

注意

关于 SVMlight 的信息，可以查看svmlight.joachims.org/。

最后，如果你是从头开始，最好从kernlab包中的 SVM 函数开始。这个包的一个有趣的优点是它是原生用 R 开发的，而不是用 C 或 C++，这使得它非常容易自定义；包的内部实现完全透明，不会隐藏任何内容。或许更重要的是，与其他选项不同，kernlab可以与caret包一起使用，这使得 SVM 模型能够使用各种自动化方法进行训练和评估（在第十一章中有介绍，提升模型性能）。

注意

如果你想更详细地了解kernlab，请参阅作者的论文，地址为www.jstatsoft.org/v11/i09/。

使用kernlab训练 SVM 分类器的语法如下。如果你确实使用的是其他包，命令大体是类似的。默认情况下，ksvm()函数使用高斯 RBF 核，但也提供了其他多种选项。

为了提供 SVM 性能的基准衡量标准，我们从训练一个简单的线性 SVM 分类器开始。如果你还没有安装，可以通过install.packages("kernlab")命令将kernlab包安装到你的库中。然后，我们可以在训练数据上调用ksvm()函数，并使用vanilladot选项指定线性（即基础）核，如下所示：

> library(kernlab)
> letter_classifier <- ksvm(letter ~ ., data = letters_train,
 kernel = "vanilladot")

根据你计算机的性能，这个操作可能需要一些时间才能完成。当完成时，输入已存储模型的名称，以查看一些关于训练参数和模型拟合的基本信息。

> letter_classifier
Support Vector Machine object of class "ksvm" 

SV type: C-svc  (classification) 
 parameter : cost C = 1 

Linear (vanilla) kernel function. 

Number of Support Vectors : 7037 

Objective Function Value : -14.1746 -20.0072 -23.5628 -6.2009 -7.5524 -32.7694 -49.9786 -18.1824 -62.1111 -32.7284 -16.2209...

Training error : 0.130062

这些信息对我们了解模型在实际应用中的表现几乎没有帮助。我们需要检查它在测试数据集上的表现，才能知道它是否能够很好地推广到未见过的数据。

第 4 步 – 评估模型性能

predict()函数允许我们使用字母分类模型对测试数据集进行预测：

> letter_predictions <- predict(letter_classifier, letters_test)

因为我们没有指定类型参数，所以使用了默认值type = "response"。这将返回一个向量，其中包含测试数据中每一行值对应的预测字母。通过使用head()函数，我们可以看到前六个预测字母分别是U、N、V、X、N和H：

> head(letter_predictions)
[1] U N V X N H
Levels: A B C D E F G H I J K L M N O P Q R S T U V W X Y Z

为了检查我们的分类器表现如何，我们需要将预测的字母与测试数据集中的真实字母进行比较。我们将使用table()函数来完成这一任务（这里只显示了完整表格的一部分）：

> table(letter_predictions, letters_test$letter)
letter_predictions   A   B   C   D   E
 A 144   0   0   0   0
 B   0 121   0   5   2
 C   0   0 120   0   4
 D   2   2   0 156   0
 E   0   0   5   0 127

144、121、120、156和127这些对角线上的值表示预测字母与真实值匹配的记录总数。类似地，错误的数量也列出了。例如，B行和D列中值为5表示有 5 个案例将字母D错误地识别为B。

分别查看每种类型的错误可能会揭示一些关于模型在处理特定字母时遇到的困难的有趣模式，但这需要大量时间。我们可以通过计算整体准确率来简化评估过程。这只考虑预测是否正确，而忽略错误的类型。

以下命令返回一个TRUE或FALSE值的向量，表示模型预测的字母是否与测试数据集中实际的字母一致（即，是否匹配）：

> agreement <- letter_predictions == letters_test$letter

使用table()函数，我们看到分类器在 4000 个测试记录中的 3357 个记录上正确识别了字母：

> table(agreement)
agreement
FALSE  TRUE
643 3357

按百分比计算，准确率大约为 84%：

> prop.table(table(agreement))
agreement
 FALSE    TRUE
0.16075 0.83925

请注意，当 Frey 和 Slate 在 1991 年发布数据集时，他们报告的识别准确率大约为 80%。仅通过几行 R 代码，我们就超越了他们的结果，尽管我们也得到了二十多年额外机器学习研究的好处。考虑到这一点，我们可能能够做得更好。

第五步 – 提高模型性能

我们之前的 SVM 模型使用了简单的线性核函数。通过使用更复杂的核函数，我们可以将数据映射到更高维的空间，并可能获得更好的模型拟合。

然而，从众多不同的核函数中选择合适的一个可能是具有挑战性的。一种常见的做法是从高斯 RBF 核开始，因为它已被证明在许多类型的数据上表现良好。我们可以像这样使用ksvm()函数训练基于 RBF 的 SVM：

> letter_classifier_rbf <- ksvm(letter ~ ., data = letters_train,
 kernel = "rbfdot")

接下来，我们像之前一样进行预测：

> letter_predictions_rbf <- predict(letter_classifier_rbf,
 letters_test)

最后，我们将与线性 SVM 的准确率进行比较：

> agreement_rbf <- letter_predictions_rbf == letters_test$letter
> table(agreement_rbf)
agreement_rbf
FALSE  TRUE

  275  3725
> prop.table(table(agreement_rbf))
agreement_rbf
 FALSE    TRUE
0.06875 0.93125

提示

由于ksvm RBF 核函数的随机性，你的结果可能与这里展示的有所不同。如果你希望它们完全匹配，可以在运行 ksvm()函数之前使用set.seed(12345)。

通过简单地更改核函数，我们成功将字符识别模型的准确率从 84%提高到了 93%。如果 OCR 程序的表现仍然不令人满意，可以尝试其他核函数，或者调整约束参数 C 的成本，以修改决策边界的宽度。作为练习，你应该尝试这些参数，看看它们如何影响最终模型的成功。

总结

在本章中，我们考察了两种机器学习方法，它们具有巨大的潜力，但由于复杂性经常被忽视。希望你现在能看到，这种声誉至少在某种程度上是不公平的。驱动 ANNs 和 SVM 的基本概念相当容易理解。

另一方面，由于人工神经网络（ANNs）和支持向量机（SVMs）已经存在了几十年，每种方法都有许多变种。本章只是简单介绍了这些方法可能实现的部分内容。通过运用你在这里学到的术语，你应该能够掌握那些每天都在发展的各种进展之间的细微差别。

现在我们已经花了一些时间学习了从简单到复杂的多种预测模型；在下一章，我们将开始考虑其他类型学习任务的方法。这些无监督学习技术将揭示数据中的迷人模式。