9.4 多元自适应回归样条

多元自适应回归样条(multivariate adaptive regression splines,MARS) 是一个回归问题的自适应过程,非常适用于高维问题中(即有大量的输入变量)。可将它视为一个对逐步线性回归的推广,或对 CART 的增强在回归问题中表现的改进。本节先从第一个角度介绍 MARS,之后再建立与 CART 的联系。

MARS 中使用的是形式为 $(x-t)_+$ 和 $(t-x)_+$ 的分段线性基函数的展开。加号 “+” 表示了正部:

$$\begin{align} (x-t)_+ &= \begin{cases} x - t & \text{若 } x > t \\ 0 & \text{其他} \end{cases} \\ (t-x)_+ &= \begin{cases} t - x & \text{若 } x < t \\ 0 & \text{其他} \end{cases} \end{align}$$

作为例子,图 9.9 展示了函数 $(x-0.5)_+$ 和 $(0.5-x)_+$。

**图 9.9**:MARS 所使用的基函数 $(x-t)\_+$(橙色实线)和 $(t-x)\_+$(蓝色虚线)。
图 9.9:MARS 所使用的基函数 $(x-t)_+$(橙色实线)和 $(t-x)_+$(蓝色虚线)。

每个函数都是分段线性,在 $t$ 值处有一个结点。以第五章中的定义,这些就是线性样条函数。在后续讨论中,将这两个函数称为 镜像对(reflected pair)。这个方法的思路是为每个输入变量 $X_j$ 构建以该输入变量每个观测值 $x_{ij}$ 为结点的镜像对。因此,基函数的集合为:

$$\mathcal{C}= \{(X_j-t)_+, (t-X_j)_+\}_{ \substack{t \in \{x_{1j}, x_{2j}, \dots, x_{Nj}\} \\ j = 1,2,\dots,p}} \tag{9.18}$$

若所有输入变量的取值都是唯一的,那么共有 $2Np$ 个基函数。注意虽然每个基函数只依赖于某一个 $X_j$,例如 $h(X)=(X_j-t)_+$,但也将其视为整个输入空间 $\mathbb{R}^p$ 上的函数。

模型的生成策略类似于一个前向逐步线性回归,但不局限于原始的输入变量,可以加入集合 $\mathcal{C}$ 中的函数及它们的乘积。因此模型的形式为:

$$f(X) = \beta_0 + \sum_{m=1}^M \beta_m h_m(X) \tag{9.19}$$

其中每个 $h_m(X)$ 为 $\mathcal{C}$ 中的一个函数,或两个或多个这样的函数的乘积。

给定了 $h_m$ 的选择后,通过最小化残差平方和估计出系数 $\beta_m$,即标准的线性回归。然而真正的技巧在于函数组 $h_m(x)$ 的构建。从只包含常数函数 $h_0(X)=1$ 的模型开始,并且集合 $\mathcal{C}$ 中的函数都是备选函数。图 9.10 描述了这个过程。

**图 9.10**:MARS 前向模型生成过程的示意图。左侧为当前模型包含的基函数:初始只有常数函数 $h(X)=1$。右侧为构建模型所考虑的所有备选基函数。这些是如图 9.9 中的分段线性函数对,结点 $t$ 的位置为每个自变量 $X_j$ 的所有观测唯一值 $x_{ij}$。在每个步骤中,考察所有的某个备选函数对与模型中某个基函数的乘积。将对残差平方和降低程度最大的乘积加入到当前模型中。图中演示了这个过程的前三个步骤,用红色标记被选中的函数。
图 9.10:MARS 前向模型生成过程的示意图。左侧为当前模型包含的基函数:初始只有常数函数 $h(X)=1$。右侧为构建模型所考虑的所有备选基函数。这些是如图 9.9 中的分段线性函数对,结点 $t$ 的位置为每个自变量 $X_j$ 的所有观测唯一值 $x_{ij}$。在每个步骤中,考察所有的某个备选函数对与模型中某个基函数的乘积。将对残差平方和降低程度最大的乘积加入到当前模型中。图中演示了这个过程的前三个步骤,用红色标记被选中的函数。

在每个步骤中,将所有的当前模型函数集 $\mathcal{M}$ 中的某个函数 $h_m$ 与 $\mathcal{C}$ 中的某个镜像对的乘积作为新的基函数对。最终将带来最大的训练误差下降的两个函数乘积添加到模型 $\mathcal{M}$ 中,其形式为:

