Matrix Calculus (1)

June 20, 2026

前言

起因是发现看dl框架的时候矩阵乘法的导数看不懂。

class MatMul(TensorOp):
  def gradient(self, out_grad, node):
        lhs,rhs=node.inputs
        lgrad=matmul(out_grad,transpose(rhs))
        rgrad=matmul(transpose(lhs),out_grad)
        return lgrad,rgrad

(引用来源：https://mp.weixin.qq.com/s/aRgYo2Bjr2k7Knl0r6kuFw) 再加上人智导期末也考了关于矩阵乘法的微分（比较基础，不涉及矩阵微积分），当时就感觉自己不是很熟悉。然而，矩阵微积分对于深度学习而言非常重要，所以自己进行了一番学习，以下是notes。(AI辅助生成，我的role是修订)

教材：Matrix Calculus (for Machine Learning and Beyond) — Alan Edelman & Steven G. Johnson (MIT 18.S096 / 18.063, arXiv:2501.14787)

核心思想

导数就是线性化 (linearization)。 对任意函数 $f$，其导数 $f'(x)$ 是一个线性算子，把输入的微小变化 $dx$ 映射到输出的微小变化： $$df = f(x+dx) - f(x) = f'(x)[dx] \quad (这里[\quad]代表输入的自变量是x)$$

第 1 章 — 概览与动机

核心洞见

  • 标量微积分"一切都简单"，但 $\frac{d}{dA}A^2 \neq 2A$（不可交换，而且不可定义除以矩阵），所以需要一套全新的微积分。
  • 把微分写成 $df = f'(x)\,dx$ — 不要除以 $dx$，因为 $dx$ 之后可能是向量/矩阵，无法做除法。

关键公式

  • 一元线性化：$ df = f'(x) dx$
  • 微分乘积法则（矩阵版，次序很重要 — 矩阵不可交换）：

$$d(AB) = (dA)\,B + A\,(dB)$$

  • 向量内积（点积可交换）：

$$d(x^Tx) = (dx)^Tx + x^T(dx) = (2x)^T dx \qquad d(u^Tv) = v^T du + u^T dv,\qquad d(uv^T) = (du)v^T + u(dv)^T$$

第 2 章 — 导数即线性算子

核心洞见

  • 导数的"形状"取决于输入/输出空间：
输入↓ \ 输出→ 标量 向量 矩阵
标量 标量 向量 矩阵
向量 梯度（列向量） 雅可比矩阵 (Jacobian) 高阶数组
矩阵 矩阵 高阶数组 高阶数组
  • 线性算子 $L$ 满足 $L[v_1+v_2]=Lv_1+Lv_2$ 且 $L[\alpha v]=\alpha L[v]$。注意 $f(x)=x+1$ 是仿射 (affine)，不是线性。有些线性算子（如 $L[A]=AB+CA$，或转置 $x\mapsto x^T$）不能写成简单的矩阵乘向量。
  • 方向导数 = 线性算子作用在某个方向上：$\left.\frac{\partial}{\partial\alpha}f(x+\alpha v)\right|_{\alpha=0} = f'(x)[v]$。这里 $v$ 可以是任意（非无穷小）向量。
  • 梯度 $\nabla f$ 是导数 $f'$ 的转置：$\nabla f = (f')^T$。梯度始终与 $x$ 同形状，所以 $df = \langle\nabla f, dx\rangle$。（注意这里的情况是输入x是vector，输出y是scalar）
  • 链式法则中，次序与结合性至关重要（线性算子不可交换）。

关键公式

  • 标量值函数 ($f:\mathbb R^n\to\mathbb R$)：$df = \nabla f \cdot dx = (\nabla f)^T dx$，其中 $f'(x)=(\nabla f)^T$ 是行向量（余向量/线性形式）。逐分量：$\nabla f=(\partial f/\partial x_1,\dots,\partial f/\partial x_n)^T$。
  • 向量值函数 ($f:\mathbb R^n\to\mathbb R^m$)：$df = J\,dx$，雅可比 $J_{ij} = \frac{\partial f_i}{\partial x_j}$（$m\times n$ 矩阵）
  • 线性映射是它自己的雅可比：$f(x)=Ax \Rightarrow df = A\,dx \Rightarrow f'(x)=A$
  • 加法法则：$f=g+h \Rightarrow f'=g'+h'$（线性算子构成向量空间 — 可以相加）
  • 乘积法则（一般情形，$g,h$ 可能不可交换）：$f=gh \Rightarrow df = (dg)\,h + g\,(dh)$
  • 二次型（这个很常用，但是理解不直观）：$f(x)=x^TAx \Rightarrow df = x^T(A+A^T)dx,\quad \nabla f = (A+A^T)x$（$A$ 对称时 $= 2Ax$）
  • 链式法则：$f=g(h(x)) \Rightarrow f'(x) = g'(h(x))\,h'(x)$（雅可比相乘；左到右 = 输出到输入）
  • 逐元素积 (Hadamard)：$d(x\,.\!*\,y) = dx\,.\!*\,y + x\,.\!*\,dy$；当 $y$ 为常量时，$g(x)=y.\!*\!x$ 的雅可比是 $\mathrm{diag}(y)$；而 $f(x)=A(x.\!*\!x)\Rightarrow J = 2A\,\mathrm{diag}(x)$

