矩阵求导（上）——标量对矩阵的求导

参考：https://zhuanlan.zhihu.com/p/24709748
这部分内容分两篇整理，上篇讲标量对矩阵的求导，下篇讲矩阵对矩阵的求导。

• 本文使用小写字母x表示标量，粗体小写字母$\mathbf{x}$表示（列）向量，大写字母X表示矩阵。

定义：标量f对矩阵X的导数，定义为$\frac{\partial f}{\partial X}=\left[\frac{\partial f}{\partial X_{ij}}\right]$，即f对X逐元素求导排成与X尺寸相同的矩阵。然而，这个定义在计算中并不好用，实用上的原因是在对较复杂的函数难以逐元素求导；哲理上的原因是逐元素求导破坏了整体性。

试想，为何要将f看做矩阵X而不是各元素$X_{ij}$的函数呢？
答案是用矩阵运算更整洁。所以在求导时不宜拆开矩阵，而是要找一个从整体出发的算法。

为此，我们来回顾，一元微积分中的导数（标量对标量的导数）与微分有联系：$df=f^{\prime }(x)dx$；
多元微积分中的梯度（标量对向量的导数）也与微分有联系：$df={\sum }_{i=1}^n\frac{\partial f}{\partial x_i}d{x}_i={\frac{\partial f}{\partial \mathbf{x}}}^{T}d\mathbf{x}$，这里第一个等号是全微分公式，第二个等号表达了梯度与微分的联系：全微分$df是n\times 1$梯度向量$\frac{\partial f}{\partial \mathbf{x}}$与$n\times 1$微分向量$d\mathbf{x}$的内积；

受此启发，我们将矩阵导数与微分建立联系：$df={\sum }_{i=1}^m{\sum }_{j=1}^n\frac{\partial f}{\partial X_{ij}}d{X}_{ij}=\text{tr}\left({\frac{\partial f}{\partial X}}^{T}dX\right)$。其中$tr$代表迹(trace)是方阵对角线元素之和，满足性质：对尺寸相同的矩阵A,B，$\text{tr}(A^{T}B)={\sum }_{i,j}{A}_{ij}{B}_{ij}，即\text{tr}(A^{T}B)$是矩阵A,B的内积。与梯度相似，这里第一个等号是全微分公式，第二个等号表达了矩阵导数与微分的联系：全微分df是$m\times n导数\frac{\partial f}{\partial X}与m\times n$微分矩阵dX的内积。

然后来建立运算法则。回想遇到较复杂的一元函数如$f=\log (2+\sin x)e^{\sqrt{x}}$，我们是如何求导的呢？通常不是从定义开始求极限，而是先建立了初等函数求导和四则运算、复合等法则，再来运用这些法则。故而，我们来创立常用的矩阵微分的运算法则：