$$\hat{\beta}_{M+1} h_\ell(X) \cdot (X_j - t)_+ + \hat{\beta}_{M+2} h_\ell(X) \cdot (t - X_j)_+, h_\ell \in \mathcal{M}$$

其中的系数 $\hat{\beta}_{M+1}$ 和 $\hat{\beta}_{M+2}$ 连同模型中其他的 $M+1$ 个系数一起通过最小二乘估计得出。将这对乘积加入到模型中,然后持续这个过程,直到模型的函数集 $\mathcal{M}$ 包含的项达到了某个预设最大个数。

举个例子,在第一步中,考虑向模型中添加的函数形式为

$$\beta_1(X_j-t)_+ + \beta_2(t-X_j)_+, t\in\{x_{ij}\}$$

当前模型中只有常数函数,而常数函数与其他函数的乘积就是这些常数的本身。假设最优的选择为

$$\hat{\beta}_1 (X_2 - x_{72})_+ + \hat{\beta}_2 (X_{72} - X_2)_+$$

那么这对基函数被添加到模型基函数集 $\mathcal{M}$。在下一步骤中,考虑如下形式的函数乘积对:

$$h_m(X) \cdot (X_j - t)_+ \text{ 和 } h_m(X) \cdot (t - X_j)_+, t \in \{x_{ij}\}$$

其中的 $h_m$ 可从下列中选择:

$$\begin{align} h_0(X) &= 1, \\ h_1(X) &= (X_2 - x_{72})_+ \text{ 或} \\ h_2(X) &= (x_{72} - X_2)_+ \end{align}$$

选择第三个函数会得到类似 $(X_1-x_{51})_+\cdot (x_{72}-X_2)_+$ 形式的函数乘积,如图 9.11 所展示。

**图 9.11**:由 MARS 分段线性基函数的乘积得出的函数 $(X_1-x_{51})\_+\cdot(x_{72}-X_2)\_+$
图 9.11:由 MARS 分段线性基函数的乘积得出的函数 $(X_1-x_{51})_+\cdot(x_{72}-X_2)_+$

这个过程最终会产生一个形式为式 9.19 的大模型。这个模型通常会对数据过拟合,所以需要进行一个后向删减过程。在每一个步骤,选择当被排除后产生最小的残差平方和增加的项从模型中移除,得到每个大小(项个数)$\lambda$ 下的估计最优模型 $\hat{f}_\lambda$。尽管可以使用交叉验证来估计最优的 $\lambda$ 取值,但为了减少计算量 MARS 过程使用了广义交叉验证。定义这个准则为:

$$\operatorname{GCV}(\lambda) = \frac {\sum_{i=1}^N (y_i- \hat{f}_\lambda(x_i)^2)} {(1 - M(\lambda) / N)^2} \tag{9.20}$$

数值 $M(\lambda)$ 是模型的有效参数个数:这包括了模型中项的个数以及用于选择最优结点位置的参数个数。一些数学推导和模拟的结果表明在分段线性回归中选择每个结点需要花费三个参数的代价。

因此如果模型中有 $r$ 个线性独立的基函数,并在前向的过程中选择了 $K$ 个结点,那么有效参数个数为 $M(\lambda)=r+cK$,其中的 $c=3$。(当模型被约束为加性时——细节见下——使用的惩罚系数是 $c=2$。)代入这个表达式后,可从最小化了 $\operatorname{GCV}(\lambda)$ 的倒序序列中选择最终模型。

为什么要选择这些分段线性基函数,并且为什么用这种特别的建模方式?图 9.9 中的函数的关键性质是它们组合成局部性函数的能力;它们在部分的定义域上取值为零。当它们相乘之后,如图 9.11 所示,得到的函数只在两者都不为零的小部分特征空间区域上是非零的。因此,回归曲面是用局部非零的成分函数一点点地构建出的——每个函数只在被需要的区域起作用。在高维问题中参数个数会随项的加入快速增加,而这个方法的重要性质可更节约地“使用”参数。比如多项式的其他基函数的乘积得到的函数会在整个空间上都不为零,也不会如这个方法一样有效。

