本节在一些较大的数据集上使用相应恰当的损失函数来演示梯度提升方法。
10.14.1 加利福尼亚州房屋数据
这个数据集(Pace and Barry, 1997)可从卡内基美隆统计实验室(Carnegie-Mellon Statlib)的资源库获取1。其中包含了加利福尼亚州的 20,460 个街区(1990 年的调查区域)的聚合数据。输出变量 $Y$ 为在每个街区内房屋价值的中位数,以十万美元为单位。其中包括一些人口统计学的自变量,比如人口数量(population)、收入中位数(MedInc)、用房屋个数来代表的房屋密度(House)和每个房屋的平均居住率(AveOccup)。其中还有表示了每个街区位置的自变量(longitude 和 latitude),以及一些反应了街区中房屋属性的数量:平均房间个数(AveRooms)和平均卧室个数(AveBedrms)。因此一共有八个自变量,全部都为数值变量。
我们用 MART 方法拟合了一个梯度提升模型,采用了 $J=6$ 个终节点、$\nu=0.1$ 的学习率(式 10.41)和预测数值输出变量的 Huber 损失准则。数据集被随机地分成训练集(80%)和测试集(20%)。
图 10.13 展示了在训练集和测试集上作为循环次数 $M$ 的函数的平均绝对误差(average absolute error):
$$\text{AAE} = \operatorname{E} |y - \hat{f}_M(x)| \tag{10.53}$$看起来测试误差随着 $M$ 的上升单调地下降,在前期阶段下降的更快速,然后随着循环次数的上升逐渐趋于平稳到接近于常数。因此,$M$ 的特定选择并不太重要,只要其取值不要太小。在很多应用中都通常如此。收缩方法(式 10.41)一般可以消除过拟合的问题,尤其是在较大的数据集中。
在 800 次循环后 AAE 的值为 0.31。与之相比,最优的常数预测为中位数 $\{y_i\}$ 的 AAE 为 0.89。以更习惯的指标来说,这个模型的平方多元相关系数为 $R^2 = 0.84$。Pace and Barry (1997) 使用了一个复杂的空间自回归方法,其中对每个街区的预测是基于附近街区的房屋价值中位数和其他的自变量。经过转换后,他们在预测 $\log Y$ 中达到了 $R^2=085$。若使用 $\log Y$ 作为输出变量,梯度提升方法可达到的值为 $R^2=0.86$。
图 10.14 展示了八个自变量的相对重要性。不出所料,街区内的收入中位数是最相关的自变量。经度、纬度和平均居住率的相关程度都大致为收入的一半,而其他的变量就不太有影响了。
图 10.15 展示了最相关的非位置自变量的单变量部分依赖图。注意到曲线图不是平滑的。这是使用了树模型所产生的结果。决策树产生的是不连续的分段常数模型(式 10.25)。这个性质在树模型的求和中(表达式 10.28)也同样存在,当然分段的个数也会更多。与本书中讨论的大多数其他方法不同,这里对结果没有平滑程度上的约束。模型中可以包含任意急剧的不连续性。这里的曲线还是呈现了大致上的平滑趋势,是因为这个问题的输出变量最佳预测的估计就是如此的。多数问题中通常会是这样的结果。
每个图底部的竖线标记描绘了相应变量在数据中的十分位数。注意这里的数据密度在边缘区域较低,尤其是在取值较大的区域这使曲线在那些区域上的准确性相对较差。所有图的纵轴尺度是一样的,所以可对不同变量的相对重要性进行视觉上的比较。
收入中位数对房屋价值中位数的部分依赖是单调递增的,在整个的数据集中几乎是线性的关系。除了在小于一的区域外,房屋价值大致是随着平均居住率的上升而单调递减的。对平均房间数量,房屋价值中位数呈现出非单调的部分依赖。房屋价值在大约三个房间处达到最小值,在增大或降低两个方向都会增加。
从图中看起来,房屋价值中位数对房龄只有很弱的部分依赖,这与其相对重要性的排序并不相符(图 10.14)。这说明了这个若主效应可能屏蔽了其与其他变量之间的强交互效应。图 10.16 展示了房屋价值中位数对房龄中位数和平均居住率联合取值的二变量部分依赖。很明显这两个变量之间存在交互作用。当平均居住率大于二时,房屋价值几乎与房龄中位数是独立的,而当小于二时,房屋价值对房龄有较强的依赖关系。
图 10.17 展示了拟合模型对经度和纬度联合取值的二变量部分依赖,表示为阴影的等值线图。很明显可看出房屋价值中位数对街区位置有较强的依赖关系。注意图 10.17 并不是不考虑其他自变量作用下的房屋价值与位置的分布图(式 10.49)。与所有的部分依赖图一样,它表示的是在已考虑了其他的街区和房屋特征的作用之后位置因素的影响(式 10.47)。可将它视为是对位置因素而付出的额外溢价。看起来这个溢价在太平洋海岸线上,尤其是在湾区(Bay Area)和洛杉矶及圣地亚哥(Los Angeles - San Diego)区域,相对更高。在加利福尼亚州北部、中央谷地和东南沙漠地区,位置因素的成本要低很多。
10.14.2 新西兰渔业数据
植物和动物生态学家使用回归模型来预测作为环境变量函数的物种存在(presence)、丰度(abundance)和多度(richness)。尽管简单的线性和参数模型流行了很多年,该领域近期开始对更复杂的模型越来越感兴趣,比如广义加性模型(第 9.1 节,GAM)、多元自适应回归样条(第 9.4 节,MARS)和提升回归树模型(Leathwick et al., 2005;Leathwick et al., 2006)。这里对黑异海鲂(Black Oreo Dory),一种在新西兰周围的海洋水域中发现的海洋鱼类,的存在和丰度进行建模2。
图 10.18 展示了 17,000 次捕捞(深海网补,最大深度为 2 千米)的位置,红色点代表了 2353 次存在黑异海鲂的捕捞。有一百多个物种会定期记录,黑异海鲂为其中之一。在每次捕捞中会记录每个物种的捕获公斤数量。在每次捕捞中,除了捕获物种数据外,还会记录一些环境变量的测量。这里包括了捕捞的平均深度(AvgDepth)和水域的温度和盐度。由于后两个变量与深度是高度相关的,Leathwick et al. (2006) 衍生出了替代的变量TempResid 和 SalResid,即(分别使用非参数回归)调整了深度后的温度和盐度的残差。SSTGrad 为海洋表面温度梯度的一个测量,Chla 为通过卫星图像测量的生态系统生产力的一个广泛指标。SusPartMatter 为在沿海岸水域的悬浮颗粒物(suspended particulate matter)的测量,也是由卫星图像得出的。
这个分析的目标是在对根据捕捞速度、捕捞距离和捕捞网目尺寸(mesh size)标准化后,估计出在一次捕捞中存在黑异海鲂的概率,以及期望的捕获量。作者在估计概率中使用了对数几率(logistic)回归。对于捕获量,可能一个自然的想法是假设其服从泊松(Poisson)分布,然后对期望捕获量的对数建模,但由于有大量的 0 取值所以这个方法通常不合适。尽管已有一些专用的方法,例如零膨胀(zero-inflated)泊松模型(Lambert, 1992),作者还是选择了更简单的方法。若 $Y$ 为(非负)捕获量,
$$\operatorname{E}(Y|X) = \operatorname{E}(Y | Y>0, X) \cdot \operatorname{Pr}(Y>0 | X) \tag{10.54}$$其中第二项由对数几率回归来估计,第一项仅用 2353 次有正捕获量的捕捞样本来估计。
在对数几率回归中,作者使用了梯度提升模型(GBM)3,采用了二项偏差损失函数,树模型的深度为 10,收缩因子为 $\nu=0.025$。在正的捕捉量回归中,作者使用了 GBM 对 $\log(Y)$ 建模,采用了平方误差损失(树模型深度也为 10,但 $\nu=0.01$),自变量为没有取对数的值。作者在两者中都使用了十折交叉验证来选择模型项的个数以及收缩因子。
图 10.19(左侧)展示了 GBM 的模型序列的二项偏差均值,包括了十折交叉验证以及在测试集上的结果。这相比于在每项 8 个自由度的平滑样条拟合的广义加性模型(GAM)有少许的改进。右侧展示了两个模型的衡量预测表现的 ROC 曲线(见第 9.2.5 节),从这个图上看,两者的表现非常相似,可能 GBM 在用 AUC 比较时有一点点优势。在灵敏度(sensitivity)和特异度(specificity)等值点处,AM 为 90%。
图 10.20 概括了对数几率 GBM 拟合中变量的贡献度。从中可见捕获到黑异海鲂有明显的深度区间,而且在较冷的水域中捕获更加频繁。这里没有给出捕获量的数量模型的细节;其中的重要变量大致相同。
这个模型中所有的自变量在一个很精细的地理网格上可以获取到;实际上,它们是从环境地图位置、卫星图像和其他类似来源衍生而来的,更多细节见 Leathwick et al. (2006)。这意味着在这个网格的位置上可以给出预测,并且输入到 GIS 地图系统中。图 10.21 展示了对捕捞条件标准化后的存在和捕获量的预测地图;由于自变量根据地理位置的变化是连续的,所以预测值也是连续变化的。
由于其可以包含交互项、自动选择变量并且对异常值和缺失值稳健,GBM 模型很快在这个数据量大而且活跃的学术群体中流行起来。
10.14.3 人口统计数据
本节在多类别的分类问题上使用 MART 来演示梯度提升方法。数据来自旧金山湾区的商场顾客填写的 9243 个调查问卷(俄亥俄州哥伦布市 Impact Resources 公司)。其中有 14 个是关于人口统计学的问题。在这个演示中,目标是用其他 13 个变量作为自变量来预测职业,也因此可识别出能区分出不同职业的人口统计学变量。数据被随机地分成训练集(80%)和测试集(20%),采用了 $J=6$ 个终节点的树以及学习率 $\nu=0.1$。
图 10.22 展示了 $K=9$ 个职位类别以及其对应的错误率。总体错误率为 42.5%,与之相比,用样本最多的类别“Prof/Man”(专家/管理)作为预测的错误率为 69%。四个预测的最好的类别为“Retired”(退休)、“Student”(学生)、“Prof/Man”(专家/管理)和“Homemaker”(料理家务)。
图 10.23 展示了在所有类别上平均的自变量相对重要性(式 10.46)。图 10.24 展示了四个预测的最好的类别各自的自变量相关重要性(式 10.45)。可见每个类别中最相关的自变量一般都不同。一个例外是“age”(年龄),在预测“Retired”(退休)、“Student”(学生)和“Prof/Man”(专家/管理)时都出现在最相关的前三个变量中。
图 10.25 展示了这三个类别的对数几率(式 10.52)对“age”(年龄)的部分依赖。横轴上的值为相应的等距年龄区间上的有序编码。可见在考虑了其他变量的贡献之后,年龄越大的人退休的几率越高,而年龄越大的人是学生的几率越低。专家/管理职业的几率在中年人中最高。显然这些结果是情理之中的。这个示例说明了对每个类别分别考察其部分依赖可得出切合实际的结论。