一、充分这两个字是什么意思

正式的说法是：找一个 d 乘 k 的矩阵 B，使得给定 B 转置乘 x 之后，y 与 x 条件独立。 翻译成人话：把 x 压到 k 维之后，关于 y 的信息一点没丢。 这个定义有两个直接后果。 第一，k 是数据自己决定的，不是你先定的。 13 章的 PCA 要你说「降到几维」，这一节可以先问「至少要几维」。 第二，这个降维是「针对这个 y」的。 换一个标签，答案就换一套——这点和判别分析一样，和 PCA 相反。 和 PCA 的关系也要说清楚。 PCA 保证的是「在 k 维里保留最多方差」，它完全不知道 y 是什么。 本节的数据里，方差最大的方向恰好是纯干扰，所以 PCA 一定会先把它捡走。

二、反过来看：切片逆回归

直接找 B 太难，因为不知道 y 依赖 x 的具体形式。 换一个角度：如果 y 只通过 B 转置乘 x 依赖 x， 那么条件均值 E[x 给定 y] 一定落在 B 张成的空间里。 道理很简单——y 含有的关于 x 的信息全在那 k 维里， 给定 y 之后 x 的期望当然也只能在那 k 维里动。 于是算法就出来了，一共四步：

第一步  把 x 标准化
第二步  按 y 的大小把样本切成 H 段，算每段里 x 的平均
第三步  这些平均值张成的空间就是候选答案
第四步  为了去掉尺度影响，在协方差度量下做主成分

第四步为什么要在协方差度量下做？ 因为切片均值本身的离散程度受 x 各列尺度影响，直接做主成分会被尺度大的列主导。 用协方差矩阵的逆去加权， 相当于先把空间「白化」再做——这和 17.1 里解广义特征值问题是同一个手法。 具体算的是 Σ 逆乘 Γ 的特征向量，其中 Σ 是 x 的协方差，Γ 是切片均值的协方差。 再强调一遍：算出方向之后要换回原始坐标，别在标准化坐标里比方向——这一步我在实现里漏过一次， 夹角算出来 86 度，全是假象。

三、结果：PCA 反而更差

10 维、400 个点，标签由两个方向决定，另有一个方差 5 倍的干扰方向。

四、切片数 H 的选择

H 是这一节最需要调的参数，效果很敏感：

五、两个真实的坑

第一个坑：协方差越病态，它越不准。 把观测噪声从 0.60 一路降到 0.05：

六、动手 Demo：充分降维：不用标签定方向，用标签定切片

这段代码造 10 维数据：标签只依赖两个方向，另有一个方差 5 倍的纯干扰方向。 先比较 PCA 和切片逆回归找回的子空间，再用降维后的坐标做近邻回归，最后扫切片数、调噪声、 换标签形式，看它在什么情况下会坏。 重点看三处。第二段两个夹角 79.9 度和 34.7 度；第三段 PCA 两维的 4.463 比不降维的 3.355 还差； 以及第六段对称依赖下那个 81.3 度。

import numpy as np


def sir_data(n=400, d=10, seed=1, big=5.0, noise=0.6):
    """标签只依赖两个方向，另有一个和标签无关的高方差方向在捣乱。"""
    rg = np.random.default_rng(seed)
    B = rg.normal(0.0, 1.0, (d, 2))
    Z = rg.normal(0.0, 1.0, (n, 2))
    X = Z @ B.T + rg.normal(0.0, noise, (n, d))
    w = rg.normal(0.0, 1.0, d)
    w = w / np.linalg.norm(w)
    X = X + big * np.outer(rg.normal(0.0, 1.0, n), w)
    y = Z[:, 0] ** 3 + Z[:, 1] + rg.normal(0.0, 0.1, n)
    return X, y, B


def pca_dir(X, k=2):
    """不用标签：按方差从大到小取方向。"""
    Xc = X - X.mean(axis=0)
    C = Xc.T @ Xc / len(X)
    vals, vecs = np.linalg.eigh(C)
    order = np.argsort(-vals)
    return vecs[:, order][:, :k], vals[order]


def sir(X, y, H=10, k=2, ridge=1e-6):
    """切片逆回归：把 y 切片，算每片的输入均值，再在协方差度量下做主成分。"""
    sd = X.std(axis=0)
    Xs = (X - X.mean(axis=0)) / sd
    n, d = Xs.shape
    Sig = Xs.T @ Xs / n
    qs = np.quantile(y, np.linspace(0.0, 1.0, H + 1))
    M = []
    for h in range(H):
        if h == H - 1:
            m = (y >= qs[h]) & (y <= qs[h + 1])
        else:
            m = (y >= qs[h]) & (y < qs[h + 1])
        M.append(Xs[m].mean(axis=0))
    M = np.array(M)
    B = M.T @ M / H
    A = np.linalg.solve(Sig + ridge * np.eye(d), B)
    vals, vecs = np.linalg.eig(A)
    order = np.argsort(-np.abs(vals))
    V = vecs[:, order].real[:, :k] / sd[:, None]
    V = V / np.linalg.norm(V, axis=0)
    return V, np.abs(vals[order])


