矩阵求导术

矩阵求导的技术，在统计学、控制论、机器学习等领域有广泛的应用。本文来做个科普，分作两篇，上篇讲标量对矩阵的求导术，下篇讲矩阵对矩阵的求导术。本文使用小写字母$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^T}{\partial \mathbf{x}}d\mathbf{x}$，这里第一个等号是全微分公式，第二个等号表达了梯度与微分的联系：全微分$df$是梯度向量$\frac{\partial f}{\partial \mathbf{x}}(n\times 1)$与微分向量$d\mathbf{x}(x\times 1)$的内积；受此启发，我们将矩阵导数与微分建立联系：$df={\sum }_{i=1}^m{\sum }_{j=1}^n\frac{\partial f}{\partial X_{ij}}d{X}_{ij}=tr\left(\frac{\partial f^T}{\partial X}dX\right)$​​​​​​​。其中tr代表迹(trace)是方阵对角线元素之和，满足性质：对尺寸相同的矩阵A,B，$tr\left(A^{T}B\right)={\sum }_{i,j}{A}_{ij}{B}_{ij}$​​​​​​​，即$tr\left(A^{T}B\right)$是矩阵A,B的内积。与梯度相似，这里第一个等号是全微分公式，第二个等号表达了矩阵导数与微分的联系：全微分$df$是导数$\frac{\partial f}{\partial X}(m\times n)$与微分矩阵$dX(m\times n)$的内积。

