一、问题长什么样

先把规模摆出来：

Y = X * A * B' + 噪声       A 是 8 乘 2，B 是 20 乘 2
所以系数矩阵 W = A * B' 的秩最多是 2

20 个样本估 160 个系数，是典型的「参数比数据多」。 这种情况最小二乘还能解（加个微小岭项防止奇异），但解出来的东西基本是噪声。

二、逐维独立回归栽在哪里

跑出来是这样：

三、奇异值给出的线索

怎么知道真秩是 2？答案在奇异值里。

四、降秩回归怎么做

思路直接：既然秩应该小，那就强行把秩卡小。 但怎么卡才不是随便砍？做法分三步：

1. 把输入白化：乘上 X'X 的负二分之一次方，让 X 的协方差变成单位阵
2. 在白化坐标里算普通最小二乘，得到 W_w
3. 对 W_w 做奇异值分解，只留前 r 个，再变换回原坐标

为什么要先白化？因为秩约束要和「预测误差」这个尺度对齐。 不白化直接截断最小二乘解的奇异值，截掉的是「系数大不大」， 而不是「对预测贡献大不大」——输入各个方向的方差不一致时，这两件事不是一回事。 白化之后，奇异值大小才真正对应「这个方向解释了多少输出变化」。 卡不同秩的结果：

五、秩选小了会更糟

秩 1 那一行是 4.1470，是秩 2 的 77.2 倍，比不卡秩还差 61 倍。 真结构有两维，硬压到一维等于把一半的信息彻底扔掉。 模型再正则化也补不回来——因为它连表达真函数的能力都没有了。 这件事说明降秩回归的风险是双向的。 秩太大 → 过拟合，测试误差慢慢变差（0.0537 到 0.0726）； 秩太小 → 欠拟合，测试误差爆炸（0.0537 到 4.1470）。 左边陡、右边缓，所以宁可秩大一点也不要太小。 还有一个数值得看：20 个输出两两之间的平均相关系数绝对值是 0.6261。 输出之间本来就有强相关，逐维独立回归完全不管这层关系， 等于把同一个信号重复估了 20 遍——每一遍都掺进独立的噪声。 降秩回归做的事情，本质上是「一次估一个共同的低维结构，再摊开」，噪声只在低维那一步掺一次。

六、动手 Demo：多次维输出：先压到低维再摊开，秩卡在断崖处

这段代码造 20 个样本、8 维输入、20 维输出，真结构是「先压到 2 维再摊开」； 然后比较逐维独立回归和不同秩的降秩回归，最后看秩选小了的代价。 重点看三处。第二段训练集 0.0204 低于噪声下限 0.04 而测试集 0.0675； 第三段奇异值从 7.986 断到 0.331；第四段测试集在秩 2 处最低（0.0537），第五段秩 1 崩到 4.1470。

# -*- coding: utf-8 -*-
import numpy as np


def rr_data(n=20, d=8, m=20, r=2, seed=5, noise=0.2):
    """输入 8 维，输出 20 维，但真正的中间变量只有 2 维。"""
    rg = np.random.default_rng(seed)
    A = rg.normal(0.0, 1.0, (d, r))
    B = rg.normal(0.0, 1.0, (m, r))
    X = rg.normal(0.0, 1.0, (n, d))
    Y = X @ A @ B.T + rg.normal(0.0, noise, (n, m))
    return X, Y, A, B


def ols_multi(X, Y):
    """每个输出维度各解一个最小二乘，等价于一次矩阵求解。"""
    return np.linalg.solve(X.T @ X + 1e-8 * np.eye(X.shape[1]), X.T @ Y)


def reduced_rank(X, Y, r):
    """降秩回归：先把输入白化，再对最小二乘解做截断奇异值分解。"""
    n, d = X.shape
    mx, my = X.mean(0), Y.mean(0)
    Xc, Yc = X - mx, Y - my
    S = Xc.T @ Xc / n
    wv, V = np.linalg.eigh(S)
    Wh = V @ np.diag(1.0 / np.sqrt(np.maximum(wv, 1e-10))) @ V.T
    Xw = Xc @ Wh
    Ww = Xw.T @ Yc / n
    U, s, Vt = np.linalg.svd(Ww, full_matrices=False)
    Wr = (U[:, :r] * s[:r]) @ Vt[:r]
    return Wh @ Wr, my - mx @ (Wh @ Wr)


