Matrix Calculus (2)

June 21, 2026

第 6 章 — 求根、优化与伴随微分

核心洞见

  • Newton 法每步大致使正确位数翻倍（二次收敛），但需要好的初值。
  • 大规模优化的关键：用反向/伴随模式；计算全部梯度的代价 ≈ 一次函数求值（与参数个数无关）。
  • 伴随法 = 把链式法则从输出往输入跑；只需要一次正向求解，反向传播+点积 即得 $g$ 和 $\nabla g$。
  • 参数很多时绝不要用正向模式或有限差分（每个参数一次求解），反向方式的优势是 一致的反向传播+不同的点积。而点积成本非常之低。

关键公式

  • Newton法求根（标量）：线性化 $f(x+\delta x)\approx f(x)+f'(x)\delta x=0\Rightarrow\delta x=-f(x)/f'(x)$，故 $x_{\text{new}} = x - \dfrac{f(x)}{f'(x)}$（二次收敛；$f'(x)=0$ 时失效）
  • Newton（多维，$f:\mathbb R^n\to\mathbb R^n$）：$x_{\text{new}} = x - f'(x)^{-1}f(x)$（每步用雅可比解一个线性系统）
  • 最速下降步：$x_{\text{new}}=x-\eta\,\nabla f$（$\eta$ = 步长 / "学习率"，由线搜索或信赖域确定）
  • 约束：最小化 $f(x)$ 满足 $g_k(x)\le0$（不等式）、$h_k(x)=0$（等式）；用 $\nabla f,\nabla g_k$ 朝最优可行点推进。（KKT条件+lagrange乘子，参考SVM求解过程）

伴随法 — $g(p)=f(x(p))$ 的完整推导（其中 $A(p)x=b$）

$$dg = f'(x)[dx] = f'(x)[d(A^{-1})b] = -\underbrace{f'(x)A^{-1}}_{v^T}\,dA\,\underbrace{A^{-1}b}_{x} = -v^T\,dA\,x$$

  • 按左到右分组：求解伴随（转置）方程 $A^Tv = f'(x)^T = \nabla_x f$。
  • 然后逐参数：$\dfrac{\partial g}{\partial p_k} = -v^T\dfrac{\partial A}{\partial p_k}x$（$\partial A/\partial p_k$ 往往很稀疏 → 廉价点积）。
  • 总共两次求解即得 $g$（解 $Ax=b$）和整个梯度 $\nabla g$（解 $A^Tv=\nabla_x f$），与参数个数无关。
  • 错误（正向）方式：$\partial g/\partial p_k = -f'(x)A^{-1}(\partial A/\partial p_k)x$ = 每个参数一次求解 — 参数多时避免。有限差分一样糟（每个参数一次求解）。

PS：详细解释 正向方式和反向方式的对比 核心矛盾：

  • 解一次方程组（如 $Ax=b$ 或 $A^Tv=f'(x)^T$）： 成本是 $O(N^3)$（直接法）或若干次矩阵-向量乘法（迭代法，如 GMRES）。
  • 显式计算 $A^{-1}$： 成本是 $O(N^3)$ 且会破坏 $A$ 的稀疏性，内存直接爆炸，实际工程中绝对禁止。
  • 我们再来看看正向和反向的实际工作量：
    • 正向模式 $$\frac{\partial g}{\partial p_k} = -f'(x) A^{-1} \left(\frac{\partial A}{\partial p_k}\right) x$$ 因为 $\frac{\partial A}{\partial p_k}$ 通常很稀疏，我们先算中间项 $y_k = \left(\frac{\partial A}{\partial p_k}\right) x$（这步很便宜）。然后，为了乘上 $A^{-1}$，你必须解一个方程组：$$A \cdot \underbrace{\frac{\partial x}{\partial p_k}}_{u_k} = -y_k$$注意！每一个参数 $p_k$，你都要解一个全新的线性方程组。
    • 反向模式 $$\frac{\partial g}{\partial p_k} = -\left(f'(x) A^{-1}\right) \left(\frac{\partial A}{\partial p_k} x\right) = -v^T \left(\frac{\partial A}{\partial p_k} x\right)$$这一次，我们通过定义 $v^T = f'(x)A^{-1}$，也就是解伴随方程：$$A^T v = \nabla_x f$$重点来了：这个伴随方程里，完全没有参数 $p_k$ 的影子！ 它只取决于当前的状态 $x$ 和目标函数 $f$。第一步： 解 $Ax=b$ 得到 $x$（第 1 次求解）。第二步： 解 $A^Tv = \nabla_x f$ 得到 $v$（第 2 次求解）。第三步： 针对每一个参数 $p_k$，计算 $-v^T \left(\frac{\partial A}{\partial p_k} x\right)$。因为 $\frac{\partial A}{\partial p_k}$ 极其稀疏，这只是一个廉价的点积（Dot Product），计算量几乎为 0。总计算量与参数个数无关。

隐函数定理（非线性约束 $h(p,x)=0\in\mathbb R^n$ 定义 $x(p)$）

$$\frac{\partial h}{\partial p}dp + \frac{\partial h}{\partial x}dx = 0 \Rightarrow dx = -\Big(\frac{\partial h}{\partial x}\Big)^{-1}\frac{\partial h}{\partial p}dp$$

$$dg = f'(x)\,dx \Rightarrow \text{伴随方程 } \Big(\frac{\partial h}{\partial x}\Big)^Tv = \nabla_x f,\quad \frac{\partial g}{\partial p_k} = -v^T\frac{\partial h}{\partial p_k}$$

仍是两次求解（一次非线性正向 + 一次线性伴随）。

第 7 章 — 行列式与逆的导数

核心洞见

  • 行列式的梯度是代数余子式矩阵 (cofactor matrix)；对数行列式的导数极其简洁，在应用数学中无处不在。
  • 逆的导数由 $A^{-1}A=I$ 一招得到；算子形式比显式雅可比更有用。

PS一些线性代数的基本知识 余子式 $M_{i,j} =$ 删掉第i行第j列的元素后形成的子矩阵的行列式 cofactor matrix(代数余子式矩阵) $C$ 伴随矩阵 adjugate matrix $\mathrm{adj}(A)=C^T$ $$ A \mathrm{adj}(A) = adj(A) A = \det(A) I $$

关键公式

  • 行列式梯度：$\nabla(\det A) = \mathrm{cofactor}(A) = (\det A)A^{-T} = \mathrm{adj}(A)^T$
    • 其中 $\mathrm{adj}(A)=\det(A)A^{-1}$（伴随/adjugate），余子式 $C_{ij}=(-1)^{i+j}\det(\text{minor}_{ij})$
  • 行列式微分：$d(\det A) = \det(A)\,\mathrm{tr}(A^{-1}dA) = \mathrm{tr}(\mathrm{adj}(A)\,dA) = \langle\,\mathrm{cofactor}(A), dA\,\rangle$
  • 两种推导：
    1. 沿第 $i$ 行的余子式展开 $\det A=\sum_j A_{ij}C_{ij}$，且 $C_{ij}$ 不含第 $i$ 行任何元素 ⇒ $\partial\det A/\partial A_{ij}=C_{ij}$，即 $\nabla\det A=C$。
    2. 在单位阵附近线性化：$\det(I+dA)-1=\mathrm{tr}(dA)$，则 $\det(A+A(A^{-1}dA))-\det A = \det A\,(\det(I+A^{-1}dA)-1)=\det A\,\mathrm{tr}(A^{-1}dA)$。
  • 对数行列式导数（很常用，如统计中）：$d(\log\det A) = \dfrac{d(\det A)}{\det A} = \mathrm{tr}(A^{-1}dA)$
  • 特征多项式 $p(x)=\det(xI-A)$：$dp = \det(xI-A)\,\mathrm{tr}((xI-A)^{-1})\,dx$（$\det M(x)=0$ 的根 → 特征值；一般 $M(x)$ → 非线性特征值问题）
  • 逆的导数：$d(A^{-1}) = -A^{-1}\,dA\,A^{-1}$；向量化雅可比：$\mathrm{vec}\,d(A^{-1}) = -(A^{-T}\otimes A^{-1})\mathrm{vec}(dA)$
  • 标量参数：$\dfrac{d(A^{-1})}{dt} = -A^{-1}\dfrac{dA}{dt}A^{-1}$

第 8 章 — 正向与反向模式自动微分

核心洞见

  • 自动微分 (AD) 不是符号微分（表达式会爆炸 — 多项式次数每步翻倍），也不是有限差分（AD 在精确算术下是精确的）。AD 更像编译器技术。
  • 正向模式 = 对偶数 (dual numbers)：把每个值扩展为（值，导数），重载运算符，导数就"魔法般"传播。
  • 计算图 + 路径乘积：导数 = 输入→输出所有路径上边权乘积之和。正向/反向只是遍历顺序不同。
  • 模式选择：正向成本 $\propto$ 输入数 $n$；反向成本 $\propto$ 输出数 $m$。
  • 二阶导用 forward-over-reverse，高效计算 Hessian–向量积。

关键公式

  • 对偶数 $a+b\epsilon$，$\epsilon^2=0$（类比复数的 $i^2=-1$）：
    • $(a+b\epsilon)\pm(c+d\epsilon)=(a\pm c)+(b\pm d)\epsilon$（加法法则）
    • $(a+b\epsilon)(c+d\epsilon) = ac + (bc+ad)\epsilon$（乘积法则）
    • $\dfrac{a+b\epsilon}{c+d\epsilon} = \dfrac{a}{c} + \dfrac{bc-ad}{c^2}\epsilon$（商法则）
  • 对偶数 = 线性化：$f(x+\epsilon) = f(x) + f'(x)\epsilon$（$\epsilon^2=0$ 规则恰好舍去高阶项）
  • 种子 (seed) 把输入的导数设为 1：输入 $D(x,1)=x+\epsilon$，读出输出的 $\epsilon$ 部分即为 $f'(x)$。

计算图 / 路径乘积

  • 给每条边 $A\to B$ 标上 $\partial B/\partial A$。导数 = 输入→输出所有路径上边权乘积之和。
  • 正向状态映射：$(\text{值},\,\text{路径积})\mapsto(f(\text{值}),\,f'\cdot\text{路径积})$，从 $(x,1)$ 开始。
  • 节点处合并多条入边：$\big(f(a,b),\ \tfrac{\partial z}{\partial a}p + \tfrac{\partial z}{\partial b}q\big)$。
  • 反向模式节点公式（对受 $a$ 影响的节点 $b_i$ 求和）：$\dfrac{\partial z}{\partial a} = \sum_i \dfrac{\partial b_i}{\partial a}\dfrac{\partial z}{\partial b_i}$，从输出 $(z,1)$ 开始。
  • 反向模式复用共享子积（如对源 $x,y$ 只算一次 $cd$ → 更少乘法）。

模式选择与二阶导

  • 正向成本 $\propto n$（输入数）；反向成本 $\propto m$（输出数）。 ML/优化：$m=1$ ⇒ 反向（反向传播）占主导。
  • 实践注意：正向更简单且与程序同序运行；反向需存储中间值的"录像带 (tape)"（内存大）。即使 $m=1$，在 $n$ 不算太大（~100）之前正向仍可能更快。
  • Forward-over-reverse 求二阶导：先用反向模式算 $\nabla f$，再用正向模式对它求导。
  • Hessian–向量积（避免构造 $n\times n$ Hessian）：$(\nabla f)'v = \left.\dfrac{\partial}{\partial\alpha}\nabla f|_{x+\alpha v}\right|_{\alpha=0}$ — 梯度的方向导数。由 Hessian 对称性，同样技巧给出 $v^T(\nabla f)'$。
  • 例：$h(x)=g(\nabla f|_x)$ ⇒ $\nabla h = \left.\dfrac{\partial}{\partial\alpha}\nabla f|_{x+\alpha\nabla g|_z}\right|_{\alpha=0}$，其中 $z=\nabla f|_x$ — 从不构造 Hessian。

速记表：值得背下来的公式

对象 微分 / 导数
$y = a^T x$ $dy = a^Tdx \quad \nabla _x y = a$
$y = x^Tx$ $dy = 2x^Tdx \quad \nabla _xy=2x$
$y = u^Tv$ $dy = v^Tdu + u^Tdv$
$y = uv^T$ $dy = (du)v^T + u(dv)^T$
$y = Ax$（对 $x$） $dy = A\,dx$, $J = \frac{\partial y}{\partial x} =A $
$y = x^TAx$ $dy=x^T(A+A^T)dx$; $\nabla = (A+A^T)x$; $H = A+A^T$
$y = A(x.\!*\!x)$ $J = 2A\,\mathrm{diag}(x)$
$y = AB$ $dy = dA\,B + A\,dB$
$y = A^2$ $dy = dA\,A + A\,dA$
$y = A^3$ $dy = dA\,A^2 + A\,dA\,A + A^2\,dA$
$y = A^{-1}$ $dy = -A^{-1}\,dA\,A^{-1}$
$y = tr(A^TX)$ $ dy = tr(A^TdX) \quad \nabla _X ​tr(A^TX)=A$
$y = \det A$ $dy = \det A\,\mathrm{tr}(A^{-1}dA)$; $\nabla = (\det A)A^{-T}$
$y = \log\det A$ $dy = \mathrm{tr}(A^{-1}dA)$
$y = \|X\|_F^2$ $\nabla_X y = 2X$
$y = \|AX - B\|_F^2$ $\nabla_X y = 2A^T(AX - B)$
$y = \|XA - B\|_F^2$ $\nabla_X y = 2(XA - B)A^T$
$y = x^TAy$（对 $A$） $\nabla = xy^T$

← Back to all posts