在机器学习中，我们有时会将损失函数表示为向量或者是矩阵的函数。由于 “learning as optimization” 的思想，对损失函数求导是一种常见的操作，而梯度下降、Newton-Raphson Methods 也要求我们计算损失函数的梯度甚至是 Hessian 矩阵。在统计中，一些多元分布，例如多元正态分布，其密度函数写成矩阵的形式更加简洁美观，而在计算参数的极大似然估计时，更是考验我们对向量或矩阵求导的功底。本文将简单介绍矩阵求导的定义与方法，并演示一些算例。
矩阵求导定义与布局
我们已经学过标量对标量的求导，比如标量 $y$ 对标量 $x$ 的求导，写成：$\frac{\partial y}{\partial x}$。
有时候，我们需要对一组标量 $y_{i}, i=1,2, \dots, m$ 对标量 $x$ 求导，此时得到一组标量求导结果：
\[\frac{\partial y_{i}}{\partial x}, i=1,2, \dots, m\]如果将这组标量 $y_i$ 表示成 $m$ 维向量 $\pmb y$ ，则求导的结果也可以用 $m$ 维向量表示，记作：$\frac{\partial \pmb{y}}{\partial x}$。
可见，向量对标量的求导，其实就是将向量中的每个元素对标量求导，再将求导得到的结果表示成向量。类似的结论可以应用到矩阵对标量的求导，标量对向量的求导，向量对向量的求导等等。
为了便于描述，若没有特别说明，则求导的自变量用 $x$ 表示标量，$\pmb x$ 表示 $n$ 维向量，$X$ 表示 $m \times n$ 维矩阵，求导的因变量用 $y$ 表示标量，用 $\pmb y$ 表示 $m$ 维向量，用 $Y$ 表示 $p \times q$ 维矩阵。
分别考虑自变量和因变量可能为标量，向量或是矩阵的情况，我们得到 9 种结果：
| 自变量\因变量 | 标量 $y$ | 向量 $\pmb y$ | 矩阵 $Y$ |
|---|---|---|---|
| 标量 $x$ | $\frac{\partial y}{\partial x}$ | $\frac{\partial \pmb y}{\partial x}$ | $\frac{\partial Y}{\partial x}$ |
| 向量 $\pmb x$ | $\frac{\partial y}{\partial \pmb x}$ | $\frac{\partial \pmb y}{\partial \pmb x}$ | $\frac{\partial Y}{\partial \pmb x}$ |
| 矩阵 $X$ | $\frac{\partial y}{\partial X}$ | $\frac{\partial \pmb y}{\partial Z}$ | $\frac{\partial Y}{\partial X}$ |
定义矩阵求导后，我们会自然而然产生一个问题。按照上述的想法，$\frac{\partial y}{\partial X}$ 可以表示成一个 $m \times n$ 的矩阵，而 $\frac{\partial Y}{\partial x}$ 可以表示成一个 $p \times q$ 的矩阵，那么 $\frac{\partial \pmb y}{\partial \pmb x}$ 应该表示成 $m \times n$ 的矩阵还是 $n \times m$ 的矩阵？
答案是：Both are OK !
因为我们关心的是这 $mn$ 个求导的结果，而不太关心这结果的形式。同理对于 $\frac{\partial Y}{\partial X}$ 来说，理论上应该表示成 $m \times n \times p \times q$ 维的张量，但事实上我们只需用一个可以囊括这 $mnpq$ 个数的矩阵来表示就可以了，这是后话。
但为了保持写法的一致，我们有必要引入求导布局的概念。最基本的求导布局有两个：分子布局(numerator layout)和分母布局(denominator layout )。
对于分子布局来说，我们求导结果的维度以分子为主，比如上面向量对标量求导的例子，结果的维度和分子的维度是一致的。即因变量 $\pmb y$ 若是 $m$ 维列向量，则求导结果 \(\frac{\partial \pmb y}{\partial x}\) 也为 $m$ 维列向量 。
同理，对于分母布局来说，我们求导结果的维度以分母为主，同样以向量对标量求导为例，因变量 $\pmb y$ 若是 $m$ 维列向量，则求导结果应为 $m$ 维行向量。
可见，对于分子布局和分母布局来说，两者相差一个转置。
但是在机器学习算法原理的资料推导里，我们通常不会看到说正在使用什么布局，这需要我们自己去判断使用的是哪种布局。一般来说我们会使用一种混合布局的思路，即分子和分母谁的维度大（这里的维度指的是标量、向量和矩阵），就用谁的布局，这大概就是降维打击吧。
但对于向量对向量求导，则有些分歧。若是采用分子布局，则结果是一个 $m \times n$ 的矩阵：
\[\frac{\partial \pmb{y}}{\partial \pmb{x}}=\left(\begin{array}{cccc}{\frac{\partial y_{1}}{\partial x_{1}}} & {\frac{\partial y_{1}}{\partial x_{2}}} & {\cdots} & {\frac{\partial y_{1}}{\partial x_{n}}} \\ {\frac{\partial y_{2}}{\partial x_{1}}} & {\frac{\partial y_{2}}{\partial x_{2}}} & {\cdots} & {\frac{\partial y_{2}}{\partial x_{n}}} \\ {\vdots} & {\vdots} & {\ddots} & {\vdots} \\ {\frac{\partial y_{m}}{\partial x_{1}}} & {\frac{\partial y_{m}}{\partial x_{2}}} & {\cdots} & {\frac{\partial y_{m}}{\partial x_{n}}}\end{array}\right)\]上边这个按分子布局的向量对向量求导的结果矩阵，我们一般叫做雅克比 (Jacobian)矩阵。有的资料上会使用 $\frac{\partial \pmb{y}}{\partial \pmb{x}^T}$ 来定义雅克比矩阵，意义是一样的。
如果是按分母布局，则求导的结果矩阵的第一维度会以分母为准，即结果是一个 $n \times m$ 的矩阵：
\[\frac{\partial \pmb{y}}{\partial \pmb{x}}=\left(\begin{array}{cccc}{\frac{\partial y_{1}}{\partial x_{1}}} & {\frac{\partial y_{2}}{\partial x_{1}}} & {\ldots} & {\frac{\partial y_{m}}{\partial x_{1}}} \\ {\frac{\partial y_{1}}{\partial x_{2}}} & {\frac{\partial y_{2}}{\partial x_{2}}} & {\cdots} & {\frac{\partial y_{m}}{\partial x_{2}}} \\ {\vdots} & {\vdots} & {\ddots} & {\vdots} \\ {\frac{\partial y_{1}}{\partial x_{n}}} & {\frac{\partial y_{2}}{\partial x_{n}}} & {\cdots} & {\frac{\partial y_{m}}{\partial x_{n}}}\end{array}\right)\]上边这个按分母布局的向量对向量求导的结果矩阵，我们一般叫做梯度矩阵。
我们暂时使用分子布局作为向量对向量求导的布局，则所用矩阵求导布局总结如下：
| 自变量\因变量 | 标量 $y$ | 向量 $\pmb y$ | 矩阵 $Y$ |
|---|---|---|---|
| 标量 $x$ | / | $\frac{\partial \pmb y}{\partial x}$ 默认布局：分子布局 | $\frac{\partial Y}{\partial x}$默认布局：分子布局 |
| 向量 $\pmb x$ | $\frac{\partial y}{\partial \pmb x}$默认布局：分母布局 | $\frac{\partial \pmb y}{\partial \pmb x}$默认布局：分子布局 | / |
| 矩阵 $X$ | $\frac{\partial y}{\partial X}$默认布局：分母布局 | / | / |
在下文中，我们主要讨论标量对向量或矩阵求导，这是因为机器学习中的目标函数值通常是标量，所以下文以分母布局为主。
矩阵求导定义法：
在我们没有学过矩阵求导之前，遇到标量对向量或是矩阵求导，一个自然而然的想法就是：将这个关于向量或是矩阵的函数分解成其各元素的函数，再分别对向量或矩阵的元素进行求导，最后找到规律，得到求导结果的向量或矩阵。
首先来看一个标量对向量求导的例子：
例1：$y = \pmb a ^T \pmb x$， 求 $\frac{\partial y}{\partial \pmb x}$。
解：首先将 $y$ 写成关于 $x_j$ 的函数，即 $y = \sum_{j=1}^{n} a_{j} x_{j}$， 则
\[\frac{\partial y}{\partial x_i} =\frac{\partial \sum_{j=1}^{n} a_{j} x_{j}}{\partial x_i} = \frac{\partial a_{i} x_{i}}{\partial x_{i}} = a_{i}\]可见，对向量的第 $i$ 个变量求导的结果就是分量 $a_i$ ，由于是分母布局，所以对向量 $\pmb x$ 求导的结果就是 $\pmb a$。
\[\frac{\partial \pmb{a}^{T} \pmb{x}}{\partial \pmb{x}}=\pmb{a}\]再看一个标量对向量求导的例子：
例2：$y=\pmb{a}^{T} {X} \pmb{b}$，求 $\frac{\partial y}{\partial X}$，其中 $\pmb a$ 是 $m$ 维向量，$\pmb b$ 是 $n$ 维向量。
解：首先将 $y$ 写成各元素的函数，即：$y = \sum_{p=1}^{m} \sum_{q=1}^{n} a_{p} A_{p q} b_{q}$，则
\[\frac{\partial y}{\partial X_{i j}}=\frac{\partial \sum_{p=1}^{m} \sum_{q=1}^{n} a_{p} A_{p q} b_{q}}{\partial X_{i j}}=\frac{\partial a_{i} A_{i j} b_{j}}{\partial X_{i j}}=a_{i} b_{j}\]对 $X_{ij}$ 求导的结果为 $\pmb a$ 的第 $i$ 个元素与 $\pmb b$ 的第 $j$ 个元素的乘积，将所有结果排列成 $m \times n$ 的矩阵，这个矩阵即为 $\pmb a \pmb b^T$，故最后求导结果为：
\[\frac{\partial y}{\partial X} =\pmb a \pmb b^T\]简单的求导的确不难，但许多用矩阵表达的函数并不容易拆分成各元素的函数，例如 $y = \lvert X\lvert$，或是用矩阵向量表达有更简洁的形式，此时再用定义法求解导数过程会很复杂，我们需要使用新的方法来求解矩阵向量求导。
矩阵求导微分法
矩阵微分
在高数中我们学习过标量导数与微分的关系，给定函数 $f = f(x) $，则有 $d f=f^{\prime}(x) d x$。
如果自变量是向量 $\pmb x$ ，则有
\[d f=\sum_{i=1}^{n} \frac{\partial f}{\partial x_{i}} d x_{i}=\left(\frac{\partial f}{\partial \pmb{x}}\right)^{T} d \pmb {x}\]我们可以发现标量对向量的求导和它的向量微分有一个转置的关系。
如果自变量是矩阵 $X$ ，$f$ 的微分可以写成：
\[d f=\sum_{i=1}^{m} \sum_{j=1}^{n} \frac{\partial f}{\partial X_{i j}} d X_{i j}=\operatorname{tr}\left(\left(\frac{\partial f}{\partial {X}}\right)^{T} d {X}\right)\]上式用到了迹函数的性质
\[\operatorname{tr}\left(A^{T} B\right) = \sum_{j=1}^n (A^T B)_{jj}= \sum_{j=1}^n\sum_{i=1}^m A_{ij}B_{ij}\]我们可以看出矩阵导数和微分也有一个转置的关系，只不过外面套了个迹函数。由于标量的迹等于自身，因此我们可以将矩阵和向量的微分统一表示成：
\[d f=\operatorname{tr}\left(\left(\frac{\partial f}{\partial \mathbf{X}}\right)^{T} d \mathbf{X}\right) \quad d f=\operatorname{tr}\left(\left(\frac{\partial f}{\partial \mathbf{x}}\right)^{T} d \mathbf{x}\right)\]于是我们得到了微分法求解矩阵求导的基本思想：给定函数 $f = f(X)$ ，若我们能将 $f$ 写成微分形式 $df = tr(AdX)$，则 $\frac{\partial f}{\partial X} =A^T$。
基于此思想，使用微分法求解矩阵求导前，我们要先介绍矩阵微分的运算法则和一些迹技巧。
矩阵微分的运算法则
对于 $m \times n$ 维的变量 $X$ ，$d X$ 仍然是 $m \times n$ 维的矩阵。
\[dX=\left(\begin{array}{cccc}{dX_{11}} & {dX_{12}} & {\ldots} & {dX_{1n}} \\ {dX_{21}} & {dX_{22}} & {\cdots} & {dX_{2n}} \\ {\vdots} & {\vdots} & {\ddots} & {\vdots} \\ {dX_{m1}} & {dX_{m2}} & {\cdots} & {dX_{mn}}\end{array}\right)\]对于 $m \times n $ 维的常量 $A$，$d A$ 则为 $m \times n$ 的零矩阵。
常见微分的运算法则如下：
- 微分加减法：$d(X+Y)=d X+d Y, d(X-Y)=d X-d Y$
- 微分乘法：$d(X Y)=(d X) Y+X(d Y)$
- 微分转置：$d\left(X^{T}\right)=(d X)^{T}$
- 微分的迹：$d t r(X)=t r(d X)$
- 微分哈达马乘积：$d(X \odot Y)=X \odot d Y+d X \odot Y$
- 逐元素求导：$d \sigma(X)=\sigma^{\prime}(X) \odot d X$
- 逆矩阵微分：$d X^{-1}=-X^{-1} d X X^{-1}$。此式可在 $X X^{-1}=I$ 两侧求微分得到。
- 行列式微分：$d\lvert X\lvert=\operatorname{tr}\left(X^{\ast} d X\right)$，其中 $X^{\ast}$ 是 $X$ 的伴随矩阵。当 $X$ 可逆时，$d\lvert X \lvert=\lvert X\lvert \operatorname{tr}\left(X^{-1} d X\right)$。
迹技巧
我们试图利用矩阵导数与微分的联系 $d f=\operatorname{tr}\left(\frac{\partial f^{T}}{\partial X} d X\right)$，在求出左侧的微分 $df$ 后，该如何写成右边的形式呢？这里还需要用到一些迹技巧：
- 标量的迹等于自己：$a=\operatorname{tr}(a)$
- 转置：$\operatorname{tr}\left(A^{T}\right)=\operatorname{tr}(A)$
- 线性：$\operatorname{tr}(A \pm B)=\operatorname{tr}(A) \pm \operatorname{tr}(B)$
- 矩阵乘法交换：$\operatorname{tr}(AB)=\operatorname{tr}(BA) $ ，其中 $A^T$$ 与 $ $B$ 尺寸相同，两侧都等于 $\sum_{i, j} A_{i j} B_{i j}$。
- 矩阵乘法/逐元素乘法交换：$\operatorname{tr}\left(A^{T}(B \odot C)\right)=\operatorname{tr}\left((A \odot B)^{T} C\right)$，其中 $A,B,C $ 尺寸相同，两侧都等于 $\sum_{i, j} A_{i j} B_{i j} C_{i j}$
有了这些性质，我们再来看看如何用微分法进行矩阵求导。
使用微分法求解矩阵向量求导
再介绍了矩阵微分与导数的联系以及上述性质后，我们可以总结出矩阵求导的基本步骤：
若标量函数 $f$ 是矩阵 $X$ 经加减乘法、逆、行列式、逐元素函数等运算构成，则使用相应的运算法则对 $f$ 求微分，给 $df$ 套上迹，再使用迹技巧将其它项交换至 $dX$ 左侧，对照导数与微分的联系 $d f=\operatorname{tr}\left(\frac{\partial f^{T}}{\partial X} d X\right)$， 即能得到导数。
特别地，若矩阵退化成向量，对照导数与微分的联系 $d f=\frac{\partial f}{\partial \boldsymbol{x}} ^{T}d \boldsymbol{x}$，即能得到导数。
接下来我们仍然以例2为例，使用微分法求解矩阵导数：
例2.2：$y=\pmb{a}^{T} {X} \pmb{b}$，求 $\frac{\partial y}{\partial X}$，其中 $\pmb a$ 是 $m$ 维向量，$\pmb b$ 是 $n$ 维向量。
解：首先使用微分的性质对 $y$ 进行微分：
\[d y=d \pmb{a}^{T} {X} \pmb{b}+\pmb{a}^{T} d {X} \pmb{b}+\pmb{a}^{T} {X} d \pmb{b}=\pmb{a}^{T} d {X} \pmb{b}\]第二步，在两边套上迹函数，并利用迹技巧将 $d X$ 移到最右侧：
\[d y=\operatorname{tr}(d y)=\operatorname{tr}\left(\mathbf{a}^{T} d \mathbf{X} \mathbf{b}\right)=\operatorname{tr}\left(\mathbf{b a}^{T} d \mathbf{X}\right)\]通过对比 $d f=\operatorname{tr}\left(\frac{\partial f^{T}}{\partial X} d X\right)$ 我们得知矩阵导数就是 $d X$ 左边部分的转置：
\[\frac{\partial f}{\partial {X}}=\left(\pmb{b} \pmb{a}^{T}\right)^{T}=\pmb a \pmb b^{T}\]这就是微分法的基本流程，比起用定义法求矩阵导数是不是要方便得多呢？
为了熟练使用这种求导方法，接下来将提供一些算例以供参考。
一些算例
例3：$f=\boldsymbol{a}^{T} \exp (X \boldsymbol{b})$，求 $\frac{\partial f}{\partial X}$。其中，$\pmb a$ 是 $m$ 维向量，$\pmb b$ 是 $n$ 维向量，$X$ 是 $m \times n$ 的矩阵。
解：使用矩阵乘法、逐元素求导进行微分：
\[\begin{align} df &= \pmb a^T d\left(\exp \left(X \boldsymbol{b}\right)\right) \\ &=\boldsymbol{a}^{T}(\exp (X \boldsymbol{b}) \odot(d X \boldsymbol{b})) \end{align}\]套上迹并作变换：
\[\begin{align} df &= \operatorname{tr}\left(\boldsymbol{a}^{T}(\exp (X \boldsymbol{b}) \odot(d X \boldsymbol{b}))\right) \\ & = \operatorname{tr}\left((\boldsymbol{a} \odot \exp (X \boldsymbol{b}))^{T} d X \boldsymbol{b}\right) \end{align}\]这里用到了公式：$\operatorname{tr}\left(A^{T}(B \odot C)\right)=\operatorname{tr}\left((A \odot B)^{T} C\right)$ ，交换了 $\pmb a, \exp (X \boldsymbol{b})$ 和 $d X \pmb b$ 。
\[\begin{align} df &= \operatorname{tr}\left(\boldsymbol{b}(\boldsymbol{a} \odot \exp (X \boldsymbol{b}))^{T} d X\right) \\ & = \operatorname{tr}\left(\left((\boldsymbol{a} \odot \exp (X \boldsymbol{b})) \boldsymbol{b}^{T}\right)^{T} d X\right) \end{align}\]则：
\[\frac{\partial f}{\partial X}=(\boldsymbol{a} \odot \exp (X \boldsymbol{b})) \boldsymbol{b}^{T}\]例4：最小二乘。$l = \Vert \pmb y - X \pmb \beta \Vert ^2$，求 $\pmb \beta$ 的最小二乘估计，即求 $\frac{\partial l}{\partial \pmb{\beta}}$ 的零点。
解：这是标量对向量的求导，可以把向量看作是矩阵的特例。向量模的平方可以改写成向量与自身的内积，即：
\[l = (\pmb y - X \pmb \beta)^T(\pmb y - X \pmb \beta)\]再求微分：
\[\begin{align} dl &= d\left((\pmb y - X \pmb \beta)^T\right)(\pmb y - X \pmb \beta) + (\pmb y - X \pmb \beta)^Td(\pmb y - X \pmb \beta) \\ & = (- X d\pmb \beta)^T(\pmb y - X \pmb \beta) - (\pmb y - X \pmb \beta)^TX d\pmb \beta \\ &=2(X\pmb \beta - \pmb y)^TXd\pmb\beta \end{align}\]对照导数与微分的联系，得到：
\[\frac{\partial l}{\partial \pmb \beta} = 2 X^T(X\pmb \beta - \pmb y)\]令导数为 0 则有：
\[2 X^T(X\hat {\pmb \beta} - \pmb y) = \pmb 0\\ X^TX\hat {\pmb \beta} = X^T\pmb y \\ \hat {\pmb \beta} = (X^TX)^{-1} X^T \pmb y\]例5：多元正态分布方差极大似然估计。样本 $\pmb x_1, \dots, \pmb x_N \sim \mathcal N(\pmb \mu , \Sigma)$ ，求 $\Sigma$ 的极大似然估计。
解：由于多元正态密度函数为：
\[p(\pmb x ; \pmb \mu, \Sigma)=\frac{1}{(2 \pi)^{n / 2}\lvert\Sigma\lvert^{1 / 2}} \exp \left(-\frac{1}{2}(\pmb x-\pmb\mu)^{T} \Sigma^{-1}(\pmb x-\pmb \mu)\right)\]最大化对数似然函数等价于最小化：
\[l=\log \lvert\Sigma\lvert+\frac{1}{N} \sum_{i=1}^{N}\left(\boldsymbol{x}_{i}-\overline{\boldsymbol{x}}\right)^{T} \Sigma^{-1}\left(\boldsymbol{x}_{i}-\overline{\boldsymbol{x}}\right)\]首先求微分，第一项是
\[d \log \lvert\Sigma\lvert=\lvert\Sigma\lvert^{-1} d\lvert\Sigma\lvert=\operatorname{tr}\left(\Sigma^{-1} d \Sigma\right)\]第二项是
\[\frac{1}{N} \sum_{i=1}^{N}\left(\boldsymbol{x}_{i}-\overline{\boldsymbol{x}}\right)^{T} d \Sigma^{-1}\left(\boldsymbol{x}_{i}-\overline{\boldsymbol{x}}\right)=-\frac{1}{N} \sum_{i=1}^{N}\left(\boldsymbol{x}_{i}-\overline{\boldsymbol{x}}\right)^{T} \Sigma^{-1} d \Sigma \Sigma^{-1}\left(\boldsymbol{x}_{i}-\overline{\boldsymbol{x}}\right)\]对第二项套上迹并作变换
\[\begin {align} &\operatorname{tr}\left(\frac{1}{N} \sum_{i=1}^{N}\left(\boldsymbol{x}_{i}-\overline{\boldsymbol{x}}\right)^{T} \Sigma^{-1} d \Sigma \Sigma^{-1}\left(\boldsymbol{x}_{i}-\overline{\boldsymbol{x}}\right)\right) \\ =&\frac{1}{N} \sum_{i=1}^{N} \operatorname{tr}\left(\left(\boldsymbol{x}_{i}-\overline{\boldsymbol{x}}\right)^{T} \Sigma^{-1} d \Sigma \Sigma^{-1}\left(\boldsymbol{x}_{i}-\overline{\boldsymbol{x}}\right)\right)\\ =&\frac{1}{N} \sum_{i=1}^{N} \operatorname{tr}\left(\Sigma^{-1}\left(\boldsymbol{x}_{i}-\overline{\boldsymbol{x}}\right)\left(\boldsymbol{x}_{i}-\overline{\boldsymbol{x}}\right)^{T} \Sigma^{-1} d \Sigma\right)\\ =&\operatorname{tr}\left(\Sigma^{-1} S \Sigma^{-1} d \Sigma\right) \end{align}\]其中
\[S=\frac{1}{N} \sum_{i=1}^{N}\left(\boldsymbol{x}_{i}-\overline{\boldsymbol{x}}\right)\left(\boldsymbol{x}_{i}-\overline{\boldsymbol{x}}\right)^{T}\]结合两项得到：
\[d l=\operatorname{tr}\left(\left(\Sigma^{-1}-\Sigma^{-1} S \Sigma^{-1}\right) d \Sigma\right)\]对照导数与微分的关系得到：
\[\frac{\partial l}{\partial \Sigma}=\left(\Sigma^{-1}-\Sigma^{-1} S \Sigma^{-1}\right)^{T}\]其零点即为 $\Sigma$ 极大似然估计 $\hat \Sigma = S$。
参考资料
Petersen K B, Pedersen M S. The matrix cookbook[J]. Technical University of Denmark, 2008, 7(15): 510.
张贤达. 矩阵分析与应用[M]. 清华大学出版社有限公司, 2004.
