本文簡單整理了以下內容：
（一）線性回歸
（二）二分類：二項Logistic回歸
（三）多分類：Softmax回歸
（四）廣義線性模型
閑話：二項Logistic回歸是我去年入門機器學習時學的第一個模型（忘記了為什么看完《統計學習方法》第一章之后直接就跳去了第六章，好像是對“邏輯斯蒂”這個名字很感興趣？。。。），對照《機器學習實戰》寫了幾行代碼敲了一個toy版本，當時覺得還是挺有意思的。我覺得這個模型很適合用來入門（但是必須注意這個模型有很多很多很多很多可以展開的地方）。更有意思的是那時候還不會矩陣微積分，推導梯度時還是把矩陣全都展開求然后再歸納成矩陣形式的（牛頓法要用的二階梯度也是）。。。
下面的文字中，“Logistic回歸”都表示用於二分類的二項Logistic回歸。
首先約定一下記號。
樣本的特征矩陣 $X=(\textbf x_1,\textbf x_2,...,\textbf x_N)=({\textbf x^{(1)}}^{\top};{\textbf x^{(2)}}^{\top};...;{\textbf x^{(d)}}^{\top})\in\mathbb R^{d\times N}$ ，$X_{ji}=x_i^{(j)}$；
$N$ 是訓練集的樣本數，每個樣本都表示成 $\textbf x_i\in\mathbb R^d$ 的列向量，真實標簽為 $y_i$ ；如果出現了 $\textbf x$ 這樣沒有上下標的記號就泛指任一樣本，相當於省略下標，真實標簽為 $y$ ；特別地，對於二分類問題，$y$ 只可能取0、1兩個值。
$d$ 是特征的維數，每維特征都表示成 $\textbf x^{(j)}\in\mathbb R^N$ 的列向量；如果出現了 $x^{(j)}$ 這樣的記號就泛指任一樣本的第 $j$ 維特征，相當於省略下標；
權重向量 $\textbf w=(w_1,w_2,...,w_d)^{\top}\in\mathbb R^d$ ，偏置 $b\in\mathbb R$ 。
$\textbf y$ 在本文可能表達兩種含義：一種是表示全部訓練樣本的真實標簽組成的列向量 $\textbf y=(y_1,y_2,...,y_N)^{\top}\in\mathbb R^N$ ；另一種含義則是表示樣本 $\textbf x$ 的one-hot表示 $\textbf y=(0,0,...,0,1,0,...,0)^{\top}\in\mathbb R^C$（只有真實類別的那一維是1，其他維均是0），相當於 $\textbf y_i$ 省略了下標。
可能看起來有點別扭，因為對於樣本來說，下標是序號索引、上標是特征索引；而對於權重來說，下標是特征索引。
（一）線性回歸
1. 概述
線性回歸（Linear regression）就是用一個超平面去擬合樣本點的標簽：
$$f(\textbf x)=\textbf w^{\top}\textbf x+b$$
對於一維特征的情況，就是用一條直線去擬合樣本點，如下圖所示。為了方便起見，將偏置也記到權重向量中並保持記號不變，同時每個樣本增加一維特征並保持記號不變：$\textbf w=(1,w_1,w_3,...,w_d)^{\top}$ ，$\textbf x=(1,x^{(1)},x^{(2)},...,x^{(d)})^{\top}$ ，$f(\textbf x)=\textbf w^{\top}\textbf x$ 。
圖片來源：[1]
2. 求解
對於回歸任務，最常使用的損失函數是平方損失函數 $L(y,f(\textbf x))=(y-f(\textbf x))^2$ ，對應的經驗風險就是均方誤差（Mean square error，MSE）：
$$\mathcal R=\frac1N\sum_{i=1}^N(y_i-f(\textbf x_i))^2=\frac1N\|X^{\top}\textbf w-\textbf y\|^2=\frac1N(X^{\top}\textbf w-\textbf y)^{\top}(X^{\top}\textbf w-\textbf y)$$
該式的 $\textbf y$ 表示全部訓練樣本的真實標簽組成的列向量 $\textbf y=(y_1,y_2,...,y_N)^{\top}\in\mathbb R^N$ 。
解一：正規方程組（Normal equations）。可以直接用 $R$ 的一階導數等於0來求極值點（省略常系數）：
$$\frac{\partial \mathcal R}{\partial\textbf w}=2X(X^{\top}\textbf w-\textbf y)=0\Rightarrow \textbf w=(XX^{\top})^{-1}X\textbf y$$
可以看出，這個不就是最小二乘法（Ordinary Least Squares，OLS）解方程 $X^{\top}\textbf w=\textbf y$ 嘛。值得注意的是 $(XX^{\top})^{-1}X$ 其實就是 $X^{\top}$ 的偽逆，計算偽逆的復雜度很高。
需要注意一個問題：$XX^{\top}$ 需要是可逆矩陣，也就是說每維特征之間線性無關，才可以求得唯一解。當其不可逆（特征個數比樣本個數還要多）時，解不唯一，需要用梯度下降（Gradient descent）來迭代求解。另外，最小二乘法的結果總是低偏差高方差的。
（注：部分求導法則，使用分母布局，維度為 $q$ 的列向量對維度為 $p$ 的列向量求導后得到的矩陣維數為 $p\times q$ 。關於矩陣求導，可以參考 [8] 。
圖片來源：[1]
$$\frac{\partial A^{\top}\textbf x}{\partial\textbf x}=\frac{\partial \textbf x^{\top}A}{\partial\textbf x}=A$$
$$\frac{\partial \textbf y^{\top}\textbf z}{\partial\textbf x}=\frac{\partial \textbf y}{\partial\textbf x}\textbf z+\frac{\partial \textbf z}{\partial\textbf x}\textbf y$$
$$\frac{\partial \textbf y^{\top}A\textbf z}{\partial\textbf x}=\frac{\partial \textbf y}{\partial\textbf x}A\textbf z+\frac{\partial \textbf z}{\partial\textbf x}A^{\top}\textbf y$$
$$\frac{\partial y\textbf z}{\partial\textbf x}=\frac{\partial y}{\partial\textbf x}\textbf z^{\top}+y\frac{\partial \textbf z}{\partial\textbf x}$$
$$\frac{\partial \text{tr}AB}{\partial A}=B^{\top}\quad\quad\frac{\partial \text{tr}AB}{\partial A^{\top}}=B$$
$$\frac{\partial f(A)}{\partial A^{\top}}=(\frac{\partial f(A)}{\partial A})^{\top}$$
）
解二：最小均方誤差（least mean squares，LMS）規則，也叫Widrow-Hoff規則，用梯度下降法求解。梯度在上面已經求出來了：
$$\frac{\partial \mathcal R}{\partial\textbf w}=2X(X^{\top}\textbf w-\textbf y)=2X(\hat{\textbf y}-\textbf y)$$
該式的 $\hat{\textbf y}$ 表示模型對全部訓練樣本的輸出標簽組成的列向量 $\hat{\textbf y}=(\hat y_1,\hat y_2,...,\hat y_N)^{\top}\in\mathbb R^N$ 。
這樣的方式是每更新一次參數就要計算整個訓練集上的梯度，是批梯度下降（batch GD）；如果把這個過程拆成 $N$ 次，也就是每次只隨機挑選一個樣本計算梯度，就是隨機梯度下降（Stochastic GD，SGD）。還有一種是mini-batch梯度下降，每次挑選一個小批量樣本計算梯度。整個訓練集計算完一次梯度稱為“一輪”。
3. 均方誤差優化目標的概率解釋
重新考慮以下問題：設樣本的特征和標簽存在關系 $y_i=\textbf w^{\top}\textbf x_i+\epsilon_i$ ，並假設每個 $\epsilon_i$ 都是服從高斯分布的隨機變量 $\epsilon\sim N(0,\sigma^2)$ 的iid樣本（之所以假設為高斯分布，是認為誤差由多個獨立的隨機因素構成，根據多個獨立隨機變量之和趨於高斯分布，所以假設 $\epsilon$ 服從高斯分布）。從而有
$$p(\epsilon_i)=\frac1{\sqrt{2\pi}\sigma}\exp(-\frac{\epsilon_i^2}{2\sigma^2})$$
也就是說
$$p(y_i|\textbf x_i;\textbf w)=\frac1{\sqrt{2\pi}\sigma}\exp(-\frac{(y_i-\textbf w^{\top}\textbf x_i)^2}{2\sigma^2})$$
用分號隔開是因為在頻率學派的觀點下 $\textbf w$ 不是隨機變量。
進一步用極大似然估計來求取參數 $\textbf w$ ：對數似然函數為
$$l(\textbf w)=\log\prod_{i=1}^Np(y_i|\textbf x_i;\textbf w)=\log\prod_{i=1}^N\frac1{\sqrt{2\pi}\sigma}\exp(-\frac{(y_i-\textbf w^{\top}\textbf x_i)^2}{2\sigma^2})$$
再往后寫一步就可以知道，極大似然估計和最小化均方誤差是一致的。
4. 局部加權回歸（Locally weighted regression，LWR）
相比於普通的線性回歸，LWR對於一個點來說使用其附近的點來做回歸（而不是全部點）。
相比於線性回歸的優化目標 $\sum_{i=1}^N(y_i-\textbf w^{\top}\textbf x_i)^2$ ，局部加權線性回歸的優化目標為
$$\sum_{i=1}^N\omega_i(y_i-\textbf w^{\top}\textbf x_i)^2$$
式中 $\omega_i$ 就是非負值的權重，一個常用的選擇為 $\omega_i=\exp(-\frac{(\textbf x_i-\textbf x)^{\top}(\textbf x_i-\textbf x)}{2\tau^2})$ ，$\tau$ 是指定的帶寬參數（bandwidth）。不難看出，LWR每預測一個點的值都要重新得到一個新的模型。
（二）二項Logistic回歸
1. 概述
對於分類任務，一個可行的思路是把基於樣本 $\textbf x$ 計算出的連續值 $z$（比如，線性加權值 $z=\textbf w^{\top}\textbf x$ ）和離散的類別標簽值聯系起來。
二項Logistic回歸（Binomial logistic regression）是工業界應用非常廣泛的一個經典的二分類模型。一般就叫邏輯回歸，這里無意爭論到底應該怎么翻譯，雖然古人雲“名不正則言不順”，但提起“邏輯回歸”大家都知道這是哪個東西，我覺得這就夠了。對Logistic回歸的歷史感興趣的朋友們可以看一下 [7] 的介紹。首先使用logistic函數 $\sigma(\cdot)$ 將 $ z$ 從實數空間 $(-\infty,+\infty)$ 映射到概率空間 $(0,1)$ 上，可以將映射之后的值 $\sigma(z)$ 解釋為樣本 $\textbf x$ 屬於正類（記類別標記為1）的可能性，也就是后驗概率的估計值：
$$\hat y=P(y=1|\textbf x)=\sigma(z)=\frac{1}{1+\exp(-z)}$$
既然解釋成后驗概率，然后就可以給出分類規則（最大后驗概率決策）：當 $P(y=1|\textbf x)>0.5$ ，認為樣本 $\textbf x$ 屬於正類；否則樣本 $\textbf x$ 屬於正類屬於負類。
下面兩個圖是一維特征、二分類的情況。大致說了線性回歸為什么不可以用來分類。因為線性回歸輸出連續值，而類別標簽只有0、1兩個，所以需要人為設定一個閾值，將輸出值與該值比較大小，從而來判斷模型將樣本分到哪一類，而這個閾值會受到離群點（outlier）的牽制，因為線性回歸的擬合曲線會因為離群點而受到較大影響，所以不好確定；相比之下，logistic回歸不會受到圖示離群點的牽制。
圖片來源：[5]、[1]
待補充：為什么使用logistic函數歸到區間 $(0,1)$ 之后就可以解釋成概率了。
2. 決策邊界
下面說一下決策邊界。當 $P(y=1|\textbf x)=0.5$ 時，意味着 $z=0$ ，這就是決策邊界的方程。換句話說，$z$ 的形式決定了邏輯回歸的決策面是線性的還是非線性的。如果 $z=\textbf w^{\top}\textbf x$ ，那決策面當然是線性的；但是如果 $z$ 的形式並不是特征的線性組合，而是非線性的形式，當然也可以是非常復雜的決策面。
圖片來源：[1]
下面我們只討論線性決策面的情況。Logistic回歸模型為：
$$\hat y=P(y=1|\textbf x)=\sigma(\textbf w^{\top}\textbf x)=\frac{\exp(\textbf w^{\top}\textbf x)}{1+\exp(\textbf w^{\top}\textbf x)}$$
$$P(y=0|\textbf x)=1-\sigma(\textbf w^{\top}\textbf x)=\frac{1}{1+\exp(\textbf w^{\top}\textbf x)}$$
稍加變換就可以看出Logistic回歸和線性回歸的區別：
線性回歸是用 $\textbf w^{\top}\textbf x$ 去擬合 $y$ ；二項Logistic回歸則是去擬合 $\ln \dfrac{\hat y}{1-\hat y}$ ，換句話說就是在擬合對數幾率（log-odds，幾率是樣本屬於正類的可能性與屬於負類的可能性的比值）。也就是說，二項Logistic回歸在對對數幾率做回歸，進而轉化為解決分類問題。
3. 求解
（1）經驗風險最小化：極大似然估計
logistic函數 $\sigma(\cdot)$ 的導函數為 $\sigma'(\cdot)=\sigma(\cdot)\odot (\textbf 1-\sigma(\cdot))$ 。也就是說當自變量為向量時，函數對逐元素進行計算，輸出同維度的向量。
首先從經驗風險最小化的角度推導參數的求解過程。使用交叉熵損失函數（單標簽情況下就是對數損失函數），模型對一個樣本 $(\textbf x,y)$ 的對數損失為：
$$\begin{aligned}\mathcal L&=-\biggl(y\ln P(y=1|\textbf x)+(1-y)\ln P(y=0|\textbf x)\biggr)\\&=-\biggl(y\ln\frac{\exp(\textbf w^{\top}\textbf x)}{1+\exp(\textbf w^{\top}\textbf x)}+(1-y)\ln\frac{1}{1+\exp(\textbf w^{\top}\textbf x)}\biggr)\\&=-\biggl(y\textbf w^{\top}\textbf x-\ln(1+\exp(\textbf w^{\top}\textbf x))\biggr)\end{aligned}$$
所以經驗風險為：
$$\mathcal R=\frac1N\sum_{i=1}^N\mathcal L_i=-\frac1N\sum_{i=1}^N\biggl(y_i\textbf w^{\top}\textbf x_i-\ln(1+\exp(\textbf w^{\top}\textbf x_i))\biggr)$$
如果不加正則的話，優化目標為上式最小化。前面的系數 $\dfrac1N$ 也就是訓練樣本數的倒數，是定值，去掉后不影響優化目標。
從極大似然估計的角度也可以推出等價的優化目標：對數似然函數為
$$\begin{aligned}l(\textbf w)&=\ln[\prod_{i=1}^NP(y_i|\textbf x_i;\textbf w)]\\&=\ln[\prod_{i=1}^NP(y_i=1|\textbf x_i)^{y_i}P(y_i=0|\textbf x_i)^{1-y_i}]\\&=\sum_{i=1}^N\biggl(y_i\ln P(y_i=1|\textbf x_i)+(1-y_i)\ln P(y_i=0|\textbf x_i)\biggr)\end{aligned}$$
該式最大化就等價於經驗風險最小化。
由於優化目標求不出解析解，但它是高階連續可微的凸函數，所以可以用迭代的方法，如梯度下降法（GD）。
因為SGD每次迭代是隨機選擇一個樣本，所以這里先求取模型對一個樣本的損失的梯度（經驗風險的梯度無非就是加個求和號再除以訓練樣本數而已）：
$$\frac{\partial \mathcal L}{\partial\textbf w}=-\biggl(y\textbf x-\textbf x\frac{\exp(\textbf w^{\top}\textbf x)}{1+\exp(\textbf w^{\top}\textbf x)}\biggr)$$
可以發現其實它就是個特別簡單的形式：
$$\frac{\partial \mathcal L}{\partial\textbf w}=-(y-\hat y)\textbf x$$
經驗風險的梯度可以寫成矩陣的形式（ $\textbf y$ 表示訓練集樣本的真實標簽組成的列向量），省略系數 $\dfrac1N$ ：
$$\begin{aligned}\frac{\partial \mathcal R}{\partial\textbf w}&=-\sum_{i=1}^N (y_i-\sigma(\textbf w^{\top}\textbf x_i))\textbf x_i\\&= -X(\textbf y-\sigma(X^{\top}\textbf w))\\&= -X(\textbf y-\hat{\textbf y})\end{aligned}$$
可以很容易把這種批處理GD的形式改寫成mini-batch SGD的形式。
不難看出，這個梯度形式和線性回歸是一樣的。（后面會知道，Softmax回歸的梯度形式和它們也是一樣的。）
這里順便把二階梯度也求一下，可以使用牛頓法或擬牛頓法來迭代求取參數：
$$\frac{\partial \mathcal L}{\partial\textbf w}=-(y-\sigma(\textbf w^{\top}\textbf x))\textbf x=\sigma(\textbf w^{\top}\textbf x)\textbf x -y\textbf x$$
$$\begin{aligned}\frac{\partial^2\mathcal L}{\partial\textbf w^2}&=\frac{\partial(\sigma(\textbf w^{\top}\textbf x)\textbf x -y\textbf x)}{\partial\textbf w}\\&=\frac{\sigma(\textbf w^{\top}\textbf x)}{\partial\textbf w}\textbf x^{\top}\\&=\sigma(\textbf w^{\top}\textbf x)(1-\sigma(\textbf w^{\top}\textbf x))\textbf x\textbf x^{\top}\\&=\hat y(1-\hat y)\textbf x\textbf x^{\top}\end{aligned}$$
$$\frac{\partial^2\mathcal R}{\partial\textbf w^2}=\sum_{i=1}^N\hat y_i(1-\hat y_i)\textbf x_i\textbf x_i^{\top}$$
（2）結構風險最小化：最大后驗概率估計
如果將正則項加上，那就是用結構風險最小化的准則來學習參數，常用的有參數的 $L_1$ 范數（LASSO）和 $L_2$ 范數（Ridge）：
$$R_{\text{srm}}=R+\lambda\|\textbf w\|_1$$
$$R_{\text{srm}}=R+\lambda\|\textbf w\|_2^2$$
從梯度的形式來看，相比於不加正則的時候，變化量為
$$\frac{\partial \|\textbf w\|_1}{\partial w_j}=\frac{\partial \sum_{j=1}^d|w_j|}{\partial w_j}=1\text{ if }w_j>0\text{ else if }w_j<0\,-1$$
$$\frac{\partial \|\textbf w\|_2^2}{\partial w_j}=\frac{\partial \sum_{j=1}^d|w_j|^2}{\partial w_j}=2w_j$$
如果 $w_j$ 為正，那么新迭代點相比之前會減去一個正數而變小；如果 $w_j$ 為負，那么新迭代點相比之前會減去一個負數而變大。也就是說避免了特別大或者特別小的權重值出現，可以使權重的絕對值變小，從而避免過於依賴某些特征的情況，減輕過擬合。
加 $L_1$ 正則時會使某些維的參數變成0，這就是所謂的稀疏解，相當於進行了一個特征選擇的過程；加 $L_2$ 正則時權重的絕對值會變小，起到平滑化的作用。更詳細地可以參考[4]。
4. 貝葉斯角度
如果從貝葉斯估計的角度來說，正則項相當於加上了先驗知識：加 $L_1$ 正則相當於是認為參數服從Laplace分布，加 $L_2$ 正則相當於是認為參數服從均值為0、協方差為 $\frac{1}{\lambda}$ 的高斯分布。此時，結構風險最小化等價於最大后驗概率估計。具體可以參考[6][9]。
5. 與其他模型的關系
最顯然的一個就是全連接的前饋神經網絡就是多層Logistic回歸模型（不知道為什么被叫成MLP，多層感知機）。其與朴素貝葉斯的聯系可以看本系列博客第二篇。更詳細地請參考 [7] ，后面有空的話會簡單談一點。
6. 並行
這里參考[7]。
（三）Softmax回歸
Softmax回歸可以用於多類分類問題，而不像Logistic回歸等二分類模型那樣需要借助One-vs-rest。設樣本 $\textbf x_i$ 的真實類別標簽 $y_i\in\{1,2,...,C\}$ ，one-hot向量為 $\textbf y_i=(0,0,...,0,1,0,...,0)^{\top}\in\mathbb R^C$ （只有真實類別的那一維是1）。
與Logistic回歸類似，Softmax回歸輸出的是樣本 $\textbf x_i$ 屬於各個類別的后驗概率的估計值 $P(y_i=c|\textbf x_i),c\in\{1,2,...,C\}$ ：
$$z_c=\textbf w_c^{\top}\textbf x_i$$
$$\begin{aligned}P(y_i=c|\textbf x_i)&=\text{softmax}(z_c)\\&=\frac{\exp(z_c)}{\sum_{j=1}^C\exp(z_j)}\\&=\frac{\exp(\textbf w_c^{\top}\textbf x_i)}{\sum_{j=1}^C\exp(\textbf w_j^{\top}\textbf x_i)},\quad c\in\{1,2,...,C\}\end{aligned}$$
將模型對一個樣本 $(\textbf x_i,\textbf y_i)$ 的后驗概率估計組成列向量
$$\hat{\textbf y}_i=(P(y_i=1|\textbf x),P(y_i=2|\textbf x),...,P(y_i=C|\textbf x))^{\top}\in\mathbb R^C$$
並將各個類別的權重向量 $\textbf w_c\in\mathbb R^d$ 組成權重矩陣 $W=(\textbf w_1,\textbf w_2,...,\textbf w_C)\in\mathbb R^{d\times C}$ ，可以寫成如下形式：
$$\textbf z_i=(z_1,z_2,...,z_C)^{\top}=W^{\top}\textbf x_i$$
$$\begin{aligned}\hat{\textbf y}_i&=\text{softmax}(\textbf z_i)=\frac{\exp(\textbf z_i)}{\sum_{j=1}^C\exp(z_j)}\\&=\frac{\exp(W^{\top}\textbf x_i)}{\sum_{j=1}^C\exp(\textbf w_j^{\top}\textbf x)}=\frac{\exp(W^{\top}\textbf x_i)}{\textbf 1^{\top}\exp(W^{\top}\textbf x_i)}\\&=\frac{\exp(\textbf z_i)}{\textbf 1^{\top}\exp(\textbf z_i)}\end{aligned}$$
使用交叉熵損失函數，模型對一個樣本 $(\textbf x_i,\textbf y_i)$ 的損失為：
$$\mathcal L_i=-\textbf y_i^{\top}\ln\hat{\textbf y}_i=-\textbf y_i^{\top}\ln\text{softmax}(\textbf z_i)=-\textbf y_i^{\top}\ln\frac{\exp(\textbf z_i)}{\textbf 1^{\top}\exp(\textbf z_i)}$$
所以經驗風險（省略系數 $\dfrac1N$ ）為
$$\mathcal R=\sum_{i=1}^N\mathcal L_i=-\sum_{i=1}^N\biggl(\textbf y_i^{\top}\ln\text{softmax}(W^{\top}\textbf x_i)\biggr)$$
下面求取對一個樣本 $(\textbf x_i,\textbf y_i)$ 的損失的梯度。先上結論：
$$\frac{\partial \mathcal L_i}{\partial\textbf w_c}=-[\textbf y_i-\hat{\textbf y}_i]_c\textbf x_i,\quad c\in\{1,2,...,C\}$$
$$\frac{\partial \mathcal L_i}{\partial W}=-\textbf x_i(\textbf y_i-\hat{\textbf y}_i)^{\top}$$
式中 $[\textbf y_i]_c$ 表示的是向量 $\textbf y_i$ 的第 $c$ 維元素的值。
經驗風險的梯度依舊可以寫成矩陣的形式：
$$\frac{\partial \mathcal R}{\partial W}=\sum_{i=1}^N\frac{\partial \mathcal L_i}{\partial W}=-X(Y-\hat Y)$$
其中 $Y\in\mathbb R^{N\times C}$ 是one-hot標簽構成的矩陣，每一行都是一個樣本的one-hot標簽；$\hat Y$ 含義類似。
Softmax回歸的 $\dfrac{\partial \mathcal L_i}{\partial\textbf w_c}$ 的形式和Logistic回歸的 $\dfrac{\partial \mathcal L_i}{\partial\textbf w}$ 是一樣的。
下面給出推導過程。
$\dfrac{\partial \mathcal L_i}{\partial\textbf w_c}$ 的求法有三種：
（1）普通方法，一步步推，可以參考我很早之前寫的一篇講word2vec的博客，我覺得寫的還挺清楚的；
（2）[1] 中的方法，在第三章；
（3）[3] 中的方法，而且用這個方法可以直接把 $\dfrac{\partial \mathcal L_i}{\partial W}$ 求出來。
$\dfrac{\partial \mathcal L_i}{\partial W}$ 推導方式可以是用 $\dfrac{\partial \mathcal L_i}{\partial\textbf w_c}$ “拼”成對矩陣 $W$ 的梯度。下面使用 [3] 里面介紹的技巧，直接求對矩陣 $W$ 的梯度。
[3] 介紹的是這樣形式的求導：已知矩陣 $X$ ，函數 $f(X)$ 的函數值為標量，求 $\dfrac{\partial f}{\partial X}$ 。一種典型的例子就是求取損失對權重矩陣的導數。
對於一元微積分，$\text{d}f=f'(x)\text{d}x$ ；多元微積分，$\text{d}f=\sum_i\dfrac{\partial f}{\partial x_i}\text{d}x_i=(\dfrac{\partial f}{\partial \textbf x})^{\top}\text{d}\textbf x$；由此建立矩陣導數和微分的聯系：
$$\text{d}f=\sum_{i,j}\frac{\partial f}{\partial X_{ij}}\text{d}X_{ij}=\text{tr}((\frac{\partial f}{\partial X})^{\top}\text{d}X)$$
上式第二個等號成立是因為對於兩個同階方陣有 $\text{tr}(A^{\top}B)=\sum_{i,j}A_{ij}B_{ij}$ 。求解的流程就是，先求微分 $\text{d}f$ 表達式，然后再套上跡（因為標量的跡等於標量本身），然后再把表達式 $\text{tr}(\text{d}f)$ 和 $\text{tr}((\dfrac{\partial f}{\partial X})^{\top}\text{d}X)$ 進行比對，進而把 $\dfrac{\partial f}{\partial X}$ 給“挖”出來。
所以，問題就從求梯度轉化成了求微分。求微分當然少不了很多法則和技巧，下面隨着講隨着介紹。接下來就來求取Softmax回歸中的 $\dfrac{\partial \mathcal L}{\partial W}$ （樣本序號 $i$ 被省略）。
首先求取 $\text{d}\mathcal L$ 。
$$\begin{aligned}\mathcal L&=-\textbf y^{\top}\ln\frac{\exp(\textbf z)}{\textbf 1^{\top}\exp(\textbf z)}\\&=-\textbf y^{\top}(\textbf z-\ln\begin{pmatrix}\textbf 1^{\top}\exp(\textbf z) \\ \textbf 1^{\top}\exp(\textbf z) \\ \vdots \\ \textbf 1^{\top}\exp(\textbf z)\end{pmatrix})\quad \textbf 1^{\top}\exp(\textbf z)\text{是標量}\\&=\ln(\textbf 1^{\top}\exp(\textbf z))-\textbf y^{\top}\textbf z\end{aligned}$$
根據法則 $\text{d}(g(X))=g'(X)\odot\text{d}X$ 、$\text{d}(XY)=(\text{d}X)Y+X(\text{d}Y)$，可得
$$\text{d}(\ln(\textbf 1^{\top}\exp(\textbf z)))=\frac{1}{\textbf 1^{\top}\exp(\textbf z)}\odot\text{d}(\textbf 1^{\top}\exp(\textbf z))$$
$$\text{d}(\textbf 1^{\top}\exp(\textbf z))=\textbf 1^{\top}\text{d}(\exp(\textbf z))=\textbf 1^{\top}(\exp(\textbf z)\odot\text{d}\textbf z)$$
所以
$$\text{d}\mathcal L=\frac{\textbf 1^{\top}(\exp(\textbf z)\odot\text{d}\textbf z)}{\textbf 1^{\top}\exp(\textbf z)}-\textbf y^{\top}\text{d}\textbf z$$
現在可以套上跡，根據恆等式 $\text{tr}(A^{\top}(B\odot C))=\text{tr}((A\odot B)^{\top}C)=\sum_{i,j}A_{ij}B_{ij}C_{ij}$ ，可得
$$\begin{aligned}\text{d}\mathcal L&=\text{tr}(\frac{(\textbf 1\odot \exp(\textbf z))^{\top}\text{d}\textbf z}{\textbf 1^{\top}\exp(\textbf z)})-\text{tr}(\textbf y^{\top}\text{d}\textbf z)\\&=\text{tr}(\biggl(\frac{(\exp(\textbf z))^{\top}}{\textbf 1^{\top}\exp(\textbf z)}-\textbf y^{\top}\biggr)\text{d}\textbf z)\\&=\text{tr}((\hat{\textbf y}-\textbf y)^{\top}\text{d}\textbf z)\\&=\text{tr}((\frac{\partial L}{\partial\textbf z})^{\top}\text{d}\textbf z)\end{aligned}$$
現在已經成功了一半，因為已經有了 $\dfrac{\partial \mathcal L}{\partial\textbf z}$ 。因為
$$\text{d}\textbf z=\text{d}(W^{\top}\textbf x)=(\text{d}W^{\top})\textbf x+W^{\top}\text{d}\textbf x=(\text{d}W^{\top})\textbf x$$
並且 $\text{tr}(ABC)=\text{tr}(BCA)=\text{tr}(CAB)$ ，所以有
$$\begin{aligned}\text{d}\mathcal L&=\text{tr}((\frac{\partial L}{\partial\textbf z})^{\top}(\text{d}W^{\top})\textbf x)\\&=\text{tr}(\textbf x(\frac{\partial L}{\partial\textbf z})^{\top}\text{d}W^{\top})\\&=\text{tr}((\frac{\partial L}{\partial W^{\top}})^{\top}\text{d}W^{\top})\end{aligned}$$
也就是說，$\dfrac{\partial \mathcal L}{\partial W^{\top}}=\dfrac{\partial \mathcal L}{\partial\textbf z}\textbf x^{\top}=(\hat{\textbf y}-\textbf y)\textbf x^{\top}$，所以
$$\frac{\partial \mathcal L}{\partial W}=-\textbf x(\textbf y-\hat{\textbf y})^{\top}$$
（四）廣義線性模型
其實上面介紹的三種模型，都屬於廣義線性模型（Generalized linear model，GLM）。
1. 指數族分布
說到GLM，就不得不說指數族分布。設有一隨機變量 $Y$ ，觀測值為 $y$ ，那么指數族分布（Exponential family distributions）的 PDF／PMF 為如下函數：
$$p(y;\boldsymbol\eta)=b(y)\exp(\boldsymbol\eta^{\top}T(y)-a(\boldsymbol\eta))$$
式中，$\boldsymbol\eta$ 被稱為nature parameter或canonical parameter，$T(y)$ 是充分統計量（通常設 $T(y)=y$ ），$a(\boldsymbol\eta)$ 是log partition function，$\exp(-a(\boldsymbol\eta))$ 用來保證PDF的積分為1（或PMF的加和為1）。把隨機變量服從指數族分布記為 $ Y\sim ExponentialFamily(\boldsymbol\eta)$ 。伯努利分布（兩點分布）、高斯分布、多項式分布、泊松分布、指數分布、伽馬分布、貝塔分布、狄利克雷分布、維希特分布……等等都屬於指數族分布。
通過選取不同的 $\boldsymbol\eta$ ，可以得到不同的分布：
例如，對於參數為 $\phi$ 兩點分布，其PMF為
$$P(y;\phi)=\phi^y(1-\phi)^{1-y}=\exp(y\ln(\frac{\phi}{1-\phi})+\ln(1-\phi))$$
所以
$$\phi=\frac{1}{1+\exp(-\eta)}$$
這正是logistic函數。
再比如參數為均值 $\mu$ 、方差1的高斯分布，其PDF為
$$p(y;\mu)=\frac{1}{\sqrt{2\pi}}\exp(-\frac12(y-\mu)^2)=\frac{1}{\sqrt{2\pi}}\exp(-\frac12y^2)\exp(\mu y-\frac12\mu^2)$$
所以
$$\mu=\eta$$
2. 廣義線性模型
通過指數族分布，可以構建廣義線性模型。設模型的參數為 $\boldsymbol\theta$ ，對於記 $X$ 、$Y$ 分別是代表特征和標簽的隨機變量，觀測值為 $\textbf x$ 、$y$ 。首先假定如下三點：
1. 條件分布服從指數族分布，即 $Y|X;\boldsymbol\theta\sim ExponentialFamily(\boldsymbol\eta)$ 。例如，
2. 給定特征 $\textbf x$ ，目標是預測 $E[T(y)|\textbf x]$ 。因為通常設 $T(y)=y$ ，所以目標就是預測 $E[y|\textbf x]$ 。
實際上，就相當於對模型輸出值進行預測。用Logistic回歸舉例：模型輸出值為 $P(y=1|\textbf x)$ ，隨機變量 $Y$ 服從兩點分布（只可能取0、1兩個值），所以 $E[y|\textbf x;\boldsymbol\theta]=0\times P(y=0|\textbf x)+1\times P(y=1|\textbf x)=P(y=1|\textbf x)$
3. 指數族分布的參數 $\boldsymbol\eta$ 和給定特征 $\textbf x$ 的關系為線性：$\boldsymbol\eta=\boldsymbol\theta^{\top}\textbf x$
下面可以開始利用不同的 $\boldsymbol\eta$ 來構建GLM。
（1）線性回歸
對於第一個假設，設指數族分布是參數為 $\mu$ 的高斯分布，即 $\mu=\eta$ ；那么對於第二個假設，可知模型輸出值為 $\mu$ ，結合第一個假設可知模型輸出值為 $\eta$ ；根據第三個假設 $\boldsymbol\eta=\boldsymbol\theta^{\top}\textbf x$ ，可知模型輸出值為 $\boldsymbol\theta^{\top}\textbf x$ 。這就推導出了線性回歸模型。
（2）Logistic回歸
對於第一個假設，設指數族分布是參數為 $\phi$ 的伯努利分布，即 $\phi=\frac{1}{1+\exp(-\eta)}$ ；那么對於第二個假設，因為伯努利分布的期望為 $\phi$ ，可知模型輸出值為 $\phi$ ；根據第三個假設 $\boldsymbol\eta=\boldsymbol\theta^{\top}\textbf x$ ，可知模型輸出值為 $\phi=\frac{1}{1+\exp(-\boldsymbol\theta^{\top}\textbf x)}$ 。這就推導出了Logistic回歸模型。
（3）Softmax回歸
相應的指數族分布是多項式分布，代表標簽的是一個隨機向量 $\boldsymbol Y$ 。詳細的推導這里就不贅述了，可以參考 [5] 的最后一部分。
參考資料：
[1] 《神經網絡與深度學習講義》
[2] 《統計學習方法》
[3] 《矩陣求導術（上）》
[5] CS229 Lecture Notes1、Lecture6 slides
[6] Regularized Regression: A Bayesian point of view
[7] 淺析Logistic Regression （寫的比我這篇真的好太多了。。。）
[8] Matrix_calculus
[9] https://www.evernote.com/shard/s146/sh/bf0d0a08-d5c5-4dc9-b70f-fa6374b9ceae/14e34bc9e8e518b1a5cc17fc585d75fc