然后来建立运算法则。回想遇到较复杂的一元函数如$f=log(2+sinx)e^{\sqrt{x}}$，我们是如何求导的呢？通常不是从定义开始求极限，而是先建立了初等函数求导和四则运算、复合等法则，再来运用这些法则。故而，我们来创立常用的矩阵微分的运算法则：
1.
加减法：$d(X\pm Y)=dX\pm dY$
矩阵乘法：$d(XY)=(dX)Y+X(dY)$
转置：$d(X^T)=(dX{)}^T$
迹：$dtr(X)=tr(dX)$
2.
逆：$d{X}^{-1}=-X^{-1}dX{X}^{-1}$。此式可在$X{X}^{-1}=I$两侧求微分来证明。
3.
行列式：$d\left|X\right|=tr(X^{\#}dX)$，其中$X^{\#}$表示$X$的伴随矩阵，在$X$可逆时又可以写作$d\left|X\right|=\left|X\right|tr\left(X^{-1}dX\right)$。此式可用Laplace展开来证明，详见张贤达《矩阵分析与应用》第279页。
4.
逐元素乘法：$d(X\odot Y)=dX\odot Y+X\odot dY$，$\odot$表示尺寸相同的矩阵X,Y逐元素相乘。
5.
逐元素函数：$d\sigma (X)={\sigma }^{\prime }(X)\odot dX,\sigma (X)=\left[\sigma (X_{ij})\right]$，是逐元素标量函数运算，${\sigma }^{\prime }(X)=\left[{\sigma }^{\prime }(X_{ij})\right]$是逐元素求导数。举个例子:
$$X=\left[\begin{array}{cc}{x}_{11}& x_{12}\\ x_{21}& x_{22}\end{array}\right],d\sin (X)=\left[\begin{array}{cc}cos{x}_{11}d{x}_{11}& cos{x}_{12}d{x}_{12}\\ cos{x}_{21}d{x}_{21}& cos{x}_{22}d{x}_{22}\end{array}\right]=cos(X)\odot dX$$
我们试图利用矩阵导数与微分的联系$df=tr\left(\frac{\partial f^T}{\partial X}dX\right)$，在求出左侧的微分$df$后，该如何写成右侧的形式并得到导数呢？这需要一些迹技巧(trace trick)：
1.标量套上迹：$a=tr(a)$
2.转置：$tr(A^T)=tr(A)$
3.线性：$tr(A\pm B)=tr(A)\pm tr(B)$
4.矩阵乘法交换：$tr(AB)=tr(BA)$，其中$A$与$B^T$尺寸相同。两侧都等于${\sum }_{ij}{A}_{ij}{B}_{ji}$。
5.矩阵乘法/逐元素乘法交换：$tr(A^T(B\odot C))=tr((A\odot B{)}^{T}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=tr(\frac{\partial f^T}{\partial Y}dY)$再将$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\times 1$列向量，$X$是$m\times n$矩阵，$\mathbf{b}$是$n\times 1$列向量, $f$是标量。

解：先使用矩阵乘法法则求微分，这里的$\mathbf{a},\mathbf{b}$是常量，$d\mathbf{a}=0,d\mathbf{b}=0$，得到：$df={\mathbf{a}}^{T}dX\mathbf{b}$，再套上迹并做矩阵乘法交换：$df=tr({\mathbf{a}}^{T}dX\mathbf{b})=tr(\mathbf{b}{\mathbf{a}}^{T}dX)$，注意这里我们根据$tr(AB)=tr(BA)$交换了${\mathbf{a}}^{T}dX$与$\mathbf{b}$。对照导数与微分的联系$df=tr\left(\frac{\partial f^T}{\partial X}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\times 1$列向量，其中$X$是$m\times n$矩阵，其中$\mathbf{b}$是$n\times 1$列向量，$exp$表示逐元素求指数，$f$是标量。

解：先使用矩阵乘法、逐元素函数法则求微分：$df={\mathbf{a}}^T(exp(X\mathbf{b}))\odot (dX\mathbf{b})$，再套上迹并做交换：$df=tr({\mathbf{a}}^T(exp(X\mathbf{b}))\odot (dX\mathbf{b})))=tr((\mathbf{a}\odot exp(X\mathbf{b}){)}^{T}dX\mathbf{b})=tr(\mathbf{b}(\mathbf{a}\odot exp(X\mathbf{b}){)}^{T}dX)$，注意这里我们先根据$tr(A^T(B\odot C))=tr((A\odot B{)}^{T}C)$交换了$\mathbf{a},exp(X\mathbf{b})$与$dX\mathbf{b}$，再根据$tr(AB)=tr(BA)$交换了$(\mathbf{a}\odot exp(X\mathbf{b}){)}^{T}dX$与$\mathbf{b}$。对照导数与微分的联系$df=tr\left(\frac{\partial f^T}{\partial X}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：$f=tr(Y^{T}MY),Y=\sigma (WX)$求$\frac{\partial f}{\partial X}$其中$W$是$l\times m$矩阵，其中$X$是$m\times n$矩阵，其中$Y$是$l\times n$矩阵，其中$M$是$l\times l$矩阵，$\sigma$是逐元素函数，$f$是标量。

解：先求$\frac{\partial f}{\partial Y}$，求微分，使用矩阵乘法、转置法则：$df=tr((dY{)}^{T}MY)+tr(Y^{T}MdY)=tr(Y^T{M}^{T}dY)+tr(Y^{T}MdY)=tr(Y^T(M+M^T)dY)$，对照导数与微分的联系，得到$\frac{\partial f}{\partial Y}=(M+M^T)Y=2MY$，这里是对称矩阵。为求$\frac{\partial f}{\partial X}$，写出$df=tr\left(\frac{\partial f^T}{\partial Y}dY\right)$，再将$dY$用$dX$表示出来代入，并使用矩阵乘法/逐元素乘法交换：$df=tr\left(\frac{\partial f^T}{\partial Y}({\sigma }^{\prime }(WX)\odot (WdX)\right)=tr\left({\left(\frac{\partial f}{\partial Y}\odot {\sigma }^{\prime }(WX)\right)}^{T}WdX\right)$，对照导数与微分的联系，得到$\frac{\partial f}{\partial X}=W^T\left(\frac{\partial f}{\partial Y}\odot {\sigma }^{\prime }(WX)\right)=W^T((2M\sigma (WXX))\odot {\sigma }^{\prime }(WX))$

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

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

例5【方差的最大似然估计】：样本${\mathbf{x}}_1,...,{\mathbf{x}}_N\sim \lambda (\mathbf{\mu },\Sigma )$，求方差$\Sigma$的最大似然估计。写成数学式是：$l=log\left|\Sigma \right|+\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$向量，$\Sigma$是$m\times m$对称正定矩阵，$\bar{\mathbf{x}}=\frac{1}{N}{\sum }_{i=1}^N{\mathbf{x}}_i$是样本均值，$l$是标量，log表示自然对数。

解：首先求微分，使用矩阵乘法、行列式、逆等运算法则，第一项是$dlog\left|\Sigma \right|={\left|\Sigma \right|}^{-1}d\left|\Sigma \right|=tr({\left|\Sigma \right|}^{-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}})$。再给第二项套上迹做交换：$tr(\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}}))=\frac{1}{N}{\sum }_{i=1}^{N}tr(({\mathbf{x}}_i-\bar{\mathbf{x}}{)}^T{\Sigma }^{-1}d\Sigma {\Sigma }^{-1}({\mathbf{x}}_i-\bar{\mathbf{x}}))=\frac{1}{N}{\sum }_{i=1}^{N}tr({\Sigma }^{-1}({\mathbf{x}}_i-\bar{\mathbf{x}})({\mathbf{x}}_i-\bar{\mathbf{x}}{)}^T{\Sigma }^{-1}d\Sigma =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=tr(({\Sigma }^{-1}-{\Sigma }^{-1}S{\Sigma }^{-1})d\Sigma )$。对照导数与微分的联系，有$\frac{\partial l}{\partial \Sigma }=({\Sigma }^{-1}-{\Sigma }^{-1}S{\Sigma }^{-1}{)}^T$，其零点即$\Sigma$的最大似然估计为$\Sigma =S$。

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

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

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

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

例7【二层神经网络】：$l=-{\mathbf{y}}^{T}logsoftmax(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\times 1$矩阵，$l$是标量；log表示自然对数，$softmax(\mathbf{a})=\frac{exp(\mathbf{a})}{1^{T}exp(\mathbf{a})}$同上，$\sigma$是逐元素sigmoid函数。

解：定义${\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}logsoftmax({\mathbf{a}}_2)$。在前例中已求出$\frac{\partial l}{\partial {\mathbf{a}}_2}=softmax({\mathbf{a}}_2)-\mathbf{y}$。使用复合法则，$dl=tr\left(\frac{\partial l^T}{\partial {\mathbf{a}}_2}d{\mathbf{a}}_2\right)=tr\left(\frac{\partial l^T}{\partial {\mathbf{a}}_2}d{W}_2{\mathbf{h}}_1\right)+\underset{d{l}_2}{\underbrace{tr\left(\frac{\partial l^T}{\partial {\mathbf{a}}_2}{W}_{2}d{\mathbf{h}}_1\right)}}$，使用矩阵乘法交换的迹技巧从第一项得到$\frac{\partial l}{\partial W_2}=\frac{\partial l}{\partial {\mathbf{a}}_2}{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}$，并利用矩阵乘法和逐元素乘法交换的迹技巧：$d{l}_2=tr\left(\frac{\partial l^T}{\partial {\mathbf{h}}_1}d{\mathbf{h}}_1\right)=tr\left(\frac{\partial l^T}{\partial {\mathbf{h}}_1}({\sigma }^{\prime }({\mathbf{a}}_1)\odot d{\mathbf{a}}_1)\right)=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}$，再用一次复合法则：$d{l}_2=tr\left(\frac{\partial l^T}{\partial {\mathbf{a}}_1}d{\mathbf{a}}_1\right)=tr\left(\frac{\partial l^T}{\partial {\mathbf{a}}_1}d{W}_1\mathbf{x}\right)=tr\left(\mathbf{x}\frac{\partial l^T}{\partial {\mathbf{a}}_1}d{W}_1\right)$，得到$\frac{\partial l}{\partial W_1}=\frac{\partial l}{\partial {\mathbf{a}}_1}{\mathbf{x}}^T$

推广：样本$({\mathbf{x}}_1,{\mathbf{y}}_1),....,({\mathbf{x}}_N,{\mathbf{y}}_N)$，$l=-{\sum }_{i=1}^N{\mathbf{y}}_i^{T}logsoftmax(W_2\sigma (W_1{\mathbf{x}}_i+{\mathbf{b}}_1)+{\mathbf{b}}_2)$，${\mathbf{b}}_1$是$p\times 1$列向量，${\mathbf{b}}_2$是$m\times 1$列向量其余定义同上。

解1：定义${\mathbf{a}}_{1,i}=W_1{\mathbf{x}}_i+{\mathbf{b}}_1,{\mathbf{h}}_{1,i}=\sigma ({\mathbf{a}}_{1,i}),{\mathbf{a}}_{2,i}=W_2{\mathbf{h}}_{1,i}+{\mathbf{b}}_2$​​​​​，则$l=-{\sum }_{i=1}^N{\mathbf{y}}_i^{T}logsoftmax({\mathbf{a}}_{2,i})$。先同上可求出$\frac{\partial l}{\partial {\mathbf{a}}_{2,i}}=softmax({\mathbf{a}}_{2,i})-{\mathbf{y}}_i$。使用复合法则，$dl=tr\left({\sum }_{i=1}^N\frac{\partial l^T}{\partial {\mathbf{a}}_{2,i}}d{\mathbf{a}}_{2,i}\right)=tr\left({\sum }_{i=1}^N\frac{\partial l^T}{\partial {\mathbf{a}}_{2,i}}d{W}_2{\mathbf{h}}_{1,i}\right)+\underset{d{l}_2}{\underbrace{tr\left(\sum _{i=1}^N\frac{\partial l^T}{\partial {\mathbf{a}}_{2,i}}{W}_{2}d{\mathbf{h}}_{1,i}\right)}}+tr\left({\sum }_{i=1}^N\frac{\partial l^T}{\partial {\mathbf{a}}_{2,i}}d{\mathbf{b}}_2\right)$，从第一项得到得到$\frac{\partial l}{\partial W_2}={\sum }_{i=1}^N\frac{\partial l}{\partial {\mathbf{a}}_{2,i}}{\mathbf{h}}_{1,i}^T$，从第二项得到$\frac{\partial l}{\partial {\mathbf{h}}_{1,i}}=W_2^T\frac{\partial l}{\partial {\mathbf{a}}_{2,i}}$，从第三项得到$\frac{\partial l}{\partial {\mathbf{b}}_2}={\sum }_{i=1}^N\frac{\partial l}{\partial {\mathbf{a}}_{2,i}}$。接下来对第二项继续使用复合法则，得到$\frac{\partial l}{\partial {\mathbf{a}}_{1,i}}=\frac{\partial l}{\partial {\mathbf{h}}_{1,i}}\odot {\sigma }^{\prime }({\mathbf{a}}_{1,i})$。为求$\frac{\partial l}{\partial W_1},\frac{\partial l}{\partial {\mathbf{b}}_1}$，再用一次复合法则：$d{l}_2=tr\left({\sum }_{i=1}^N\frac{\partial l^T}{\partial a_{1,i}}d{\mathbf{a}}_{1,i}\right)=tr\left({\sum }_{i=1}^N\frac{\partial l^T}{\partial {\mathbf{a}}_{1,i}}d{W}_1{\mathbf{x}}_i\right)+tr\left({\sum }_{i=1}^N\frac{\partial l^T}{\partial {\mathbf{a}}_{1,i}}d{\mathbf{b}}_1\right)$得到$\frac{\partial l}{\partial W_1}={\sum }_{i=1}^N\frac{\partial l}{\partial {\mathbf{a}}_{1,i}}{x}_i^T,\frac{\partial l}{\partial {\mathbf{b}}_1}={\sum }_{i=1}^N\frac{\partial l}{\partial {\mathbf{a}}_{1,i}}$

解2：可以用矩阵来表示N个样本，以简化形式。定义$X=[{\mathbf{x}}_1,...,{\mathbf{x}}_N],A_1=[{\mathbf{a}}_{1,1},..,{\mathbf{a}}_{1,N}]=W_{1}X+{\mathbf{b}}_1{1}^T,H_1=[{\mathbf{h}}_{1,1},...,{\mathbf{h}}_{1,N}]=\sigma (A_1),A_2=[{\mathbf{a}}_{2,1},..,{\mathbf{a}}_{2,N}]=W_2{H}_1+{\mathbf{b}}_2{1}^T$，注意这里使用全1向量来扩展维度。先同上求出$\frac{\partial l}{\partial A_2}=[softmax({\mathbf{a}}_{2,1})-{\mathbf{y}}_1,...,softmax({\mathbf{a}}_{2,N})-{\mathbf{y}}_N]$。使用复合法则，$dl=tr\left(\frac{\partial l^T}{\partial A_2}d{A}_2\right)=tr\left(\frac{\partial l^T}{\partial A_2}d{W}_2{H}_1\right)+\underset{d{l}_2}{\underbrace{tr\left(\frac{\partial l^T}{\partial A_2}{W}_{2}d{H}_1\right)}}+tr\left(\frac{\partial l^T}{\partial A_2}d{\mathbf{b}}_2{1}^T\right)$，从第一项得到$\frac{\partial l}{\partial W_2}=\frac{\partial l}{\partial A_2}{H}_1^T$，从第二项得到$\frac{\partial l}{\partial W_1}=W_2^T\frac{\partial l}{\partial A_2}$，从第三项得到到$\frac{\partial l}{\partial {\mathbf{b}}_2}=\frac{\partial l}{\partial A_2}1$。接下来对第二项继续使用复合法则，得到$\frac{\partial l}{\partial A_1}=\frac{\partial l}{\partial H_1}\odot {\sigma }^{\prime }(A_1)$。为求$\frac{\partial l}{\partial W_1},\frac{\partial l}{\partial {\mathbf{b}}_1}$，再用一次复合法则：$d{l}_2=tr\left(\frac{\partial l^T}{\partial A_1}d{A}_1\right)=tr\left(\frac{\partial l^T}{\partial A_1}d{W}_{1}X\right)+tr\left(\frac{\partial l^T}{\partial A_1}d{\mathbf{b}}_1{1}^T\right)$，得到$\frac{\partial l}{\partial W_1}=\frac{\partial l}{\partial A_1}{X}^T,\frac{\partial l}{\partial {\mathbf{b}}_1}=\frac{\partial l}{\partial A_1}1$