分段线性基函数的第二个重要优势在计算量上。考虑 $\mathcal{M}$ 中的一个函数与输入变量 $X_j$ 的所有 $N$ 个镜像对函数的乘积。这需要拟合 $N$ 个单输入变量的线性回归模型,每个拟合需要 $O(N)$ 次运算,所以一共需要 $O(N^2)$ 次运算。然而可以充分利用分段线性函数的简单结构。首先用最右一个位置结点的镜像对进行拟合。然后随着结点向左侧逐步地移动,基函数在定义域的左侧部分保持不变,在右侧部分有一个常数的改动。因此在每个移动后更新拟合只需要 $O(1)$ 次运算。这样仅用 $O(N)$ 次运算就可以拟合所有的结点。

MARS 的前向建模方式是层级式的,也就是说在每一步中构建所考察的多个函数乘积项会使用到当前步骤模型中已存在的函数项。例如在向模型中添加四函数乘积项时,只考虑由模型中已有的三函数乘积构建出的项。这样做的道理是一般只有当某个高阶的交叉项中的低阶的成分存在于模型中时,它才可能也会存在于模型中。现实并不一定如此,但这是一个可操作的合理假设,并且避免了在指数级增长的备选函数空间进行搜寻最优。

对模型中项的形式有一个约束条件:每个输入变量在一个乘积项中只能出现最多一次。这将从模型中排除了在特征空间边界附近剧烈上升或下降的输入变量高阶次幂项。但可用分段线性函数以更稳健的方式来近似这种高次幂项。

MARS 过程中一个有用的选项是设置交叉项阶数的上限。例如可设置上限为 2,只考虑分段线性函数的两两乘积,但不考虑三个或更多函数的乘积。这样可以改进最终模型的可解释性。设置上限为 1 就会得到一个加性模型。

9.4.1 垃圾邮件例子(续)

下面在本章之前的“垃圾邮件”数据上应用 MARS 方法。为加强可解释性,将 MARS 限制在二阶交叉项以内。虽然目标变量是二分类变量,但仍使用平方误差损失函数(见第 9.4.3 节)。

**图 9.12**:垃圾邮件数据:MARS 过程的测试误差误分类率,横轴为模型的秩(独立的基函数个数)。
图 9.12:垃圾邮件数据:MARS 过程的测试误差误分类率,横轴为模型的秩(独立的基函数个数)。

图 9.12 展示了测试误差误分类率,横轴为模型的秩(独立的基函数个数)。误差率稳定在大约 5.5%,比之前介绍的广义线性模型稍微高一点(5.3%)。根据 GCV 选择的模型大小为 60,大约为可得到最优表现的最小模型。MARS 过程得出的主要的交叉项由输入变量对 $(\text{ch\$},\text{remove})$、$(\text{ch\$},\text{free})$ 和 $(\text{hp},\text{CAPTOT})$。不过这些交叉项并不会带来比广义加性模型效果上的提升。

9.4.2 示例(模拟数据)

下面在三个不同的场景中比较 MARS 的表现。场景中有 $N=100$ 个样本,自变量 $X_1$、$X_2$、……、$X_p$ 和误差 $\varepsilon$ 都服从独立的标准正态分布。

  • 场景 1:数据生成模型为 $$Y = (X_1 - 1)_+ + (X_1 - 1)_+ \cdot (X_2 - .8)_+ + 0.12 \cdot \varepsilon \tag{9.21}$$ 噪声项的标准差设为 0.12,使得信号噪声比大概为 5。将这称为张量积场景;乘积项的曲面类似于图 9.11 中所示。
  • 场景 2:与场景 1 一样,但共有 $p=20$ 个自变量;即有 18 个与输出变量独立的输入变量。
  • 场景 3:这是一个神经网络的结构: $$\begin{align} \ell_1 &= X_1 + X_2 + X_3 + X_4 + X_5 \\ \ell_2 &= X_6 - X_7 + X_8 - X_9 + X_{10} \\ \sigma(t) &= 1 / (1 + e^{-t}) \\ Y &= \sigma(\ell_1) + \sigma(\ell_2) + 0.12 \cdot \varepsilon \end{align}\tag{9.22}$$

场景 1 和 2 最适合使用 MARS,而场景 3 包含了高阶交叉项,可能难以通过 MARS 来近似。每个模型进行了五次模拟,并记录了结果。