def angle_between(U, V):
    """两个子空间的主角度（度），越小越接近。"""
    Q1, _ = np.linalg.qr(U)
    Q2, _ = np.linalg.qr(V)
    s = np.clip(np.linalg.svd(Q1.T @ Q2, compute_uv=False), -1.0, 1.0)
    return float(np.degrees(np.arccos(s.min())))


def knn_rmse(Z, y, k=5):
    """在降维后的坐标上做 k 近邻回归，返回留一均方根误差。"""
    D2 = ((Z[:, None, :] - Z[None, :, :]) ** 2).sum(-1)
    np.fill_diagonal(D2, np.inf)
    idx = np.argsort(D2, axis=1)[:, :k]
    return float(np.sqrt(np.mean((y[idx].mean(axis=1) - y) ** 2)))


X, Y, B = sir_data(400, 10, 1)
VP, lam_p = pca_dir(X, 2)
VS, vals = sir(X, Y, H=10, k=2)

print('[1] 10 维、400 个点，标签只依赖两个方向')
print('  相关系数矩阵的最小特征值 %.2e' % np.linalg.eigvalsh(np.corrcoef(X.T)).min())
print('  方差最大的那个方向是纯干扰，和标签毫无关系')

print('')
print('[2] 谁把真正的低维结构找回来了')
print('  PCA 前两方向与真子空间夹角 %.1f 度' % angle_between(VP, B))
print('  SIR 前两方向与真子空间夹角 %.1f 度' % angle_between(VS, B))
print('  SIR 的十个广义特征值 ' + '  '.join('%.3f' % v for v in vals[:4]) + '  ...')

print('')
print('[3] 降到两维之后做 5 近邻回归')
print('  原始 10 维   %.3f' % knn_rmse(X, Y))
print('  PCA 两维     %.3f' % knn_rmse(X @ VP, Y))
print('  SIR 两维     %.3f' % knn_rmse(X @ VS, Y))

print('')
print('[4] 切片数 H 的影响')
for H in (2, 5, 10, 20, 40):
    V, _ = sir(X, Y, H=H, k=2)
    print('  H=%-3d 夹角 %5.1f 度  回归误差 %.3f' % (H, angle_between(V, B), knn_rmse(X @ V, Y)))

print('')
print('[5] 协方差越病态，SIR 越不准')
for noise in (0.60, 0.30, 0.05):
    Xn, Yn, Bn = sir_data(400, 10, 1, noise=noise)
    Vn, _ = sir(Xn, Yn, H=10, k=2)
    print('  噪声 %.2f  最小特征值 %.2e  夹角 %5.1f 度'
          % (noise, np.linalg.eigvalsh(np.corrcoef(Xn.T)).min(), angle_between(Vn, Bn)))

print('')
print('[6] 它抓不到什么：把标签换成对称依赖')
rg = np.random.default_rng(1)
Bq = rg.normal(0.0, 1.0, (10, 2))
Zq = rg.normal(0.0, 1.0, (400, 2))
Xq = Zq @ Bq.T + rg.normal(0.0, 0.6, (400, 10))
wq = rg.normal(0.0, 1.0, 10)
wq = wq / np.linalg.norm(wq)
Xq = Xq + 5.0 * np.outer(rg.normal(0.0, 1.0, 400), wq)
Yq = Zq[:, 0] ** 2 + Zq[:, 1] + rg.normal(0.0, 0.1, 400)
Vq, vq = sir(Xq, Yq, H=10, k=2)
print('  y 等于 z0 的平方加 z1 时，SIR 前两方向与真子空间夹角 %.1f 度'
      % angle_between(Vq, Bq))
print('  两个广义特征值 %.3f 和 %.3f，一大一小' % (vq[0], vq[1]))

实跑输出：

--------------------------------------------------------------------
充分降维：不用标签定方向，用标签定切片
--------------------------------------------------------------------
[1] 10 维、400 个点，标签只依赖两个方向
  相关系数矩阵的最小特征值 3.74e-02
  方差最大的那个方向是纯干扰，和标签毫无关系

[2] 谁把真正的低维结构找回来了
  PCA 前两方向与真子空间夹角 79.9 度
  SIR 前两方向与真子空间夹角 34.7 度
  SIR 的十个广义特征值 0.770  0.195  0.072  0.057  ...

