本节中介绍一个通过对输出变量衍生变量进行线性回归从而实现 LDA 的方法。这个方法又可以将 LDA 扩展为非参数以及更灵活的模型。与第四章一样地，假设观测样本为分类输出变量 $G$，取值范围是 $K$ 个类别的集合 $\mathcal{G}=\{1,\dots,K\}$，输入特征变量为 $X\in\mathbb{R}^p$。令 $\theta:\mathcal{G}\mapsto\mathbb{R}^1$ 为一个给类别“评分”（score）的函数，然后用 $X$ 的线性回归来预测这个类别标签变换后的输出变量：若将训练集写为 $(g_i,x_i),i=1,2,\dots,N$，则需要求解
$$\min_{\beta,\theta} \sum_{i=1}^N (\theta(g_i) - x_i^T\beta)^2 \tag{12.52}$$为避免无效解，需要对 $\theta$ 加上某些约束条件（例如在样本集上均值为零方差为一）。这样就得到了类别之间在一维空间上的分离。
更一般地，我们可以为类别标签寻找 $L\leq K-1$ 组不同（或独立）的评分函数，$\theta_1,\theta_2,\dots,\theta_L$，以及 $L$ 个相应地在 $\mathbb{R}$ 上的回归解出的线性函数 $\eta_\ell(X)=X^T\beta_\ell, \ell=1,\dots,L$。选择合适的评分函数 $\theta_\ell(g)$ 和线性函数 $\beta_\ell$，使得平均残差平方（average squared residual）最小：
$$\text{ASR} = \frac{1}{N}\sum_{\ell=1}^L \left [ \sum_{i=1}^N (\theta_\ell(g_i) - x_i^T \beta_\ell)^2 \right ] \tag{12.53}$$这些评分函数需要满足相互正交（mutually orthogonal）并且按一个恰当的内积标准化，以避免无效的零解。
这些计算过程的目的是什么？可以证明第 4.3.3 节中推导的判别向量或典型（canonical）向量 $\nu_\ell$ 即为 $\beta_\ell$ 加上一个常数（Mardia et al.，1979；Hastie et al.，1995）。并且，一个（测试集或新的）点与类别 $k$ 的中心点 $\hat{\mu}_k$ 的马氏距离（Mahalanobis distance）为
$$\delta_J(x, \hat{\mu}_k) = \sum_{\ell=1}^{K-1} w_\ell(\hat{\eta}_\ell(x) - \bar{\eta}_\ell^k)^2 + D(x) \tag{12.54}$$其中 $\bar{\eta}_\ell^k$ 为类别 $k$ 中 $\hat{\eta}_\ell(x_i)$ 的均值，$D(x)$ 不依赖于 $k$。这里的坐标权重 $w_\ell$ 的定义依赖于对第 $\ell$ 个评分函数最优拟合的平均平方残差 $r_\ell^2$。
$$w_\ell = \frac{1}{r_\ell^2(1 - r_\ell^2)} \tag{12.55}$$第 4.3.2 节中说明了在高斯分布假设下的分类问题中，当每个类别的协方差矩阵相同时，模型只需要这些典型（canonical）距离即可。结论可总结为：
LDA 可以通过一系列线性回归来实现，分类规则为在拟合结果的空间中与样本点距离最近的中心点对应的类别。这种实现方式在降秩（$L<K-1$）和满秩（$L=K-1$）的情况下都适用。
这个结果真正的有用之处是可以基于此对模型推广。我们可以用更加灵活的以及非参数的拟合函数替代线性回归拟合 $\eta_\ell(x)=x^T\beta_\ell$，而照此类推可得到一个比 LDA 更灵活的分类器。我们考虑中的替代拟合函数有一般加性模型、样条函数、MARS 模型和其他类似的模型。在更一般的设定下，回归问题的准则函数可定义为
$$\text{ASR}(\{\theta_\ell, \eta_\ell\}_{\ell=1}^L) = \frac{1}{N} \sum_{\ell=1}^L \left[ \sum_{i=1}^N (\theta_\ell(g_i) - \eta_\ell(x_i))^2 + \lambda J(\eta_\ell) \right]$$ $$\tag{12.56}$$其中的 $J$ 是那些非参数回归对应的恰当的正则项，例如平滑样条、加性样条、低阶 ANOVA 样条模型等等。并且也可能是由核函数（如同第 12.3.3）生成的一系列函数以及相应的惩罚项。
在进入这个推广方法的计算细节之前，我们先看一个很简单的例子。假设我们对每个 $\eta_\ell$ 选用的是二阶多项式回归。则每个拟合函数是二次的，并且在比较距离的时候它们的平方项会和 LDA 中一样被消掉，因此 12.54 所隐含的决策边界是以后二次曲面。我们也可以用更传统的方法得出一致的二次边界，即将原始预测变量的平方和交叉项添加到特征变量中。在扩大的特征空间中应用 LDA，扩大空间上的线性边界映射到原空间上就会是一个二次边界。一个经典的例子是在中心在原点的两个多元高斯分布，一个协方差矩阵为 $I$，另一个的为 $cI， c>1$，如图 12.9 所示。贝叶斯决策边界是球形 $\|x\|=\frac{pc\log(c)}{2(c-1)}$，这在扩大特征空间上是一个线性边界。
很多非参数回归方法的流程是先生成一组衍生变量作为基扩展，然后在扩大的特征空间上进行线性回归。MARS 方法（第 9 章）就是一个典型。平滑样条和加性样条模型会生成非常大的基扩展（加性样条有 $N\times p$ 个基函数），但之后在扩大特征空间上进行了带惩罚项的回归拟合。SVM 原理也是如此，参考第 12.3.7 节中基函数生成回归的例子。而本节的 FDA 可以看成是在扩大特征空间上进行的带惩罚项的线性判别分析。第 12.6 节 会更详细说明。扩大特征空间中的线性边界退化映射到原空间上成为非线性边界。这与 SVM 所使用的范式也是完全相同的（第 12.3 节）。
我们用语音识别的实例来演示 FDA，其中有 $K=11$ 个类别和 $p=10$ 个预测变量。 数据的类别对应了 11 个元音，分别出现在 11 个不同的词汇中。下表列出了这些词汇以及代表其中元音的音标符号。
| Vowel | Word | Vowel | Word | Vowel | Word | Vowel | Word |
|---|---|---|---|---|---|---|---|
| i: | heed | O | hod | I | hid | C: | hoard |
| E | head | U | hood | A | had | u: | who’d |
| a: | hard | 3: | heard | Y | hud |
在训练集中，八位发音者每个人把每个单词说六次；测试集与之相似，但有七位发音者。十个预测变量是从语音数字化后通过一个语音识别领域标准而复杂的方法衍生出的。共有 528 个训练观测样本和 462 各测试观测样本。图 12.10 展示了LDA 和 FDA 的结果在二维平面的投影。FDA 模型使用了自适应加性样条回归函数来拟合 $\eta_\ell(x)$，右图中的点的坐标为 $\hat{\eta}_1(x_i)$ 和 $\hat{\eta}_2(x_i)$。在 S-PLUS 中所使用的程序是 bruto，标记在了图和表格 12.3 的标题中。可见在这个例子中（对 $\eta$ 的）灵活的拟合函数在分离类别上有一定作用。表 12.3 展示了一些分类方法的训练和测试误差率。FDA/MARS 指的是 Friedman 的多元自适应回归样条；degree = 2 意思是允许两两乘积。注意到 FDA/MARS 模型的最好分类结果出现在一个降秩子空间上。
| Error Rates | Error Rates | ||
|---|---|---|---|
| Technique | Training | Test | |
| (1) | LDA | 0.32 | 0.56 |
| Softmax | 0.48 | 0.67 | |
| (2) | QDA | 0.01 | 0.53 |
| (3) | CART | 0.05 | 0.56 |
| (4) | CART (linear combination splits) | 0.05 | 0.54 |
| (5) | Single-layer perceptron | 0.67 | |
| (6) | Multi-layer perceptron (88 hidden units) | 0.49 | |
| (7) | Gaussian node network (528 hidden units) | 0.45 | |
| (8) | Nearest neighbor | 0.44 | |
| (9) | FDA/BRUTO | 0.06 | 0.44 |
| Softmax | 0.11 | 0.50 | |
| (10) | FDA/MARS (degree = 1) | 0.09 | 0.45 |
| Best reduced dimension (=2) | 0.18 | 0.42 | |
| Softmax | 0.14 | 0.48 | |
| (11) | FDA/MARS (degree = 2) | 0.02 | 0.42 |
| Best reduced dimension (=6) | 0.13 | 0.39 | |
| Softmax | 0.10 | 0.50 |
表 12.3：元音识别数据的表现结果。神经网络的文献结果是在更大范围的比较中，神经网络的效果最好。FDA/BRUTO 标记是指 FDA 中所使用的回归方法。
12.5.1 FDA 估计的计算
在很多重要的场景中，尤其是当非参数回归方法可以用一个线性运算表述时，FDA 坐标的计算可以被简化。我们用 $\mathbf{S}_\lambda$ 来代表这个线性运算，即 $\hat{\mathbf{y}}=\mathbf{S}_\lambd\mathbf{y}$，其中的 $\mathbf{y}$ 为输出变量的向量，而 $\hat{\mathbf{y}}$ 是拟合值的向量。当平滑参数固定是，加性样条可满足这个条件；当选定了基函数时，MARS 也满足这个条件。下标 $\lambda$ 代表了所有的平滑参数。这里的最优评分函数就等价于一个典型关联（canonical correlation）问题，它的求解只需要进行一次特征分解。练习 12.6 会详细说明这个过程，这里只介绍最终的算法。
我们从输出变量 $g_i$ 构建一个 $N\times K$ 的指示变量矩阵 $\mathbf{Y}$，使得当 $g_i=k$ 时 $yh_{ik}=1$，否则 $y_{ik}=0$。一个五个类别的 $\mathbf{Y}$ 的结构如下：
$$\begin{matrix} & \begin{matrix}C_1&C_2&C_3&C_4&C_5\end{matrix} \\ \begin{matrix}g_1=2\\g_2=1\\g_3=1\\g_4=5\\g_5=4\\\vdots\\g_N=3\end{matrix} & \begin{pmatrix} 0 & 1 & 0 & 0 & 0 \\ 1 & 0 & 0 & 0 & 0 \\ 1 & 0 & 0 & 0 & 0 \\ 0 & 0 & 0 & 0 & 1 \\ 0 & 0 & 0 & 1 & 0 \\ \vdots & \vdots & \vdots & \vdots & \vdots \\ 0 & 0 & 1 & 0 & 0 \\ \end{pmatrix} \end{matrix}$$计算步骤如下：
- 多元非参数回归。
用 $\mathbf{X}$ 对 $\mathbf{Y}$ 进行一个多输出变量的、自适应非参数回归，得到拟合结果 $\hat{\mathbf{Y}}$。令 $\mathbf{S}_\lambda$ 为最终选定模型的拟合结果对应的线性运算，并令 $\eta^*(x)$ 为回归函数拟合值的向量。 - 最优评分函数。
计算 $\mathbf{Y}^T\hat{\mathbf{Y}}=\mathbf{Y}^T\mathbf{S}_\lambda\mathbf{Y}$ 的特征分解，对特征向量 $\mathbf{\Theta}$ 进行标准化 $\mathbf{\Theta}^T\mathbf{D}_\pi\mathbf{\Theta}=\mathbf{I}$。其中的 $\mathbf{D}_\pi=\mathbf{Y}^T\mathbf{Y}/N$ 是类别先验概率估计的对角线矩阵。 - 用最优评分函数 $\eta(x)=\mathbf{\Theta}^T\eta^*(x)$ 更新步骤（1）中的模型。
$K$ 个 $\eta(x)$ 函数中的第一个是常数函数，一个零解（trivial solution），其他 $K-1$ 个函数是判别函数。这个常数函数，以及进行的标准化，使所有其他函数中心化（centered）。
再次说明一下，$\mathbf{S}_\lambda$ 可以代表任意回归方法。如果 $\mathbf{S}_\lambda=\mathbf{H}_X$，即线性回归投影算子，则 FDA 就变成了线性判别分析。在关于计算量 中提到的软件就利用了这个模块性质，fda 函数有一个 method= 参数可以让用户提供任意的回归函数，只需要满足一些正常的规范。我们提供的回归函数有多项式回归、自适应加性模型、和 MARS。这些函数都可以有效地处理多个输出变量，因此步骤（1）只需要调用一次回归函数。步骤（2）中的特征分解会同时计算出所有的最优评分函数。
第 4.2 节提到了在分类问题中使用对指示输出变量矩阵的线性回归时存在的问题。尤其当有三个或更多类别时可能会出现屏蔽（masking）的问题。FDA 先从步骤 1 中的回归得到拟合值，但然后将他们进一步变换得到不再受这些问题影响的有效的判别函数。练习 12.9 从另一个角度解释了这个现象。
本节练习
练习 12.6
Suppose that the regression procedure used in FDA (Section 12.5.1) is a linear expansion of basis functions $h_m(x)$ $m=1,\dots,M$. Let $\mathbf{D}_pi=\mathbf{Y}T\mathbf{Y}/N$ be the diagonal matrix of class proportions.
- Show that the optimal scoring problem (12.52) can be written in vector notation as $$\min_{\theta,\beta} \|\mathbf{Y}\theta - \mathbf{H}\beta\|^2 \tag{12.65}$$ where $\theta$ is a vector of $K$ real numbers, and $\mathbf{H}$ is the $N\times M$ matrix of evaluations $h_j(x_i)$.
- Suppose that the normalization on $\theta$ is $\theta^T\mathbf{D}_\pi\mathbf{1}=0$ and $\theta^T\mathbf{D}_\pi\theta=1$. Interpret these normalizations in terms of the original scored $\theta(g_i)$.
- Show that, with this normalization, (12.65) can be partially optimized w.r.t. $\beta$, and leads to $$\max_\theta \theta^T \mathbf{Y}^T\mathbf{S}\mathbf{Y} \theta \tag{12.66}$$ subject to the normalization constraints, where $\mathbf{S}$ is the projection operator corresponding to the basis matrix $\mathbf{H}$.
- Suppose that the $h_j$ include the constant function. Show that the largest eigenvalue of $\mathbf{S}$ is 1.
- Let $\mathbf{\Theta}$ be a $K\times K$ matrix of scores (in columns), and suppose the normalization is $\mathbf{\Theta}^T\mathbf{D}_pi\mathbf{\Theta}=\mathbf{I}$. Show that the solution to (12.53) is given by the complete set of eigenvectors of $\mathbf{S}$; the first eigenvector is trivial, and takes care of the centering of the scores. The remainder characterize the optimal scoring solution.
练习 12.9
Let $\hat{\mathbf{Y}}=\mathbf{X}\hat{\mathbf{B}}$ be the fitted $N\times K$ indicator response matrix after linear regression on the $N\times p$ matrix $\mathbf{X}$, where $p>K$. Consider the reduced features $x_i^*=\hat{\mathbf{B}}^Tx_i$. Show that LDA using $x_i^*$ is equivalent to LDA in the original space.