场景 1 中,MARS 通常可以几乎完美地得出正确的模型。场景 2 中,MARS 在正确的模型中又添加了若干涉及其他自变量的外部项。

令 $\mu(x)$ 为 $Y$ 的真实均值,令

$$\begin{align} \operatorname{MSE}_0 &= \operatorname{ave}_{x \in \operatorname{Test}} (\bar{y} - \mu(x))^2 \\ \operatorname{MSE} &= \operatorname{ave}_{x \in \operatorname{Test}} (\hat{f}(x) - \mu(x))^2 \end{align}\tag{9.23}$$

这两者分别为常数模型和 MARS 模型拟合在 $x$ 的 1000 个测试集样本上估计出的均方误差。表 9.4 列出了每个场景的模型误差下降比例,或 $R^2$:

$$R^2 = \frac{\text{MSE}_0 - \text{MSE}}{\text{MSE}_0} \tag{9.24}$$
场景 均值 (标准误差)
1:张量积 $p=2$ 0.97 (0.01)
2:张量积 $p=20$ 0.96 (0.01)
3:神经网络 0.79 (0.01)

表9.4:MARS 应用在三个不同场景中的模型误差下降比例($R^2$)。

表中的值为五次模拟的均值和标准误差。场景 2 引入无效的变量后,MARS 的表现只有稍微的下降;但在场景 3 中它的表现明显地下降。

9.4.3 其他问题

分类问题中的 MARS

MARS 方法经推广后可以处理分类问题。有如下的若干拓展方式。

对于二分类问题,可将输出变量编码为 0/1,并如同回归问题一样训练模型;在“垃圾邮件”例子中就是如此操作的。如果有多于两个类别,可以用第 4.2 节中介绍的指示输出变量的方法。将 $K$ 个输出类别通过 0/1 的指示变量编码,然后进行一个多输出变量的 MARS 回归。在这个回归中对所有输出变量使用的是共同的一组基函数。最大的输出变量预测值的类别就是模型给出的分类结果。然而如同第 4.2 节中所述,这个方法可能会存在屏蔽(masking)的问题。第 12.5 节会介绍的“最优得分”(optimal scoring)方法通常是一个更优秀的方法。

Stone et al. (1997) 设计了专门处理分类问题的 MARS 的混合模型,称为“多项式 MARS”(PolyMARS)1。这个方法使用了第 4.4 节中介绍的多个对数几率回归的结构。它如 MARS 一样以前向分段的方式训练模型,但在每一步中使用多项分布对数似然函数的二次近似来寻找下一个基函数对。找到之后,用最大似然方法拟合添加了新基函数的模型,并重复这个过程。

MARS 与 CART 的关系

尽管看似非常不同的方法,但 MARS 和 CART 实际上非常相似。假设将 MARS 的过程进行如下的改动:

  • 用分步函数 $I(x-t>1)$ 和 $I(x-t\leq 0)$ 替换分段线性基函数。
  • 当候选的乘积项中包含了模型中的某个原有项,则用该交叉项替换原有项,因此在接下来的交叉项中不再会出现原有项。

这样的改动后,MARS 的前向过程与 CART 的树生成算法是一样的。将一个分步函数与一个分步函数镜像对的乘积,就等价于在这一步骤对节点的分割。第二个约束使得一个节点不可以被多次分割,因而得出 CART 模型所希望的二叉树结构。另一方面,这个约束也正是 CART 难以对加性结构建模的原因。MARS 舍弃了树的结构从而获得了捕捉加性效应的能力。

混合自变量

如 CART 一样,MARS 可自然地处理混合自变量(既有量化也有分类变量)。MARS 会考察所有将分类变量分为两组的可能分割方式。每个如此的分割都会产生一对分段常数基函数,即两组类别集合的指示函数。接下来可将这个基函数对与其他函数对一样对待,并可用于与模型中已有的其他基函数生成张量积。


本节练习

练习 9.6

使用图 6.9 中的臭氧浓度数据。

  1. 利用温度、风速、和辐射水平,对臭氧浓度的三次方根拟合一个加性模型。将拟合结果与图 6.9 中的格架(trellis)图结果进行比较。
  2. 在同样的数据上拟合树模型、MARS、和 PRIM,并于图 6.9 和加性模型的结果进行比较。

  1. R 程序包的文档中的介绍是“multivariate adaptive polynomial splines regression”。 ↩︎

上一页
下一页