Ramaswamy et al. (2001) 提出了一个更困难的微阵列分类问题,其中有 14 个不同的癌症类别、144 个病人的训练集、和 54 个病人的测试集。共有 16,063 个基因的表达水平测量值。
方法 | CV erros (SE) Out of 144 | Test errors Out of 54 | Number of Genes Used |
---|---|---|---|
1. 最近收缩中心 | 35 (5.0) | 17 | 6,520 |
2. $L_2$ 惩罚的判别分析 | 25 (4.1) | 12 | 16,063 |
3. 支持向量分类器 | 26 (4.2) | 14 | 16,063 |
4. 套索回归(OVA) | 30.7 (1.8) | 12.5 | 1,429 |
5. k-最近邻 | 41 (4.6) | 26 | 16,063 |
6. $L_2$ 惩罚的多项分布 | 26 (4.2) | 15 | 16,063 |
7. $L_1$ 惩罚的多项分布 | 17 (2.8) | 13 | 269 |
8. 弹性网惩罚的多项分布 | 22 (3.7) | 11 | 384.8 |
表 18.1: 14 个癌症类别在微阵列数据上的预测结果。第 18.2 节介绍了方法 1。第 18.3 节介绍了方法 2、3、和 6。第 18.4 节介绍了方法 4、7、和 8。第 13.3 节介绍了方法 5。弹性网惩罚的多项分布在测试数据上表现最好,不过每个测试误差估计的标准差都大约为 3,所以这个对比也无定论。
表 18.1 展示了八个不同分类方法的预测结果。首先将患者的数据进行均值为 0 和方差为 1 的标准化。在这个例子中,数据标准化貌似改进了预测的准确率,这说明了重要的是每个基因表达分布的“形状”,而不是表达水平的绝对值。在每个模型中,通过对交叉验证误差最小化来选择正则参数,表中展示了所选参数对应的测试误差。如果交叉验证误差最小化的正则参数有多个取值,则展示这些取值对应的测试误差的平均数。
正则判别分析(RDA)、正则多项对数几率回归、以及支持向量机都是从多元数据集中提取信息的更复杂的方法。我们依次介绍每个方法,以及一些正则化的方法,包括了 $L_1$、$L_2$、以及介于之间的正则化方法。
18.3.1 正则判别分析
第 4.3.1 节中介绍了正则判别分析(RDA)。线性判别分析中涉及到对一个 $p\times p$ 的组内协方差矩阵求逆运算。当 $p\gg N$ 时,这个矩阵可能很大,它的秩最多为 $N<p$,所以它是奇异矩阵。RDA 通过对组内协方差的估计 $\hat{\mathbf{\Sigma}}$ 进行正则化,来克服奇异性的问题。这里我们使用了 RDA 的一个版本,它将 $\hat{\mathbf{\Sigma}}$ 向它的对角矩阵收缩:
$$\hat{\mathbf{\Sigma}}(\gamma) = \gamma \hat{\mathbf{\Sigma}} + (1 - \gamma) \operatorname{diag}(\hat{\mathbf{\Sigma}}), \text{ with } \gamma \in [0,1] \tag{18.9}$$注意到 $\gamma=0$ 就对应了对角 LDA,也就是“无收缩”版本的最小收缩中心。式 18.9 中的收缩形式与岭回归(第 3.4.1 节)非常相似,后者是将特征的总体协方差矩阵向一个对角(常数)矩阵收缩。实际上,如果将线性判别分析看作是带有分类输出变量的最优得分函数(optimal scoring)的线性回归,就能更清楚地看出他们的等价性。
第 18.3.5 节将介绍一些克服对 $p\times p$ 大矩阵求逆的计算负担的方法。在表 18.1 的第二行中使用了交叉验证来选择 $\gamma$ 的值;当 $\gamma\in(0.002,0.550)$ 区间上时,可以得到同样的 CV 和测试误差。Guo et al. (2006) 介绍了 RDA 进一步的进展,包括除了对协方差矩阵外增加了对中心点的收缩。
18.3.2 二次正则的对数几率回归
对数几率回归(第 4.4 节)可以通过相似的方式进行改造,来应对 $p\gg N$ 的情况。假设有 $K$ 个类别,我们使用一个对称版本的多类别对数几率模型(119 页的式 4.17):
$$\operatorname{Pr}(G=k | X=x) = \frac{\exp(\beta_{k0} + x^T \beta_k)} {\sum_{\ell=1}^K \exp( \beta_{\ell 0} + x^T \beta_\ell)} \tag{18.10}$$这里有 $K$ 个对数几率参数的系数向量 $\beta_1$、$\beta_2$、……、$\beta_K$。我们通过带惩罚项的对数似然函数最大化来对拟合正则化:
$$\max_{\{\beta_{0k}, \beta_k\}_1^K} \left[ \sum_{i=1}^N \log \operatorname{Pr}(g_i | x_i) - \frac{\lambda}{2} \sum_{k=1}^K \|\beta_k\|_2^2 \right] \tag{18.11}$$这个正则化会自动地解决参数的冗余,并且会控制了 $\sum_{k=1}^K\hat{\beta}_{kj}=0$,$j=1,\dots,p$(练习 18.3)。需要注意常数项 $\beta_{k0}$ 没有被正则化(其中之一应被设置为零)。这样得出的最优化问题是一个凸问题,可以用牛顿算法或其他的数值方法求解。更多细节可参考 Zhu and Hastie (2004)。Friedman et al. (2010) 给出了计算两个和多个类别对数几率模型的正则化路径的软件。表 18.1 的第 6 行中的“多项分布”就代表了多类别的对数几率回归模型。可证明(Rosset et al., 2004a)当数据可分(separable)时,随着 $\lambda\to0$,正则化(两类别)对数几率回归估计(重新标准化后)收敛到最大间隔(maximal margin)分类器(第 12.2 节)。这是对下一段介绍的支持向量机的一个很有吸引力的替代方法,尤其是在多类别的情况下。
18.3.3 支持向量分类器
第 12.2 节介绍了两个类别情况下的支持向量分类器。当 $p>N$ 时,它非常有优势,因为除非在不同类别中有完全一样的特征向量,类别一般都是完全可分的。在不加任何正则条件时,支持向量分类器得出可使得间隔(margin)最大的分离超平面,也就是在训练数据集上两个类别之间的间隙最大的超平面。有点出乎意料的是当 $p\gg N$ 时,无正则条件的支持向量分类器通常与最佳正则条件版本的表现几乎一样好。过拟合(在这个情景中)通常不会出现,误分类损失的不敏感性也是一部分原因。
将两个类别的支持向量分类器推广到 $K>2$ 个类别,有很多不同的方法。在“一对一(one versus one,OVO)”方法中,我们要计算所有的 $\binom{K}{2}$ 个两两分类器。对每一个测试点,它的预测类别就是在所有的两两比较中胜出次数最多的类别。在“一对所有(one versus all,OVA)”方法中,在两个类别的 $K$ 个分类器中,将每个类别与其他所有类别进行比较。要对一个测试点进行分类,我们要计算 $K$ 个分类器每一个的置信度(距离超平面的有向距离)。置信度最高的类别就是这个点的预测类别。最后,Vapnik (1998) 以及 Weston and Watkins (1999) 提出了推广两个类别准则(式 12.7)的一个有些复杂的多类别准则。
Tibshirani and Hastie (2007) 提出了 间隔树(margin tree) 分类器,在一个二叉树模型中使用支持向量分类器,它与 CART(第九章)很相似。类别是以分层的方式组织起来的,例如在本节中它可以用来将患者分类到不同的癌症类型中。
表 18.1 的第三行展示了使用 OVA 方法的支持向量分类器的结果;Ramaswamy et al. (2001) 调查得出这个方法对这个问题是最有效的(我们也证实了)。误分类个数与第六行的非常相似,这从上一节最后一段的论述中可能猜想得到。当 $C>0.001$ 时,误差率对 $C$ 的选择(420 页式 12.8 中的正则参数)不敏感。既然 $p>N$,支持向量超平面可以令 $C=\infty$ 来完全地分离训练数据集。
18.3.4 特征选择
当 $p$ 比较大时,特征选择是分类器的一个重要的科学理解上的要求。由于使用了二次正则条件,判别分析、对数几率回归、和支持向量分类器都不会自动进行特征选择。在这些模型中,所有的特征都有非零的权重。有一些设计出的特征选择方法,例如去除掉系数小的基因,然后再次拟合分类器。这是以一个后向分布的方式进行的,先从(绝对值)最小的权重开始,然后再到更大的权重。这也被称为 递归特征消除(recursive feature elimination,RFE)(Guyon et al., 2002)。但在这个例子中这个方法效果不好;例如 Ramaswamy et al. (2001) 的研究发现,当(使用的)基因数量从全部的 16,603 个开始缩减时,支持向量分类器的准确度就开始下降了。这是非常奇怪的现象,因为训练样本只有 144 个。对这个现象我们也尚且无法解释。
本节中的三个方法(RDA、LR、和 SVM)都可以使用核函数来拟合非线性的决策边界。这样的做法一般是为了要增加模型的复杂度。但是当 $p\gg N$ 时,模型已经足够复杂了,而且总会有过拟合的危险。不过尽管维度很高,径向核函数(第 12.3.3 节)在这些高维度问题中有时也能得到更好的结果。径向基往往会削弱彼此距离较远的点之间的内积,从而使其对离群值比较稳健。这在高维度问题中会经常出现,所以可能是它表现良好的原因。我们在表 18.1 中的支持向量机中尝试了径向基,但在这个例子中使用后的表现变差了。
18.3.5 $p\gg N$ 时的计算捷径
本节介绍的计算方法可用于任意的对系数进行二次正则的线性模型的拟合中。这包含了本节上文中介绍的模型,以及其他很多模型。当 $p>N$ 时,是在一个 $N$ 维空间上通过第 14.5 节介绍的奇异值分解来进行计算,而不是 $p$ 维空间。从代数几何上可以理解为:就如同三维空间上的两个点总会形成一条直线,$p$ 维空间上的 $N$ 个点会形成一个 $(N-1)$ 维的仿射子空间。
给定 $N\times p$ 的矩阵 $\mathbf{X}$:
$$\begin{align} \mathbf{X} &= \mathbf{U} \mathbf{D} \mathbf{V}^T \tag{18.12}\\ &= \mathbf{R} \mathbf{V}^T \tag{18.13} \end{align}$$令以上为 $\mathbf{X}$ 的奇异值分解(SVD);即 $\mathbf{V}$ 是 $p\times N$ 的列标准正交的矩阵,$\mathbf{U}$ 是 $N\times N$ 的正交矩阵,$\mathbf{D}$ 是一个对角矩阵,它的对角元素 $d_1\geq d_2\geq\dots\geq d_N\geq0$。$\mathbf{R}$ 是 $N\times N$ 的矩阵,它的行为 $r_i^T$。
首先,用一个岭回归的估计作为一个简单的例子:
$$\hat{\beta} = (\mathbf{X}^T \mathbf{X} + \lambda \mathbf{I})^{-1} \mathbf{X}^T \mathbf{y} \tag{18.14}$$用 $\mathbf{R}\mathbf{V}^T$ 代替 $\mathbf{X}$,经过一些变换操作后,可改写为(练习 18.4):
$$\hat{\beta} = \mathbf{V} (\mathbf{R}^T \mathbf{R} + \lambda \mathbf{I})^{-1} \mathbf{R}^T \mathbf{y} \tag{18.15}$$因此 $\hat{\beta}=\mathbf{V}\hat{\theta}$,其中的 $\hat{\theta}$ 是使用 $N$ 个观测样本 $(r_i,y_i)$,$i=1,2,\dots,N$ 的岭回归估计。也就是说,我们可以简单地将数据矩阵从 $\mathbf{X}$ 缩减到 $\mathbf{R}$,然后基于 $\mathbf{R}$ 的行进行计算。当 $p>N$ 时,这个技巧将计算量从 $O(p^3)$ 降低至 $O(pN^2)$。
这个结果可以推广到所有的对参数呈线性并且有二次惩罚项的模型。在任意的监督学习问题中,我们用线性函数 $f(X)=\beta_0+X^T\beta$ 来建立条件概率 $Y|X$ 的参数模型。我们通过对某个损失函数 $\sum_{i=1}^NL(y_i,f(x_i))$ 以及一个对 $\beta$ 的二次惩罚项在数据集上的最小化,来拟合参数 $\beta$。对数几率回归就是一个具体的例子。则我们可得出以下这个简单的定理:
令 $f^*(r_i)=\theta_0+r_i^T\theta$,其中的 $r_i$ 如式 18.13 所定义,并有以下两个最优化问题:
$$\begin{align} (\hat{\beta}_0, \hat{\beta}) &= \underset{\beta_0,\beta \in \mathbb{R}^p}{\arg\min} \sum_{i=1}^N L(y_i, \beta_0 + x_i^T \beta) + \lambda \beta^T \beta \tag{18.16}\\ (\hat{\theta}_0, \hat{\theta}) &= \underset{\theta_0,\theta \in \mathbb{R}^N}{\arg\min} \sum_{i=1}^N L(y_i, \theta_0 + r_i^T \theta) + \lambda \theta^T \theta\tag{18.17} \end{align}$$则 $\hat{\beta}_0=\hat{\theta}_0$,并且 $\hat{\beta}=\mathbf{V}\hat{\theta}$。
这个定理说明可以简单地用 $N$ 维向量 $r_i$ 代替 $p$ 维向量 $x_i$,和之前一样进行带惩罚项的拟合,但自变量的数量要少得多。然后可以通过一个简单的矩阵乘法将 $N$ 维向量的解 $\hat{\theta}$ 转化回 $p$ 维向量的解。这个结论在统计学的坊间传播(statistics folklore),它应该更广泛地被了解。更多细节可参考 Hastie and Tibshirani (2004)。
从代数几何的角度看,我们将特征变量旋转到一个坐标系中,其中除前 $N$ 个坐标外其他坐标都是零。由于二次惩罚项对旋转是不变的,而且线性模型是等变的(equivariant),所以这种旋转变换是可以操作的。
这个结果可以被应用在本章介绍的很多学习方法中,例如正则化(多类别)对数几率回归、线性判别分析(练习 18.6)、支持向量机。它也可以被应用在二次正则化的神经网络中(第 11.5.2 节)。不过需要注意,它不可应用在使用对系数的非二次惩罚项(例如 $L_1$)的方法中,比如套索方法。
我们通常使用交叉验证来选择参数 $\lambda$。可证明(练习 18.12)我们只需要在原始数据集上构建一次 $\mathbf{R}$,然后在每一次的交叉验证中都可以将它用作为数据集。
第 12.3.7 节中支持向量的“核函数技巧”利用了和本节一样的化简,只是场景稍微不同。假设我们已经得到了 $N\times N$ 的格拉姆矩阵(gram matrix,内积矩阵)$\mathbf{K}=\mathbf{X}\mathbf{X}^T$。从式 18.12 可得 $\mathbf{K}=\mathbf{U}\mathbf{D}^2\mathbf{U}^T$,所以 $\mathbf{K}$ 与 $\mathbf{R}$ 包含着同样的信息。练习 18.13 会展示如何利用本节中的思路来用 $\mathbf{K}$ 和它的奇异值分解(SVD)来拟合一个岭式对数几率回归(ridged logistic regression)。
本节练习
练习 18.3
证明正则化的多类别对数几率回归问题(式 18.10)的拟合参数满足:
$$ \sum_{k=1}^K \hat{\beta}_{kj} = 0, j = 1, \dots, p $$$\hat{\beta}_{k0}$ 是否会满足?考虑这些常数参数的问题,以及如何解决这些问题。
练习 18.4
推导岭回归的计算公式(式 18.15)。
【提示】:使用带惩罚项的平方和准则的一阶导数,证明如果 $\lambda>0$ 则存在某个 $s\in\mathbb{R}^N$ 使得 $\hat{\beta}=\mathbf{X}^Ts$。
练习 18.5
Prove the theorem (18.16)–(18.17) in Section 18.3.5, by decomposing $\beta$ and the rows of $\mathbf{X}$ into their projections into the column space of $\mathbf{V}$ and its complement in $\mathbb{R}^p$.
练习 18.6
说明第 18.3.5 节中的定理可以如何被应用在正则判别分析中(式 4.14 和 18.9)。
练习 18.12
假设在一个 $p\gg N$ 的(任意的线性模型)情境中,我们想要通过 10-折交叉验证来选择岭回归参数 $\lambda$。我们想使用第 18.3.5 节中介绍的计算捷径。试说明我们只需要进行一次从 $N\times p$ 矩阵 $\mathbf{X}$ 到 $N\times N$ 矩阵 $\mathbf{R}$ 的简化过程,然后在所有的交叉验证中都可以使用这个结果。
练习 18.13
假设 $p>N$ 个自变量形成了一个 $N\times N$ 的内积矩阵 $\mathbf{K}=\mathbf{X}\mathbf{X}^T$,我们想要进行一个与原特征变量上的二次正则化线性对数几率回归模型等价的拟合。(对测试点的)预测也是使用内积来得出的,即一个新的数据点 $x_0$ 会被表述为 $k_0=\mathbf{X}x_0$。令 $\mathbf{K}=\mathbf{U}\mathbf{D}^2\mathbf{U}^T$ 为矩阵 $\mathbf{K}$ 的特征分解。试说明预测结果是 $\hat{f}_0=k_0^T\hat{\alpha}$,其中的
- $\hat{\alpha}=\mathbf{U}\mathbf{D}^{-1}\hat{\beta}$。
- $\hat{\beta}$ 是使用输入变量矩阵 $\mathbf{R}=\mathbf{U}\mathbf{D}$ 的岭式对数几率回归的估计。
试说明对任一恰当的核矩阵 $\mathbf{K}$ 都可以使用同样的方法。