def err(W, b, X, Y):
    return float(np.mean((X @ W + b - Y) ** 2))


X, Y, A, B = rr_data()

# 测试集必须用同一组 A、B 生成，否则学的和考的不是同一个任务
rg = np.random.default_rng(555)
Xt = rg.normal(0.0, 1.0, (400, 8))
Yt = Xt @ A @ B.T + rg.normal(0.0, 0.2, (400, 20))
print('')

print('[1] 20 个样本，输入 8 维，输出 20 维')
print('  真结构是「输入先压到 2 维，再摊到 20 维」，也就是 8 乘 20 的系数矩阵秩只有 2')
print('  20 个样本要估 160 个系数，逐维独立回归必然过拟合')
print('  噪声标准差 0.20，预测误差的下限是 0.04')

print('')
print('[2] 逐维独立回归：训练集几乎完美，换一批数据就露馅')
Wo = ols_multi(X, Y)
bo = Y.mean(0) - X.mean(0) @ Wo
print('  训练集 %.4f   测试集 %.4f' % (err(Wo, bo, X, Y), err(Wo, bo, Xt, Yt)))
print('  系数矩阵的秩 %d，该是 2' % int(np.linalg.matrix_rank(Wo, tol=1e-8)))

print('')
print('[3] 看奇异值，断崖在哪里')
sv = np.linalg.svd(Wo, compute_uv=False)
print('  前八个奇异值 ' + '  '.join('%.3f' % v for v in sv[:8]))
print('  第三个开始断崖式下跌，说明真正有用的方向只有两个')

print('')
print('[4] 把秩卡到 2：训练集略差，测试集明显好')
print('  秩 r   训练集   测试集')
for r in (1, 2, 3, 4, 8):
    Wr, br = reduced_rank(X, Y, r)
    print('  %d      %.4f   %.4f' % (r, err(Wr, br, X, Y), err(Wr, br, Xt, Yt)))

W2, b2 = reduced_rank(X, Y, 2)
print('  秩卡到 2 之后，测试集从 %.4f 降到 %.4f，改善 %.2f 倍'
      % (err(Wo, bo, Xt, Yt), err(W2, b2, Xt, Yt), err(Wo, bo, Xt, Yt) / err(W2, b2, Xt, Yt)))

print('')
print('[5] 秩选错的代价：卡到 1 会比不卡还差得多')
W1, b1 = reduced_rank(X, Y, 1)
print('  r=1 测试集 %.4f，是 r=2 的 %.1f 倍' % (err(W1, b1, Xt, Yt), err(W1, b1, Xt, Yt) / err(W2, b2, Xt, Yt)))
print('  真结构有两维，硬压到一维等于把一半的信息扔掉')

print('')
print('[6] 输出之间本来就有强相关')
Yc = Y - Y.mean(0)
cc = np.corrcoef(Yc.T)
print('  20 个输出两两之间的平均相关系数绝对值 %.4f' % np.mean(np.abs(cc[np.triu_indices(20, 1)])))
print('  逐维独立回归完全不管这层相关，等于把同一个信号重复估了 20 遍')

实跑输出：

--------------------------------------------------------------------
多次维输出：先压到低维再摊开，秩卡在断崖处
--------------------------------------------------------------------

[1] 20 个样本，输入 8 维，输出 20 维
  真结构是「输入先压到 2 维，再摊到 20 维」，也就是 8 乘 20 的系数矩阵秩只有 2
  20 个样本要估 160 个系数，逐维独立回归必然过拟合
  噪声标准差 0.20，预测误差的下限是 0.04

[2] 逐维独立回归：训练集几乎完美，换一批数据就露馅
  训练集 0.0204   测试集 0.0675
  系数矩阵的秩 8，该是 2

[3] 看奇异值，断崖在哪里
  前八个奇异值 13.385  7.986  0.331  0.289  0.243  0.157  0.138  0.116
  第三个开始断崖式下跌，说明真正有用的方向只有两个

[4] 把秩卡到 2：训练集略差，测试集明显好
  秩 r   训练集   测试集
  1      2.6993   4.1470
  2      0.0320   0.0537
  3      0.0279   0.0602
  4      0.0248   0.0665
  8      0.0197   0.0726
  秩卡到 2 之后，测试集从 0.0675 降到 0.0537，改善 1.26 倍