[3] 降到两维之后做 5 近邻回归
  原始 10 维   3.355
  PCA 两维     4.463
  SIR 两维     2.414

[4] 切片数 H 的影响
  H=2   夹角  85.2 度  回归误差 3.265
  H=5   夹角  71.6 度  回归误差 2.717
  H=10  夹角  34.7 度  回归误差 2.414
  H=20  夹角  37.3 度  回归误差 2.447
  H=40  夹角  40.8 度  回归误差 2.497

[5] 协方差越病态，SIR 越不准
  噪声 0.60  最小特征值 3.74e-02  夹角  34.7 度
  噪声 0.30  最小特征值 9.79e-03  夹角  51.9 度
  噪声 0.05  最小特征值 2.78e-04  夹角  82.4 度

[6] 它抓不到什么：把标签换成对称依赖
  y 等于 z0 的平方加 z1 时，SIR 前两方向与真子空间夹角 81.3 度
  两个广义特征值 0.456 和 0.085，一大一小

第一段只有两行，但把问题的性质定住了：最小特征值 3.74e-02，说明协方差矩阵可逆但不算健康。 后面第五段会看到，这个数降到 2.78e-04 时算法就崩了。 第二段是全节的核心对比：PCA 79.9 度，切片逆回归 34.7 度。 79.9 度意味着 PCA 找到的两维和真子空间几乎正交——它把干扰方向当成了主要结构。 十个广义特征值 0.770、0.195、0.072、0.057 里，前两个明显大于后面， 这就是「结构维数等于 2」的经验判据。 第三段三个误差 3.355、4.463、2.414 要连起来读。 不降维 3.355，PCA 降到两维变成 4.463（亏了），切片逆回归降到两维 2.414（赚了）。 同样是降维，一个亏一个赚，区别只在「降掉的是噪声还是信息」。 第四段 H 的扫描是典型的「中间好两头坏」：85.2、71.6、34.7、37.3、40.8。 H 等于 2 之所以灾难，是因为两段只能表达一个方向，而真实结构是两维。 H 大于 20 之后缓慢变差，是每段样本太少导致的估计噪声。 第五段三行很有教育意义：噪声 0.60 时夹角 34.7 度，噪声 0.05 时反而 82.4 度。 噪声越小越不准，因为协方差矩阵越来越病态，它的逆把切片均值的采样噪声放大了。 代码里那个 ridge 参数就是干这个的——不能省。 第六段是「它抓不到什么」：把标签换成 z0 的平方，夹角从 34.7 度变成 81.3 度。 因为 z0 对称，给定 y 之后它的条件均值恒为零，一阶方法看不见。 知道算法的盲区，比记住它的优点更有用。

30 秒实验

1. 把 sir_data 的 big 从 5.0 改成 0.0，去掉那个干扰方向，重跑第二段。 这时 PCA 也该找到真子空间了，两者夹角都会很小。 这个对照说明 PCA 的失败完全来自「方差最大的方向没用」这一个原因。
2. 把 sir 的 ridge 从 1e-6 改成 1e-1，重跑第五段。看噪声 0.05 那一行的夹角会不会降下来。 如果降了，就证明了「病态」这个诊断是对的。
3. 把 [2] 里的 k=2 改成 k=1，看一维时两边各是什么夹角。 真实结构维数是 2，降到 1 维一定会丢信息，但丢多丢少差别很大。

常见误区

1. 「降维总能提升近邻方法的效果」。 取决于降掉什么。实测 PCA 降到两维，回归误差从 3.355 涨到 4.463； 切片逆回归降到两维才降到 2.414。 降对方向是去噪，降错方向是毁数据。
2. 「切片越多，对条件期望的估计越准，结果越好」。 H 从 10 加到 40，夹角从 34.7 度涨到 40.8 度。 因为每段的样本数从 40 掉到 10，均值估计的噪声变大，而噪声会经过协方差逆被放大。
3. 「数据越干净，这个方法越准」。 正好相反：噪声从 0.60 降到 0.05，夹角从 34.7 度涨到 82.4 度。 数据越接近低维子空间，协方差矩阵越病态，它的逆越会把切片均值的采样噪声放大。 必须加岭项。
4. 「它和 PCA 一样，都是找方差结构」。 它找的是「条件均值变化的方向」，和方差结构不是一回事。 本节里方差最大的方向恰好是干扰，所以两者的答案几乎正交（79.9 度）。

17 章两节都在假设一件事：训练和测试的数据来自同一个分布。 下一章把这个假设也拆掉——如果训练时见的和上线后遇到的分布不一样，该怎么办。

• ☐ （本章完，下一节 18.1 18.1）

相关笔记

• 与分类问题相对应的判别分析
• 线性降维的原理
• 灵活应用输入数据的流形构造
• 🤖 图解机器学习 · 索引