---- 手造数据: 5 个样本, 2 个特征; x2 几乎 = 2*x1 (偏差手工指定, 无随机) ----

岭回归：会”缩水”的最小二乘

小象实战讲义 · 机器学习算法精讲

上一节结尾留了个悬念：特征一多、$X^TX$ 求不出逆，正规方程就翻车了——讲师说，需要往里面加一些更多的……加的是什么？就是 regularization（正则化）。这个词在汉语里经常翻译成”正则化”，英文里就是 regulate：给更多的规则、更多的法则——对参数给一些约束。今天我们仍然在研究线性模型、线性拟合，唯一的新东西，就是在损失函数上多加一个”限制”。第一个登场的方法叫 ridge regression（岭回归）——ridge，就是一个山脊、一道山岭。

💡 核心导读

• 动机：$X^TX$ 不可逆时正规方程失守——往损失函数里加东西，把问题救回来。
• 岭回归损失：最小二乘损失 + 参数的 L2 范数平方惩罚（$\theta_1$ 到 $\theta_n$，不含偏置 $\theta_0$）。
• 两个视角：带约束的优化（拉格朗日乘子）与梯度下降的”先缩水、再更新”。
• 解析解：$\hat{\theta} = (X^TX+\lambda I')^{-1}X^Ty$，其中 $I'$ 是把单位矩阵左上角的 1 换成 0——千万不要写成对角线全 1。
• 实战口径：真正使用的最小二乘，一般都带着 Ridge。

回顾：记号与那个求不出的逆

每次课都把记号重新钉一遍，不要忘了。训练数据 $D = \lbrace(x^{(1)}, y^{(1)}), \dots, (x^{(N)}, y^{(N)})\rbrace$ 是监督学习的数据：大 $N$ 是数据个数，每个 $x^{(i)}$ 属于小 $n$ 维。给每个数据前头添一个恒为 1 的 $x_0$，把所有数据按行堆成大 $X$——上角标表示第几个数据，下角标表示第几个特征。所以 $X$ 的形状是 $N \times (n+1)$：行有 $N$ 个数据，列有 $n+1$ 维（每行开头那个 1）。参数是 $\theta = (\theta_0, \theta_1, \dots, \theta_n)$，输出堆成向量 $y$。因为我们默认每个向量都是列向量，所以把 $x^{(i)}$ 写成横排时要加转置：$x^{(i)T} = [x_0^{(i)}, x_1^{(i)}, \dots, x_n^{(i)}]$。

我们希望 $X\theta = y$，但实际上经常得不到，只能约等于。那就选一个让”最接近”的 $\theta$：让预测 $\hat{y} = X\theta$ 与真实 $y$ 之间的距离最小（几何上就是第 8 章讲过的投影）。距离用 L2 范数，损失写成

$$J(\theta) = \frac{1}{2N}\Vert X\theta - y\Vert _2^2$$

（前面乘 $1/2$ 还是 $1/N$，并不重要。）对 $\theta$ 求导置零，得到上次的最小二乘解：$\hat{\theta} = (X^TX)^{-1}X^Ty$。

问题就出在那个求逆上。求一个矩阵的逆叫 inverse，inverse 的形容词是 invertible，non-invertible 就是不能求逆的。上一节说过：特征数 $n$ 大于样本数 $N$ 时，$X^TX$ 多半是奇异矩阵；特征之间高度相关（近似线性相关）时也一样——第 8 章练一练里那道 $A^TA$ 行列式恰为 0 的题，就是这种局面。如果它不能求逆，怎么办？ 我们在原有的基础之上加东西。

Ridge：在损失函数上加限制

Ridge regression 做的事，跟以前没有任何区别——唯一的区别是在损失函数上加一些限制条件：希望每个系数不要特别大，把它限制到某个具体范围之内。

具体做法：原来那部分不动，加上一个惩罚项（正则项）：

$$J(\theta) = \frac{1}{2N}\Big(\Vert X\theta - y\Vert _2^2 + \lambda\Vert \theta_{1:n}\Vert _2^2\Big)$$

逐块读：$\Vert X\theta-y\Vert _2^2$ 还是原来的损失；$\Vert \theta_{1:n}\Vert _2^2 = \theta_1^2 + \theta_2^2 + \cdots + \theta_n^2$ 是所有参数的平方和。注意这里的 $\theta$ 从 $\theta_1$ 加到 $\theta_n$，不包含 $\theta_0$——$\theta_0$ 是偏置（bias），先不考虑它，不惩罚。范数家族（L1/L2 的老朋友们）在第 5 章《距离测算》里见过。

为什么是平方？因为不管 $\theta_j$ 是正是负，平方都是一个正数——我们实际上希望它的绝对值比较小。如果其中某个 $\theta_j$ 特别大，整个损失立刻变大。所以新的损失函数要同时照顾两件事：一方面选取合适的 $\theta$，让预测与期待的 $y$ 距离最小；另一方面不能让 $\theta$ 本身的值过大——找一个不是很大、但还能尽量满足你所需条件的 $\theta$。

一个记号说明：课件与板书把惩罚项的系数记作 $\alpha$，更新公式里又改用 $\lambda$——因为学习率也叫 $\alpha$，讲师特意换了字母：“这个阿尔法跟刚才的阿尔法不是一个，它是学习速率。“本讲义统一用 $\lambda$；$\lambda$ 放在哪里（除不除 $N$、乘不乘 $1/2$）只改变 $\lambda$ 的读数尺度，不改变任何结论——以讲义为准。

$\lambda$ 起的作用，是调整两边的相对权重：$\lambda$ 特别大时，$\theta$ 一变大整个函数就变得非常大，于是更偏重”别把 $\theta$ 变大”；$\lambda$ 小的时候，惩罚项没那么重要，更关注前面那部分——是否把 $y$ 拟合好。课件把这件事画成一条 U 形曲线：横轴是 $\lambda$ 增大，纵轴是误差——最左端 $\lambda=0$，只顾拟合损失，过拟合；最右端 $\lambda \to \infty$，$\theta \to 0$，只顾惩罚，欠拟合；中间有个最低点，就是 $\lambda$ 该待的地方。图上还画着一条水平线，标着 Linear Regression——普通最小二乘所在的误差水平。

换个写法：带约束的优化

同一件事还可以写成两个部分：先优化原来的损失，同时满足一个条件——

$$\min_{\theta}\ \frac{1}{2N}\Vert X\theta-y\Vert _2^2 \quad \text{s.t.}\quad \Vert \theta_{1:n}\Vert _2^2 \le t$$

即”参数的平方和不许超过某个数 $t$“。这种有条件的优化可以构造拉格朗日函数：给定 $\min f(x)$、约束 $g(x)=0$，就构造 $L(x,\lambda) = f(x) + \lambda\, g(x)$——把约束乘上一个乘子加进目标，变成无约束的问题。两相对照：惩罚项里那个 $\lambda$，扮演的正是拉格朗日乘子的角色。拉格朗日函数本身我们在第 1 章就打过照面。

梯度下降：先缩水，再走老路

怎么求解？还是老朋友梯度下降。第 4 章的更新公式：每个分量

$$\theta_j := \theta_j - \alpha\cdot\frac{1}{N}\sum_{i=1}^{N}\big(h_\theta(x^{(i)}) - y^{(i)}\big)x_j^{(i)}$$

其中 $\theta_0$ 是特例——$x_0$ 恒为 1，所以不乘后面任何东西。现在多了惩罚项，对 $\frac{\lambda}{2}\theta_j^2$ 求导，平方正好约掉，正好剩下一个 $\lambda\theta_j$（带着 $1/N$ 的口径就是 $\frac{\lambda}{N}\theta_j$）。于是对 $j \ge 1$：

$$\theta_j := \theta_j - \alpha\Big[\frac{1}{N}\sum_{i=1}^{N}\big(h_\theta(x^{(i)})-y^{(i)}\big)x_j^{(i)} + \frac{\lambda}{N}\theta_j\Big] = \Big(1-\frac{\alpha\lambda}{N}\Big)\theta_j - \alpha\cdot\frac{1}{N}\sum_{i=1}^{N}\big(h_\theta(x^{(i)})-y^{(i)}\big)x_j^{(i)}$$

$\theta_0$ 的更新公式完全不变——它没有被惩罚。

把缩水因子单拎出来看，讲师给了一段数值直觉：$\alpha$ 是学习速率，一般是 0.1、0.5 这样的小数；$\lambda$ 代表两者相对权重，也是个不太大的正数；而 $N$ 是数据个数——成千上万，多的还有几百万。所以 $\alpha\lambda/N$ 经常是一个很小的数，比如 0.01、0.02 这个量级。也就是说，每次更新其实是：先把原有的 $\theta_j$ 乘上 0.99 或 0.98 或 0.97 缩小一下，再减去与从前完全一样的更新量。

“缩水”这个词就来自这里——讲师的原话：“shrink，就是缩水，叫 shrinkage。“（教科书里叫参数收缩。）每轮迭代 = 缩一次水 + 做一次老更新，这就是 ridge regression 在梯度下降里做的事。

解析解：留空的矩阵填上了

不用梯度下降、按上次最小二乘的方式直接推导解析表达，也能得到解。推导过程上次讲过，这次不从头推了——“否则这节课也不够”，了解背后的道理即可。上次的结果是 $\hat{\theta} = (X^TX)^{-1}X^Ty$；讲师板书时在这个公式中间留了一个空——填进去的正是那个惩罚：

$$\hat{\theta}_{\text{ridge}} = (X^TX + \lambda I')^{-1}X^Ty$$

这里有一个必须注意的小 trick，先做个维度检查：$X$ 是 $N\times(n+1)$，所以 $X^TX$ 是 $(n+1)\times(n+1)$；那这个 $I'$ 也必须是 $(n+1)\times(n+1)$。但它跟普通的单位矩阵稍稍有点区别——最上面（左上角）是 0，其余对角线是 1：

$$I' = \begin{bmatrix} 0 & 0 & \cdots & 0 \\ 0 & 1 & \cdots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & \cdots & 1 \end{bmatrix}$$

它并不完全是一个 identity matrix（单位矩阵：对角线是 1、其他地方都是零），相当于把单位矩阵左上角的 1 换成了 0。为什么？因为 $\theta_0$ 跟其他参数不一样——$\theta_0$ 不需要加惩罚项。讲师在这里特意加重语气：真正做计算的时候，千万不要把 $I'$ 写成对角线全 1——那是不对的。

这个公式还有一层特别重要的好处：以前那些求不出逆的矩阵——非满秩的、奇异的——加上这个 $\lambda I'$ 之后，很多矩阵就重新有逆了，可以直接计算。这就是 regularization 一个特别重要的好处。

下面把三件事跑给你看：$\lambda$ 怎么让系数缩水、$\lambda$ 怎么救活奇异矩阵、$I'$ 的左上角写错会怎样。

代码一：$\lambda$ 扫描，看系数缩水

手造一份近共线数据：$x_2$ 几乎等于 $2x_1$（偏差只有百分之几，手工指定、无随机）；两批 $y$ 的主趋势完全相同，只是”毛刺”方向相反——毛刺的形状恰好长在 $x_2 - 2x_1$ 这条近零方向上：

# ---- 手造数据: 5 个样本, 2 个特征; x2 几乎 = 2*x1 (偏差手工指定, 无随机) ----
X = [
    [1.0, 0.0, 0.02],
    [1.0, 0.1, 0.19],
    [1.0, 0.2, 0.41],
    [1.0, 0.3, 0.58],
    [1.0, 0.4, 0.82],
]
# 两批 y: 主趋势都是 2.5 + 0.2*x2, 只是"毛刺"方向相反 (毛刺形状 = x2 - 2*x1 的形状)
y_A = [2.604, 2.488, 2.632, 2.516, 2.764]
y_B = [2.404, 2.588, 2.532, 2.716, 2.564]

def transpose(A):
    return [list(col) for col in zip(*A)]

def matmul(A, B):
    return [[sum(A[i][k] * B[k][j] for k in range(len(B))) for j in range(len(B[0]))] for i in range(len(A))]

def matvec(A, v):
    return [sum(A[i][k] * v[k] for k in range(len(v))) for i in range(len(A))]

def det3(M):
    return (M[0][0] * (M[1][1] * M[2][2] - M[1][2] * M[2][1])
          - M[0][1] * (M[1][0] * M[2][2] - M[1][2] * M[2][0])
          + M[0][2] * (M[1][0] * M[2][1] - M[1][1] * M[2][0]))

def solve(A, b):
    # 高斯消元(列主元); 零主元 => 矩阵奇异, 返回 None
    n = len(A)
    aug = [row[:] + [b[i]] for i, row in enumerate(A)]
    for i in range(n):
        piv = max(range(i, n), key=lambda r: abs(aug[r][i]))
        if abs(aug[piv][i]) < 1e-12:
            return None
        aug[i], aug[piv] = aug[piv], aug[i]
        for r in range(i + 1, n):
            f = aug[r][i] / aug[i][i]
            for c in range(i, n + 1):
                aug[r][c] -= f * aug[i][c]
    x = [0.0] * n
    for i in range(n - 1, -1, -1):
        x[i] = (aug[i][n] - sum(aug[i][j] * x[j] for j in range(i + 1, n))) / aug[i][i]
    return x

XT = transpose(X)
XTX = matmul(XT, X)
XTy_A = matvec(XT, y_A)
XTy_B = matvec(XT, y_B)

print("X^T X (对称):", [[round(v, 4) for v in r] for r in XTX])
print(f"det(X^T X) = {det3(XTX):.6f}  -- 三个特征值相乘, 有一个被压到接近 0 (x2 几乎 = 2*x1)")

print("\nlambda 扫描 (I' = diag(0,1,1), 不惩罚 theta0):")
for lam in [0.0, 0.5, 2.0, 10.0]:
    A = [[XTX[i][j] + (lam if i == j and i > 0 else 0.0) for j in range(3)] for i in range(3)]
    tA = solve(A, XTy_A)
    tB = solve(A, XTy_B)
    print(f"  lambda = {lam:4.1f}: 批A theta = [{tA[0]:7.4f}, {tA[1]:8.4f}, {tA[2]:8.4f}]"
          f"   批B theta = [{tB[0]:7.4f}, {tB[1]:8.4f}, {tB[2]:8.4f}]")

# 验证 1: lambda=2 的解确实站在正规方程的零点
A2 = [[XTX[i][j] + (2.0 if i == j and i > 0 else 0.0) for j in range(3)] for i in range(3)]
theta2 = solve(A2, XTy_A)
resid = [a - b for a, b in zip(matvec(A2, theta2), XTy_A)]
print("\nlambda=2 验证 (X^T X + 2I')theta - X^T y =", [f"{v:.1e}" for v in resid])

# 验证 2: lambda=0 那组离谱参数, 在训练点上竟是"完美拟合"
theta0 = solve(XTX, XTy_A)
diff = [a - b for a, b in zip(matvec(X, theta0), y_A)]
print(f"lambda=0 批A 在 5 个训练点上的预测残差: {[f'{v:.1e}' for v in diff]}")
print("=> 两批数据各被拟合到机器零, theta1 却相差 20: 这就是近零方向上的自由漂移")
print("   加了 lambda 之后, 对冲消失, 两批的系数都缩回到千分位量级")

X^T X (对称): [[5.0, 1.0, 2.02], [1.0, 0.3, 0.603], [2.02, 0.603, 1.2134]]
det(X^T X) = 0.000655  -- 三个特征值相乘, 有一个被压到接近 0 (x2 几乎 = 2*x1)

lambda 扫描 (I' = diag(0,1,1), 不惩罚 theta0):
  lambda =  0.0: 批A theta = [ 2.5000, -10.0000,   5.2000]   批B theta = [ 2.5000,  10.0000,  -4.8000]
  lambda =  0.5: 批A theta = [ 2.5630,   0.0323,   0.0776]   批B theta = [ 2.5182,   0.0475,   0.0818]
  lambda =  2.0: 批A theta = [ 2.5857,   0.0137,   0.0306]   批B theta = [ 2.5438,   0.0182,   0.0331]
  lambda = 10.0: 批A theta = [ 2.5972,   0.0033,   0.0073]   批B theta = [ 2.5568,   0.0043,   0.0079]

lambda=2 验证 (X^T X + 2I')theta - X^T y = ['0.0e+00', '0.0e+00', '0.0e+00']
lambda=0 批A 在 5 个训练点上的预测残差: ['-1.5e-14', '1.3e-14', '-5.8e-15', '2.2e-14', '-1.6e-14']
=> 两批数据各被拟合到机器零, theta1 却相差 20: 这就是近零方向上的自由漂移
   加了 lambda 之后, 对冲消失, 两批的系数都缩回到千分位量级

三个观察。$\lambda=0$ 的疯狂对冲：两批 $y$ 只是毛刺方向相反，$\theta_1$ 却从 $-10$ 跳到 $+10$、$\theta_2$ 从 $+5.2$ 跳到 $-4.8$——系数沿着 $x_2 \approx 2x_1$ 这条近零方向自由漂移（$\det(X^TX)=0.000655$，有一个特征值被压到接近零），而两批数据都被拟合到 $10^{-14}$ 的机器零。这就是”希望系数不要特别大”的真实含义：不约束，解就在悬崖边上跳舞。 $\lambda$ 一开，全缩回来：0.5、2、10 一路扫上去，两批的系数都缩到千分位，而且两批之间的差异也跟着消失——毛刺方向那点自由度被惩罚买断了。解验证：$\lambda=2$ 的解代回正规方程 $(X^TX+2I')\theta - X^Ty$，残差精确为 0——它确实站在新损失函数的零点上。

代码二：奇异矩阵的复活

把偏差去掉，让 $x_2 = 2x_1$ 严格成立——$X$ 的第三列恰好等于第二列乘 2，$X^TX$ 彻底奇异：

# ---- 完全共线数据: x2 = 2*x1 严格成立, X 的第三列 = 2 * 第二列 ----
X = [
    [1.0, 0.0, 0.0],
    [1.0, 0.1, 0.2],
    [1.0, 0.2, 0.4],
    [1.0, 0.3, 0.6],
    [1.0, 0.4, 0.8],
]
y = [2.604, 2.488, 2.632, 2.516, 2.764]   # 沿用批A

def transpose(A):
    return [list(col) for col in zip(*A)]

def matmul(A, B):
    return [[sum(A[i][k] * B[k][j] for k in range(len(B))) for j in range(len(B[0]))] for i in range(len(A))]

def matvec(A, v):
    return [sum(A[i][k] * v[k] for k in range(len(v))) for i in range(len(A))]

def det3(M):
    return (M[0][0] * (M[1][1] * M[2][2] - M[1][2] * M[2][1])
          - M[0][1] * (M[1][0] * M[2][2] - M[1][2] * M[2][0])
          + M[0][2] * (M[1][0] * M[2][1] - M[1][1] * M[2][0]))

def solve(A, b):
    # 高斯消元(列主元); 零主元 => 矩阵奇异, 返回 None
    n = len(A)
    aug = [row[:] + [b[i]] for i, row in enumerate(A)]
    for i in range(n):
        piv = max(range(i, n), key=lambda r: abs(aug[r][i]))
        if abs(aug[piv][i]) < 1e-12:
            return None
        aug[i], aug[piv] = aug[piv], aug[i]
        for r in range(i + 1, n):
            f = aug[r][i] / aug[i][i]
            for c in range(i, n + 1):
                aug[r][c] -= f * aug[i][c]
    x = [0.0] * n
    for i in range(n - 1, -1, -1):
        x[i] = (aug[i][n] - sum(aug[i][j] * x[j] for j in range(i + 1, n))) / aug[i][i]
    return x

XT = transpose(X)
XTX = matmul(XT, X)
XTy = matvec(XT, y)

print("X^T X (对称):", [[round(v, 4) for v in r] for r in XTX])
print(f"det(X^T X) = {det3(XTX):.6f}  -- 第三行 = 2 * 第二行, 矩阵奇异")

print("\nlambda 扫描:")
for lam in [0.0, 0.1, 1.0]:
    A = [[XTX[i][j] + (lam if i == j and i > 0 else 0.0) for j in range(3)] for i in range(3)]
    theta = solve(A, XTy)
    if theta is None:
        print(f"  lambda = {lam:3.1f}: 高斯消元撞上零主元 -> (X^T X)^(-1) 不存在, 正规方程失守")
    else:
        print(f"  lambda = {lam:3.1f}: theta = [{theta[0]:7.4f}, {theta[1]:8.4f}, {theta[2]:8.4f}]"
              f"   det(X^T X + lambda I') = {det3(A):.4f}")

# lambda=0.1 的解在训练点上的拟合 (x2 方向的毛刺没人能拟合, 残差诚实留在那里)
theta = solve([[XTX[i][j] + (0.1 if i == j and i > 0 else 0.0) for j in range(3)] for i in range(3)], XTy)
pred = matvec(X, theta)
mse = sum((a - b) ** 2 for a, b in zip(pred, y)) / len(y)
print(f"\nlambda=0.1 的训练 MSE = {mse:.4f} (x1 与 x2 完全共线, 毛刺方向无人能拟合, 残差留在原地)")
print("=> lambda 买回来的是: 唯一、有界的解 + 一个可以求逆的矩阵")

X^T X (对称): [[5.0, 1.0, 2.0], [1.0, 0.3, 0.6], [2.0, 0.6, 1.2]]
det(X^T X) = 0.000000  -- 第三行 = 2 * 第二行, 矩阵奇异

lambda 扫描:
  lambda = 0.0: 高斯消元撞上零主元 -> (X^T X)^(-1) 不存在, 正规方程失守
  lambda = 0.1: theta = [ 2.5428,   0.0580,   0.1160]   det(X^T X + lambda I') = 0.3000
  lambda = 1.0: theta = [ 2.5776,   0.0232,   0.0464]   det(X^T X + lambda I') = 7.5000

lambda=0.1 的训练 MSE = 0.0072 (x1 与 x2 完全共线, 毛刺方向无人能拟合, 残差留在原地)
=> lambda 买回来的是: 唯一、有界的解 + 一个可以求逆的矩阵

$\det(X^TX) = 0$ 精确成立（第三行 = 2×第二行），$\lambda=0$ 时高斯消元撞上零主元——$(X^TX)^{-1}$ 不存在，这就是上一节”翻车现场”的完整复现。$\lambda=0.1$ 一加上，$\det = 0.3000$；$\lambda=1$ 时 $\det = 7.5$——矩阵活了，解唯一且有界。这正对应讲师说的”以前那些非满秩的、奇异的矩阵，加上这个控制项就可以直接求逆”。还有一处诚实的细节：$x_1$ 与 $x_2$ 完全共线，$y$ 里那个沿 $x_2-2x_1$ 方向的毛刺谁也拟合不了，训练 MSE = 0.0072 的残差就留在原地——$\lambda$ 买回来的是可解性，不是魔法。顺带一提，解里 $\theta_2 = 2\theta_1$：完全共线时，ridge 把拟合任务沿惩罚最小的方向分摊给两个系数。

代码三：千万别把左上角写成 1

讲师那句”千万不要把 $I'$ 写成对角线全 1”值得单独跑一次——同一份数据、同一个 $\lambda$，只对比两种写法：

# ---- 同一份近共线数据, 对比两种 I 矩阵: 正确的 I' (左上角 0) vs 错误的 I (对角全 1) ----
X = [
    [1.0, 0.0, 0.02],
    [1.0, 0.1, 0.19],
    [1.0, 0.2, 0.41],
    [1.0, 0.3, 0.58],
    [1.0, 0.4, 0.82],
]
y = [2.604, 2.488, 2.632, 2.516, 2.764]

def transpose(A):
    return [list(col) for col in zip(*A)]

def matmul(A, B):
    return [[sum(A[i][k] * B[k][j] for k in range(len(B))) for j in range(len(B[0]))] for i in range(len(A))]

def matvec(A, v):
    return [sum(A[i][k] * v[k] for k in range(len(v))) for i in range(len(A))]

def solve(A, b):
    n = len(A)
    aug = [row[:] + [b[i]] for i, row in enumerate(A)]
    for i in range(n):
        piv = max(range(i, n), key=lambda r: abs(aug[r][i]))
        if abs(aug[piv][i]) < 1e-12:
            return None
        aug[i], aug[piv] = aug[piv], aug[i]
        for r in range(i + 1, n):
            f = aug[r][i] / aug[i][i]
            for c in range(i, n + 1):
                aug[r][c] -= f * aug[i][c]
    x = [0.0] * n
    for i in range(n - 1, -1, -1):
        x[i] = (aug[i][n] - sum(aug[i][j] * x[j] for j in range(i + 1, n))) / aug[i][i]
    return x

XT = transpose(X)
XTX = matmul(XT, X)
XTy = matvec(XT, y)
lam = 2.0

# 正确: I' = diag(0, 1, 1) -- theta0 不惩罚
A_right = [[XTX[i][j] + (lam if i == j and i > 0 else 0.0) for j in range(3)] for i in range(3)]
t_right = solve(A_right, XTy)

# 错误: I = diag(1, 1, 1) -- 连 theta0 一起惩罚
A_wrong = [[XTX[i][j] + (lam if i == j else 0.0) for j in range(3)] for i in range(3)]
t_wrong = solve(A_wrong, XTy)

print(f"lambda = {lam}:")
print(f"  正确 I' = diag(0,1,1): theta = [{t_right[0]:7.4f}, {t_right[1]:8.4f}, {t_right[2]:8.4f}]")
print(f"  错误 I  = diag(1,1,1): theta = [{t_wrong[0]:7.4f}, {t_wrong[1]:8.4f}, {t_wrong[2]:8.4f}]")

# 偏置被压小 => 整条预测线整体下沉
mean_y = sum(y) / len(y)
mean_right = sum(matvec(X, t_right)) / len(y)
mean_wrong = sum(matvec(X, t_wrong)) / len(y)
print(f"\n训练集均值: y 平均 {mean_y:.4f} | 正确 I' 预测平均 {mean_right:.4f} | 错误 I 预测平均 {mean_wrong:.4f}")
print(f"=> 错误写法把 theta0 从 {t_right[0]:.4f} 压到 {t_wrong[0]:.4f}, 全部预测值整体偏移 {mean_wrong - mean_right:.4f}")

lambda = 2.0:
  正确 I' = diag(0,1,1): theta = [ 2.5857,   0.0137,   0.0306]
  错误 I  = diag(1,1,1): theta = [ 1.6546,   0.2781,   0.5663]

训练集均值: y 平均 2.6008 | 正确 I' 预测平均 2.6008 | 错误 I 预测平均 1.9390
=> 错误写法把 theta0 从 2.5857 压到 1.6546, 全部预测值整体偏移 -0.6618

两个观察。错误写法的代价：把 $\theta_0$ 也惩罚，$\theta_0$ 从 2.5857 被压到 1.6546，全部预测值整体下沉 0.6618——模型整体偏离了数据的中心，这就是”对角线全写 1”的系统性偏移。正确写法的红利：$\theta_0$ 不受惩罚，它可以自由移动去吸收数据的均值——正确 $I'$ 的预测平均 2.6008，与 $y$ 的平均 2.6008 精确相等。偏置的职责就是整体平移，惩罚它没有任何好处。

实战口径：实战的最小二乘，一般都带着 Ridge

最后是讲师的实战经验之谈，值得原样记住：真正使用最小二乘法的时候，往往并不是真的只用前面那部分损失——一般都是加上这个 L2 惩罚项的。 换句话说，我们说的”最小二乘法”，实战里一般就是加了 Ridge Regression 的最小二乘。多花一个 $\lambda$ 的心思，换回可逆的矩阵、有界的解、不过分的系数——这笔账在实战里几乎总是划算的。

📝 动手练一练

1. 更新公式的两部分：岭回归的梯度下降更新 $\theta_j := (1-\frac{\alpha\lambda}{N})\theta_j - \alpha\cdot\frac{1}{N}\sum_i (h_\theta(x^{(i)})-y^{(i)})x_j^{(i)}$（$j \ge 1$）可以拆成两部分。请分别说明：这两部分各来自代价函数的哪一项？它们在优化过程中各自起什么作用？$\theta_0$ 的更新为什么没有第一部分（缩水因子）？

2. 手算一个 2×2 的岭回归：两个样本 $x = 0, 1$，$y = 0, 2$（添上 $x_0=1$ 后 $X = \begin{bmatrix}1 & 0\\ 1 & 1\end{bmatrix}$）。先验证普通最小二乘解是 $\theta = (0, 2)$（直线恰好穿过两点）；再取 $\lambda = 1$，手算 $(X^TX + \lambda I')\theta = X^Ty$ 的解（注意 $I'$ 左上角是 0）；最后算 $\lambda = 100$ 的解，观察 $\theta_0$、$\theta_1$ 各自奔向哪里。

本章小结

• 动机：$X^TX$ 求不出逆（特征数大于样本数、特征近似线性相关）时正规方程失守；解法是 regularization——在损失函数上加限制（英文 regulate：给更多规则和法则）。
• 岭回归损失：$J(\theta) = \frac{1}{2N}(\Vert X\theta-y\Vert _2^2 + \lambda\Vert \theta_{1:n}\Vert _2^2)$，惩罚只加 $\theta_1 \sim \theta_n$、不含偏置 $\theta_0$；平方恒正，逼着绝对值小。
• 两个视角：约束形式 $\Vert \theta_{1:n}\Vert _2^2 \le t$ 配拉格朗日乘子；$\lambda$ 是权衡两边的旋钮（$\lambda=0$ 过拟合、$\lambda\to\infty$ 欠拟合、中间有最低点）。
• 梯度下降：$\theta_j := (1-\frac{\alpha\lambda}{N})\theta_j - \alpha\frac{1}{N}\sum_i(h_\theta(x^{(i)})-y^{(i)})x_j^{(i)}$——先缩水、再走老路；$\alpha\lambda/N$ 通常只有 0.01 量级（$N$ 成千上万），每轮把 $\theta_j$ 乘上 0.99/0.98/0.97。
• 解析解：$\hat{\theta} = (X^TX+\lambda I')^{-1}X^Ty$，$I'$ 是单位矩阵把左上角换成 0（$\theta_0$ 不惩罚）；奇异矩阵加上 $\lambda I'$ 后重新可逆。
• 实战口径：实战使用的最小二乘，一般都带着 Ridge。

📋 行动清单

• 默写岭回归的损失函数与解析解，说出 $I'$ 左上角为什么是 0
• 从 $J(\theta)$ 出发推出梯度下降更新公式，整理出 $(1-\frac{\alpha\lambda}{N})$ 缩水因子
• 跑代码一：把两批 $y$ 的毛刺再放大一倍，看 $\lambda=0$ 的 $\theta_1$ 对冲是否跟着翻倍、$\lambda=2$ 时又如何
• 跑代码三：把 lambda 改成 10，观察错误 $I$ 的整体偏移是变大还是变小，并解释
• 预习下一节《LASSO》：把惩罚里的 L2 范数换成 L1 范数，会发生什么

—— 小象教研组