[5] 秩选错的代价：卡到 1 会比不卡还差得多
  r=1 测试集 4.1470，是 r=2 的 77.2 倍
  真结构有两维，硬压到一维等于把一半的信息扔掉

[6] 输出之间本来就有强相关
  20 个输出两两之间的平均相关系数绝对值 0.6261
  逐维独立回归完全不管这层相关，等于把同一个信号重复估了 20 遍

第一段那句「20 个样本要估 160 个系数」是整节的引子。 它把「参数比数据多」这件事量化了：每个样本要支撑 8 个系数。 这种比例下，任何不带结构约束的方法都注定失败。 第二段训练集 0.0204 低于噪声 0.04，这个反差要停一下。 它说明模型把噪声也当成信号学进去了。 判断过拟合最硬的标准不是「训练和测试差多少」，而是「训练误差是不是低于噪声水平」。 第三段那八个奇异值，前两个和后面差了 24 倍。 在真实数据里，这个比值往往只有两三倍，断崖要靠经验或者交叉验证去认。 本节数据造得干净，所以断崖这么明显——这是教学上的便利，不是普遍现象。 第四段那张表最该横着读一遍。 训练集 0.0320、0.0279、0.0248、0.0197 一路降，测试集 0.0537、0.0602、0.0665、0.0726 一路升。 两条线反着走，中间那个交点就是该取的秩。 第五段秩 1 的 4.1470 是这一节最有冲击力的数字。 它不是「差一点」，是「彻底废掉」——比不卡秩的 0.0675 差 61 倍。 正则化不是「越强越好」，结构约束一旦强过真实结构，代价是断崖式的。 代码里 reduced_rank 的白化那一步用的是 eigh 而不是 svd。 因为协方差矩阵是对称的，eigh 更快更稳，而且特征值天然有序。 那个 np.maximum(wv, 1e-10) 是防止输入某个方向方差恰好为零时除零。

30 秒实验

1. 把 rr_data 的 r 从 2 改成 3，重跑。奇异值的断崖会移到第三个位置，最优秩也跟着变成 3。 这说明「奇异值断崖」这个判据是真的在工作，不是碰巧。
2. 把 rr_data 的 n 从 20 改成 60，重跑。逐维独立回归会追上降秩回归。 样本一多，秩约束就不必要了——它本来就是为「参数比数据多」这种情况准备的。
3. 把 rr_data 的 noise 从 0.2 改成 0.6，重跑。奇异值的断崖会变钝，最优秩变得难认。 噪声越大，断崖越模糊，这是这类方法在真实数据上最常遇到的困难。

常见误区

1. 「多个输出就逐个独立回归，简单可靠」。 实测 20 个输出、20 个样本时，独立回归测试误差 0.0675，降秩到 2 之后 0.0537。 输出之间的平均相关系数绝对值是 0.6261，这层关系被白白丢掉了。 参数越多，丢得越狠。
2. 「秩越小，模型越简单越安全」。 实测秩 1 的测试误差 4.1470，是秩 2 的 77.2 倍，比不卡秩还差 61 倍。 结构约束强过真实结构时，模型连表达真函数的能力都没有，再正则化也救不回来。 左边陡右边缓，宁可秩偏大。
3. 「降秩之前不用白化」。 白化是为了让奇异值的大小对应「解释了多少输出变化」，而不是「系数大不大」。 输入各方向方差不一致时，这两件事完全不同。 不白化直接截断，截掉的可能是有用的方向。
4. 「训练误差低于测试误差就是过拟合了」。 这个说法太粗。 更准的判据是「训练误差是不是低于噪声水平」。 本节训练集 0.0204 低于噪声 0.04，这才是过拟合的铁证。 而如果训练误差是 0.045、测试 0.05，两者都高于噪声，那只是正常的学习误差。

第 19 章把「多个任务」的几种形态都过了一遍。 到这里，这本书的主要方法就讲完了。 最后一章不再介绍新算法，而是回头看一眼：这些方法之间的关系是什么，以及它们各自成立的前提。

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

相关笔记

• 使用最小二乘概率分类器的多任务学习
• 最小二乘解的性质
• 主成分分析
• 🤖 图解机器学习 · 索引