整体推导：$\nabla(x^TAx)$（例 10）

$$df=(x+dx)^TA(x+dx)-x^TAx = dx^TAx + x^TA\,dx + \underbrace{dx^TA\,dx}_{\text{舍去}}$$

计算成本（§2.5.1，全书最重要的洞见之一）

  • $m\times q$ 乘 $q\times p$ 矩阵 ≈ $2mpq$ 次运算，即 $\Theta(mpq)$。
  • 对 $f=a(b(c(x)))$，链式法则 $f'=a'b'c'$ 满足结合律：$(a'b')c' = a'(b'c')$。结果相同，但成本天差地别。
  • 反向模式 = 左到右（链式法则里面从左边开始计算）；正向模式 = 右到左。
  • 例 16：$n\gg1$ 个输入，$m=1$ 个输出，中间维度 $q=p=n$。反向模式成本 $\Theta(n^2)$；正向模式成本 $\Theta(n^3)$。（反向/反向传播胜出。）若角色互换（$m\gg1$, $n=1$），则正向胜出。
  • 结论：多输入/少输出 → 反向（这是神经网络里面常用的，反向传播）；多输出/少输入 → 正向。

多元导数（矩阵进 / 矩阵出）

  • $f(A)=A^2 \Rightarrow df = dA\,A + A\,dA$（展开 $(A+dA)^2-A^2$，舍去 $(dA)^2$）
  • $f(A)=A^3 \Rightarrow df = dA\,A^2 + A\,dA\,A + A^2\,dA$（不是 $3A^2\,dA$，除非 $A,dA$ 可交换）
  • $f(A)=A^{-1} \Rightarrow df = -A^{-1}\,dA\,A^{-1}$。推导：$d(AA^{-1})=dA\,A^{-1}+A\,d(A^{-1})=dI=0$。
  • 从刚才提到的“倒数作为线性算子”以及“方向导数角度理解”，对任意 $X$ 有 $f'(A)[X]=XA+AX$ — 写成"规则"容易，写成单个矩阵则困难。

第 3 章 — 矩阵函数的雅可比

核心洞见

  • 对矩阵→矩阵的函数，把 $f'(A)[dA]$ 写成"规则"（如 $dA\,A + A\,dA$）往往比显式雅可比更清晰。
  • 若硬要写成雅可比矩阵，需用向量化 (vec) + Kronecker 积。
  • 成本问题：用 Kronecker 积显式计算通常会让成本从 $\sim m^3$ 爆炸到 $\sim m^4$（存储与运算）。把它当作概念工具，或用于稀疏矩阵。

关键公式

  • vec：把 $m\times n$ 矩阵按列堆叠成 $mn$ 维列向量（列主序；Julia/Fortran/Matlab 默认）。$\mathrm{vec}\,A$ = $A$ 在矩阵的"笛卡尔"基下的坐标。
  • Kronecker 积：对 $A$（$m\times n$）和 $B$（$p\times q$），$A\otimes B$ 是 $mp\times nq$ 矩阵，把 $B$ 的副本按 $a_{ij}$ 缩放后"贴入"：

$$A\otimes B = \begin{pmatrix} a_{11}B & \cdots & a_{1n}B\\ \vdots & & \vdots\\ a_{m1}B & \cdots & a_{mn}B\end{pmatrix},\qquad A\otimes B \neq B\otimes A \text{（一般情况下）}$$

  • 矩阵平方：$f(A)=A^2 \Rightarrow f'(A)[dA] = dA\,A + A\,dA$。向量化后的 $m^2\times m^2$ 雅可比为 $\tilde f' = I\otimes A + A^T\otimes I$。
  • Kronecker 关键恒等式（命题 27） — 连接算子形式与向量化形式的桥梁：

$$(A\otimes B)\,\mathrm{vec}(C) = \mathrm{vec}(BCA^T)$$

（口诀：$B$ 在 $C$ 左边，$A$ 绕到右边并转置。）两个构成它的特例：