• 加减法：$d(X\pm Y)=dX\pm dY$；
• 矩阵乘法：$d(XY)=(dX)Y+XdY$ ；
• 转置：$d(X^T)=(dX{)}^T$；
• 迹：$d\text{tr}(X)=\text{tr}(dX)$。
• 逆：$d{X}^{-1}=-X^{-1}dX{X}^{-1}$。此式可在$X{X}^{-1}=I$两侧求微分来证明。
• 行列式：$d|X|=\text{tr}(X^{\#}dX)$ ，其中$X^{\#}$表示X的伴随矩阵，在X可逆时又可以写作$d|X|=|X|\text{tr}(X^{-1}dX)$。此式可用Laplace展开来证明，详见张贤达《矩阵分析与应用》第279页。
• 逐元素乘法：$d(X\odot Y)=dX\odot Y+X\odot dY$，$\odot$表示尺寸相同的矩阵X,Y逐元素相乘。
• 逐元素函数：$d\sigma (X)={\sigma }^{\prime }(X)\odot dX$ ，$\sigma (X)=\left[\sigma (X_{ij})\right]$是逐元素标量函数运算， ${\sigma }^{\prime }(X)=[{\sigma }^{\prime }(X_{ij})]$是逐元素求导数。举个例子，$X=[x_1,x_2],d\sin (X)=[\cos x_{1}d{x}_1,\cos x_{2}d{x}_2]=\cos (X)\odot dX$。

我们试图利用矩阵导数与微分的联系$df=\text{tr}\left({\frac{\partial f}{\partial X}}^{T}dX\right)$，在求出左侧的微分$df$后，该如何写成右侧的形式并得到导数呢？这需要一些迹技巧(trace trick)：

1. 标量套上迹：$a=\text{tr}(a)$
2. 转置：$\mathrm{t}\mathrm{r}(A^T)=\mathrm{t}\mathrm{r}(A)$
3. 线性：$\text{tr}(A\pm B)=\text{tr}(A)\pm \text{tr}(B)$。
4. 矩阵乘法交换：$\text{tr}(AB)=\text{tr}(BA)$，其中A与$B^T$尺寸相同。两侧都等于${\sum }_{i,j}{A}_{ij}{B}_{ji}$。
5. 矩阵乘法/逐元素乘法交换：$\text{tr}(A^T(B\odot C))=\text{tr}((A\odot B{)}^{T}C)$，其中A, B, C尺寸相同。两侧都等于${\sum }_{i,j}{A}_{ij}{B}_{ij}{C}_{ij}$。

观察一下可以断言，若标量函数f是矩阵X经加减乘法、行列式、逆、逐元素函数等运算构成，则使用相应的运算法则对f求微分，再使用迹技巧给$df$套上迹并将其它项交换至$dX$左侧，即能得到导数。

在建立法则的最后，来谈一谈复合：假设已求得$\frac{\partial f}{\partial Y}$，而$Y是X$的函数，如何求$\frac{\partial f}{\partial X}$呢？在微积分中有标量求导的链式法则$\frac{\partial f}{\partial x}=\frac{\partial f}{\partial y}\frac{\partial y}{\partial x}$，但这里我们不能沿用链式法则，因为矩阵对矩阵的导数$\frac{\partial Y}{\partial X}$截至目前仍是未定义的。于是我们继续追本溯源，链式法则是从何而来？源头仍然是微分。我们直接从微分入手建立复合法则：先写出$df=\text{tr}\left({\frac{\partial f}{\partial Y}}^{T}dY\right)$，再将$dY用dX$表示出来代入，并使用迹技巧将其他项交换至$dX$左侧，即可得到$\frac{\partial f}{\partial X}$。

接下来演示一些算例。特别提醒要依据已经建立的运算法则来计算，不能随意套用微积分中标量导数的结论，比如认为AX对X的导数为A，这是没有根据、意义不明的。

• 例1：$f={\mathbf{a}}^{T}X\mathbf{b}$，求$\frac{\partial f}{\partial X}$。其中$\mathbf{a}$是m×1列向量，X是$m\times n$矩阵，$\mathbf{b}$是n×1列向量，f是标量。

  解：先使用矩阵乘法法则求微分，这里的$\mathbf{a}$, $\mathbf{b}$是常量，$d\mathbf{a}=0,d\mathbf{b}=0$，得到：$df={\mathbf{a}}^{T}dX\mathbf{b}$ ，再套上迹并做矩阵乘法交换：$df=\text{tr}({\mathbf{a}}^{T}dX\mathbf{b})=\text{tr}(\mathbf{b}{\mathbf{a}}^{T}dX)$，注意这里我们根据$\text{tr}(AB)=\text{tr}(BA)交换了{\mathbf{a}}^{T}dX与\mathbf{b}$。对照导数与微分的联系$df=\text{tr}\left({\frac{\partial f}{\partial X}}^{T}dX\right)$，得到$\frac{\partial f}{\partial X}=(\mathbf{b}{\mathbf{a}}^T{)}^T=\mathbf{a}{\mathbf{b}}^T$。

• 注意：这里不能用$\frac{\partial f}{\partial X}={\mathbf{a}}^T\frac{\partial X}{\partial X}\mathbf{b}=?$，导数与乘常数矩阵的交换是不合法则的运算（而微分是合法的）。有些资料在计算矩阵导数时，会略过求微分这一步，这是逻辑上解释不通的。

• 例2：$f={\mathbf{a}}^T\exp (X\mathbf{b})$，求$\frac{\partial f}{\partial X}$。其中$\mathbf{a}$是m×1列向量，X是$m\times n$矩阵，$\mathbf{b}$是n×1列向量，exp表示逐元素求指数，$f$是标量。
  解：先使用矩阵乘法、逐元素函数法则求微分：$df={\mathbf{a}}^T(\exp (X\mathbf{b})\odot (dX\mathbf{b}))$，再套上迹并做矩阵乘法/逐元素乘法交换、矩阵乘法交换：$df=\text{tr}({\mathbf{a}}^T(\exp (X\mathbf{b})\odot (dX\mathbf{b})))=\text{tr}((\mathbf{a}\odot \exp (X\mathbf{b}){)}^{T}dX\mathbf{b})=\text{tr}(\mathbf{b}(\mathbf{a}\odot \exp (X\mathbf{b}){)}^{T}dX)$，注意这里我们先根据$\text{tr}(A^T(B\odot C))=\text{tr}((A\odot B{)}^{T}C)$交换了$\mathbf{a}、\exp (X\mathbf{b})与dX\mathbf{b}$，再根据$\text{tr}(AB)=\text{tr}(BA)$交换了$(\mathbf{a}\odot \exp (X\mathbf{b}){)}^{T}dX$与$\mathbf{b}$。对照导数与微分的联系$df=\text{tr}\left({\frac{\partial f}{\partial X}}^{T}dX\right)$，得到$\frac{\partial f}{\partial X}=(\mathbf{b}(\mathbf{a}\odot \exp (X\mathbf{b}){)}^T{)}^T=(\mathbf{a}\odot \exp (X\mathbf{b})){\mathbf{b}}^T$。

• 例3.【线性回归】：$l=\Vert X\mathbf{w}-\mathbf{y}{\Vert }^2$， 求$\mathbf{w}$的最小二乘估计，即求$\frac{\partial l}{\partial \mathbf{w}}$的零点。其中$\mathbf{y}$是$m\times 1$列向量，X是$m\times n$矩阵，$\mathbf{w}$是n×1列向量，$l$是标量。

  解：严格来说这是标量对向量的导数，不过可以把向量看做矩阵的特例。
  先将向量模平方改写成向量与自身的内积：$l=(X\mathbf{w}-\mathbf{y}{)}^T(X\mathbf{w}-\mathbf{y})$，求微分，使用矩阵乘法、转置等法则：$dl=(Xd\mathbf{w}{)}^T(X\mathbf{w}-\mathbf{y})+(X\mathbf{w}-\mathbf{y}{)}^T(Xd\mathbf{w})=2(X\mathbf{w}-\mathbf{y}{)}^{T}Xd\mathbf{w}$。
  对照导数与微分的联系$dl={\frac{\partial l}{\partial \mathbf{w}}}^{T}d\mathbf{w}，得到\frac{\partial l}{\partial \mathbf{w}}=(2(X\mathbf{w}-\mathbf{y}{)}^{T}X{)}^T=2X^T(X\mathbf{w}-\mathbf{y})$。$\frac{\partial l}{\partial \mathbf{w}}的零点即\mathbf{w}的最小二乘估计为\mathbf{w}=(X^{T}X{)}^{-1}{X}^T\mathbf{y}$。

• 例4.【方差的最大似然估计】：样本${\mathbf{x}}_1,\dots ,{\mathbf{x}}_n\sim N(\mathbf{\mu },\Sigma )$，求方差$\Sigma$的最大似然估计。写成数学式是：$l=\log |\Sigma |+\frac{1}{n}{\sum }_{i=1}^n({\mathbf{x}}_i-\bar{\mathbf{x}}{)}^T{\Sigma }^{-1}({\mathbf{x}}_i-\bar{\mathbf{x}})$，求$\frac{\partial l}{\partial \Sigma }$的零点。其中${\mathbf{x}}_i是m\times 1$列向量，$\overline{\mathbf{x}}=\frac{1}{n}{\sum }_{i=1}^n{\mathbf{x}}_i$是样本均值，$\Sigma 是m\times m对称正定矩阵，l是标量$.

  解：首先求微分，使用矩阵乘法、行列式、逆等运算法则，
  第一项是$d\log |\Sigma |=|\Sigma {|}^{-1}d|\Sigma |=\text{tr}({\Sigma }^{-1}d\Sigma )$，
  第二项是$\frac{1}{n}{\sum }_{i=1}^n({\mathbf{x}}_i-\bar{\mathbf{x}}{)}^{T}d{\Sigma }^{-1}({\mathbf{x}}_i-\bar{\mathbf{x}})=-\frac{1}{n}{\sum }_{i=1}^n({\mathbf{x}}_i-\bar{\mathbf{x}}{)}^T{\Sigma }^{-1}d\Sigma {\Sigma }^{-1}({\mathbf{x}}_i-\bar{\mathbf{x}})$。
  再给第二项套上迹做交换：
  $\text{tr}\left(\frac{1}{n}{\sum }_{i=1}^n({\mathbf{x}}_i-\bar{\mathbf{x}}{)}^T{\Sigma }^{-1}d\Sigma {\Sigma }^{-1}({\mathbf{x}}_i-\bar{\mathbf{x}})\right)=\frac{1}{n}{\sum }_{i=1}^n\text{tr}\left(({\mathbf{x}}_i-\bar{\mathbf{x}}{)}^T{\Sigma }^{-1}d\Sigma {\Sigma }^{-1}({\mathbf{x}}_i-\bar{\mathbf{x}})\right)=\frac{1}{n}{\sum }_{i=1}^n\text{tr}\left({\Sigma }^{-1}({\mathbf{x}}_i-\bar{\mathbf{x}})({\mathbf{x}}_i-\bar{\mathbf{x}}{)}^T{\Sigma }^{-1}d\Sigma \right)=\text{tr}({\Sigma }^{-1}S{\Sigma }^{-1}d\Sigma )$，
  其中先交换迹与求和，然后将 ${\Sigma }^{-1}({\mathbf{x}}_i-\bar{\mathbf{x}})$交换到左边，最后再交换迹与求和，
  并定义$S=\frac{1}{n}{\sum }_{i=1}^n({\mathbf{x}}_i-\bar{\mathbf{x}})({\mathbf{x}}_i-\bar{\mathbf{x}}{)}^T$为样本方差矩阵。得到$dl=\text{tr}\left(\left({\Sigma }^{-1}-{\Sigma }^{-1}S{\Sigma }^{-1}\right)d\Sigma \right)$。
  对照导数与微分的联系，有$\frac{\partial l}{\partial \Sigma }=({\Sigma }^{-1}-{\Sigma }^{-1}S{\Sigma }^{-1}{)}^T$，其零点即$\Sigma$的最大似然估计为$\Sigma =S$。

• 例5【多元logistic回归】：$l=-{\mathbf{y}}^T\log \text{softmax}(W\mathbf{x})，求\frac{\partial l}{\partial W}$。其中$\mathbf{y}$是除一个元素为1外其它元素为0的m×1列向量，$W是m\times n$矩阵，$\mathbf{x}$是n×1列向量，$l$是标量；$\text{softmax}(\mathbf{a})$ = $\frac{\exp (\mathbf{a})}{1^T\exp (\mathbf{a})}$，其中$\exp (\mathbf{a})$表示逐元素求指数，$1$代表全1向量。

  解：首先将softmax函数代入并写成
  $l=-{\mathbf{y}}^T\left(\log (\exp (W\mathbf{x}))-1\log (1^T\exp (W\mathbf{x}))\right)$ = $-{\mathbf{y}}^{T}W\mathbf{x}+\log (1^T\exp (W\mathbf{x}))$，这里要注意逐元素log满足等式$\log (\mathbf{u}/c)=\log (\mathbf{u})-1\log (c)$，以及$\mathbf{y}满足{\mathbf{y}}^{T}1=1$。
  求微分，使用矩阵乘法、逐元素函数等法则：$dl=-{\mathbf{y}}^{T}dW\mathbf{x}+\frac{1^T\left(\exp (W\mathbf{x})\odot (dW\mathbf{x})\right)}{1^T\exp (W\mathbf{x})}$。再套上迹并做交换，注意可化简$1^T\left(\exp (W\mathbf{x})\odot (dW\mathbf{x})\right)=\exp (W\mathbf{x}{)}^{T}dW\mathbf{x}$，
  这是根据等式$1^T(\mathbf{u}\odot \mathbf{v})={\mathbf{u}}^T\mathbf{v}$，故$dl=\text{tr}\left(-{\mathbf{y}}^{T}dW\mathbf{x}+\frac{\exp (W\mathbf{x}{)}^{T}dW\mathbf{x}}{1^T\exp (W\mathbf{x})}\right)=\text{tr}(\mathbf{x}(\text{softmax}(W\mathbf{x})-\mathbf{y}{)}^{T}dW)$。
  对照导数与微分的联系，得到$\frac{\partial l}{\partial W}=(\text{softmax}(W\mathbf{x})-\mathbf{y}){\mathbf{x}}^T$。

  另解：定义$\mathbf{a}=W\mathbf{x}$，则$l=-{\mathbf{y}}^T\log \text{softmax}(\mathbf{a})，先如上求出\frac{\partial l}{\partial \mathbf{a}}=\text{softmax}(\mathbf{a})-\mathbf{y}$ ，
  再利用复合法则：$dl=\text{tr}\left({\frac{\partial l}{\partial \mathbf{a}}}^{T}d\mathbf{a}\right)=\text{tr}\left({\frac{\partial l}{\partial \mathbf{a}}}^{T}dW\mathbf{x}\right)=\text{tr}\left(\mathbf{x}{\frac{\partial l}{\partial \mathbf{a}}}^{T}dW\right)，得到\frac{\partial l}{\partial W}=\frac{\partial l}{\partial \mathbf{a}}{\mathbf{x}}^T$。

• 最后一例留给经典的神经网络。神经网络的求导术是学术史上的重要成果，还有个专门的名字叫做BP算法， 我相信如今很多人在初次推导BP算法时也会颇费一番脑筋，事实上使用矩阵求导术来推导并不复杂。
  为简化起见，我们推导二层神经网络的BP算法。

• 例6.【二层神经网络】：$l=-{\mathbf{y}}^T\log \text{softmax}(W_2\sigma (W_1\mathbf{x}))，求\frac{\partial l}{\partial W_1}和\frac{\partial l}{\partial W_2}$。其中$\mathbf{y}$是除一个元素为1外其它元素为0的的$m\times 1$列向量，$W_2是m\times p$矩阵，$W_1是p\times n$矩阵，$\mathbf{x}$是n×1列向量，$l$是标量；$\text{softmax}(\mathbf{a})=\frac{\exp (\mathbf{a})}{1^T\exp (\mathbf{a})}$,同例3，$\sigma (\cdot )$是逐元素$sigmoid函数\sigma (a)=\frac{1}{1+\exp (-a)}$。

  解：定义${\mathbf{a}}_1=W_1\mathbf{x}，{\mathbf{h}}_1=\sigma ({\mathbf{a}}_1)，{\mathbf{a}}_2=W_2{\mathbf{h}}_1，则l=-{\mathbf{y}}^T\log \text{softmax}({\mathbf{a}}_2)$。
  在前例中已求出$\frac{\partial l}{\partial {\mathbf{a}}_2}=\text{softmax}({\mathbf{a}}_2)-\mathbf{y}$ 。使用复合法则，注意此处${\mathbf{h}}_1,W_2$都是变量：
  $dl=\text{tr}\left({\frac{\partial l}{\partial {\mathbf{a}}_2}}^{T}d{\mathbf{a}}_2\right)=\text{tr}\left({\frac{\partial l}{\partial {\mathbf{a}}_2}}^{T}d{W}_2{\mathbf{h}}_1\right)+\text{tr}\left({\frac{\partial l}{\partial {\mathbf{a}}_2}}^T{W}_{2}d{\mathbf{h}}_1\right)$，
  使用矩阵乘法交换的迹技巧从第一项得到$\frac{\partial l}{\partial W_2}=\frac{\partial l}{\partial {\mathbf{a}}_2}{\mathbf{h}}_1^T$，从第二项得到$\frac{\partial l}{\partial {\mathbf{h}}_1}=W_2^T\frac{\partial l}{\partial {\mathbf{a}}_2}$。
  $接下来求\frac{\partial l}{\partial {\mathbf{a}}_1}$，继续使用复合法则，并利用矩阵乘法和逐元素乘法交换的迹技巧：
  $\text{tr}\left({\frac{\partial l}{\partial {\mathbf{h}}_1}}^{T}d{\mathbf{h}}_1\right)=\text{tr}\left({\frac{\partial l}{\partial {\mathbf{h}}_1}}^T({\sigma }^{\prime }({\mathbf{a}}_1)\odot d{\mathbf{a}}_1)\right)=\text{tr}\left({\left(\frac{\partial l}{\partial {\mathbf{h}}_1}\odot {\sigma }^{\prime }({\mathbf{a}}_1)\right)}^{T}d{\mathbf{a}}_1\right)$，得到$\frac{\partial l}{\partial {\mathbf{a}}_1}=\frac{\partial l}{\partial {\mathbf{h}}_1}\odot {\sigma }^{\prime }({\mathbf{a}}_1)$.
  $为求\frac{\partial l}{\partial W_1}$，再用一次复合法则：$\text{tr}\left({\frac{\partial l}{\partial {\mathbf{a}}_1}}^{T}d{\mathbf{a}}_1\right)=\text{tr}\left({\frac{\partial l}{\partial {\mathbf{a}}_1}}^{T}d{W}_1\mathbf{x}\right)=\text{tr}\left(\mathbf{x}{\frac{\partial l}{\partial {\mathbf{a}}_1}}^{T}d{W}_1\right)，得到\frac{\partial l}{\partial W_1}=\frac{\partial l}{\partial {\mathbf{a}}_1}{\mathbf{x}}^T$。