$$(I\otimes B)\,\mathrm{vec}\,C = \mathrm{vec}(BC),\qquad (A\otimes I)\,\mathrm{vec}\,C = \mathrm{vec}(CA^T)$$

  • $A^3$ 的向量化雅可比：$\tilde f' = (A^2)^T\otimes I + A^T\otimes A + I\otimes A^2$（对 $dA\,A^2 + A\,dA\,A + A^2\,dA$ 三次应用命题 27）。
  • Kronecker 代数（很好的练习）：
    • $(A\otimes B)^T = A^T\otimes B^T$
    • $(A\otimes B)(C\otimes D) = (AC)\otimes(BD)$
    • $(A\otimes B)^{-1} = A^{-1}\otimes B^{-1}$；若 $A,B$ 正交则 $A\otimes B$ 正交
    • $\det(A\otimes B) = (\det A)^m(\det B)^n$（$A\in\mathbb R^{n\times n},B\in\mathbb R^{m\times m}$）；$\mathrm{tr}(A\otimes B) = \mathrm{tr}(A)\,\mathrm{tr}(B)$
    • 特征值：$Au=\lambda u,\ Bv=\mu v \Rightarrow (A\otimes B)(u\otimes v)=\lambda\mu\,(u\otimes v)$
  • 2×2 雅可比示例：$A^2$，$A=\left(\begin{smallmatrix}p&r\\q&s\end{smallmatrix}\right)$（行/列按矩阵元素列主序编号）：

$$\tilde f' = I\otimes A + A^T\otimes I = \begin{pmatrix} 2p & r & q & 0\\ q & p+s & 0 & q\\ r & 0 & p+s & r\\ 0 & r & q & 2s\end{pmatrix}$$

  • Sylvester 方程 $AX+XB=C$ 向量化：$(I\otimes A + B^T\otimes I)\mathrm{vec}\,X = \mathrm{vec}\,C$（朴素高斯消元 $\sim (m^2)^3 = m^6$；专门算法 $\sim m^3$）。

成本警告（量化）

  • 直接计算 $C\mapsto BCA^T$：两次 $m\times m$ 乘法 $\sim m^3$ 运算，$\sim m^2$ 存储。
  • 通过 $(A\otimes B)\mathrm{vec}\,C$：构造 $A\otimes B$ 需 $\sim m^4$ 运算和 $\sim m^4$ 存储 — 对 $m=1000$ 是 TB 级 vs MB 级。Kronecker 积应作为概念工具使用（或用于稀疏矩阵，此时 $A\otimes B$ 仍稀疏）。

第 4 章 — 有限差分近似

核心洞见

  • 手算导数极易出错 — 永远要用有限差分做合理性检验（这才是它的主要用途，而非真正用来计算导数）。
  • 截断误差（$\delta x$ 太大）与舍入误差（$\delta x$ 太小 → 灾难性抵消）之间的拉锯。
  • 高维下完整雅可比需要 $n$ 次差分，故对大问题不实用。

关键公式

  • 底层展开：$\delta f = f(x+\delta x)-f(x) = f'(x)\,\delta x + o(\|\delta x\|)$
  • 前向差分（一阶）：$f'(x)\,\delta x \approx f(x+\delta x) - f(x)$；后向：$f'(x)\,\delta x\approx f(x)-f(x-\delta x)$
  • 标量形式（仅当 $x$ 为标量）：$f'(x)\approx \dfrac{f(x+\delta x)-f(x)}{\delta x}$
  • 相对误差：$\dfrac{\|\text{近似}-\text{精确}\|}{\|\text{精确}\|}$（相对于正确答案而言要小）
  • 中心差分（二阶）：$f'(x)\,\delta x \approx \dfrac{f(x+\delta x)-f(x-\delta x)}{2}$（误差 $\propto\|\delta x\|^2$ 。理解：使用tayler展开发现可以消掉h一次方项 $\frac{f(x+h) - f(x-h)}{2h} = f'(x) + \frac{h^2}{6}f'''(x) + \dots$）
  • Frobenius 范数：$\|A\| = \sqrt{\sum_{ij}|A_{ij}|^2} = \sqrt{\mathrm{tr}(A^TA)}$
  • 机器精度 $\epsilon = 2^{-52}\approx 2.22\times10^{-16}$；存储数的舍入：$|\tilde y - y|\le\epsilon|y|$（约 15–16 位十进制有效数字）
  • 精度阶：误差 $\propto\|\delta x\|^n$ ⇒ "$n$ 阶精度"；前向 = 1 阶，中心 = 2 阶。
  • 经验法则：最佳步长 $\|\delta x\| \approx \sqrt{\epsilon}\,\|x\| \approx 10^{-8}\|x\|$（保留约一半有效数字）

实例检验：$f(A)=A^2$

  • 精确：$df = A\,dA + dA\,A$。取随机 $A$ 和 $dA = 10^{-8}A$。
  • $\|f(A+dA)-f(A) - (A\,dA+dA\,A)\|/\|\cdots\| \approx 10^{-8}$（✓ 吻合），而与错误的 $2A\,dA$ 比较给出相对误差 $\approx 10^0$（✗）— 立刻抓住 bug。这就是经典的"永远用有限差分检查你的导数"工作流。

为什么太小也不好：若 $\delta x$ 极小，$f(x+\delta x)-f(x)$ 会因灾难性抵消（计算机计算的精度限制）丢失全部有效数字，相对误差反而升高。整体呈 U 形误差曲线，在 $\sqrt\epsilon$ 附近最小。

第 5 章 — 一般向量空间中的导数

核心洞见

  • 把"点积"推广为内积、把"长度"推广为范数，导数即可推广到矩阵、函数等（Fréchet 导数）。
  • 梯度依赖于内积的选择 — 换内积就换梯度。
  • 对矩阵用 Frobenius 内积，得到 $(\nabla f)_{ij} = \partial f/\partial A_{ij}$，但整体求导胜过逐元素。

关键公式

  • 内积公理：(1) 对称 $\langle x,y\rangle=\langle y,x\rangle$，(2) 线性 $\langle x,\alpha y+\beta z\rangle=\alpha\langle x,y\rangle+\beta\langle x,z\rangle$，(3) 正定 $\langle x,x\rangle=\|x\|^2\ge0$（$=0$ 当且仅当 $x=0$）；Cauchy–Schwarz：$|\langle x,y\rangle|\le\|x\|\,\|y\|$
  • 标准欧氏内积：$\langle x,y\rangle = x^Ty$
  • 加权内积：$\langle x,y\rangle_W = x^TWy$（$W$ 对称正定）。换内积就换梯度：例如对 $f=x^TAx$，欧氏梯度是 $(A+A^T)x$，但 $\nabla^{(W)}f = W^{-1}(A+A^T)x$。
  • Frobenius 内积（矩阵默认）：$\langle A,B\rangle_F = \sum_{ij}A_{ij}B_{ij} = \mathrm{tr}(A^TB) = \mathrm{vec}(A)^T\mathrm{vec}(B)$，范数 $\|A\|_F=\sqrt{\mathrm{tr}(A^TA)}$
  • Riesz 表示定理：任意线性形式都等于与某向量的内积 — 那个向量就是梯度：$f'(x)[dx] = \langle\nabla f, dx\rangle = df$。所以 $\nabla f$ 始终与 $x$ 同形状。
  • 一般矩阵梯度（经 Frobenius）：$df = \langle\nabla f, dA\rangle \Rightarrow (\nabla f)_{ij}=\partial f/\partial A_{ij}$
  • 范数梯度：$\nabla\|A\|_F = \dfrac{A}{\|A\|_F}$；类似地 $\nabla\|x\| = x/\|x\|$
  • $f(A) = x^TAy \Rightarrow \nabla f = xy^T$（迹技巧：$df = x^T dA\,y = \mathrm{tr}(yx^T dA) = \langle xy^T, dA\rangle$，矩阵trace的性质-循环移位不变性）
  • 范数公理：非负、齐次 $\|\alpha v\|=|\alpha|\|v\|$、三角不等式 $\|u+v\|\le\|u\|+\|v\|$
  • Hilbert 空间（完备 + 内积）⊂ Banach 空间（完备 + 范数）；$o(\delta x)$ 需要范数：$\lim_{\delta x\to0}\frac{\|o(\delta x)\|}{\|\delta x\|}=0$。Banach 空间上的导数 = Fréchet 导数。

整体推导：$\nabla\|A\|_F$（例 35）

$$df = \tfrac{1}{2\|A\|_F}\,d(\mathrm{tr}\,A^TA) = \tfrac{1}{2\|A\|_F}\big(\mathrm{tr}(dA^TA)+\mathrm{tr}(A^TdA)\big) = \tfrac{1}{\|A\|_F}\mathrm{tr}(A^TdA) = \Big\langle \tfrac{A}{\|A\|_F}, dA\Big\rangle$$

其中用到 $\mathrm{tr}(B)=\mathrm{tr}(B^T)$ 和 $d(\mathrm{tr}\,B)=\mathrm{tr}(dB)$。

← Back to all posts