機器學習-softmax 回歸原理與實作

文章目錄

• 一、什麼是 softmax 回歸？
• 二、原理
• • 2.1 梯度下降法參數求解
  • 2.2 模型參數特點
  • 2.3 正則化
  • 2.4 softmax 與 logistic 回歸的關系
• 三、實作
• • 3.1 python 手動實作
  • 3.2 sklearn 算法包實作
• 四、參考

一、什麼是 softmax 回歸？

softmax 回歸(softmax regression)其實是 logistic 回歸的一般形式，logistic 回歸用于二分類，而 softmax 回歸用于多分類，關于 logistic 回歸可以看我的這篇部落格👉[機器學習-logistic回歸原理與實作]。

對于輸入資料 { ( x 1 , y 1 ) , ( x 2 , y 2 ) , … , ( x m , y m ) } \{(x_1,y_1),(x_2,y_2),\ldots,(x_m,y_m)\} {(x1​,y1​),(x2​,y2​),…,(xm​,ym​)}有 k k k 個類别，即 y i ∈ { 1 , 2 , … , k } y_i \in \{1,2,\ldots,k\} yi​∈{1,2,…,k}，那麼 softmax 回歸主要估算輸入資料 x i x_i xi​ 歸屬于每一類的機率，即

h θ ( x i ) = [ p ( y i = 1 ∣ x i ; θ ) p ( y i = 2 ∣ x i ; θ ) ⋮ p ( y i = k ∣ x i ; θ ) ] = 1 ∑ j = 1 k e θ j T x i [ e θ 1 T x i e θ 2 T x i ⋮ e θ k T x i ] (1) h_{\theta}\left(x_i\right)=\left[\begin{array}{c}{p\left(y_i=1 | x_i ; \theta\right)} \\ {p\left(y_i=2 | x_i ; \theta\right)} \\ {\vdots} \\ {p\left(y_i=k | x_i ; \theta\right)}\end{array}\right]=\frac{1}{\sum_{j=1}^{k} e^{\theta_{j}^{T} x_i}}\left[\begin{array}{c}{e^{\theta_{1}^{T} x_i}} \\ {e^{\theta_{2}^{T} x_i}} \\ {\vdots} \\ {e^{\theta_{k}^{T} x_i}}\end{array}\right]\tag{1} hθ​(xi​)=⎣⎢⎢⎢⎡​p(yi​=1∣xi​;θ)p(yi​=2∣xi​;θ)⋮p(yi​=k∣xi​;θ)​⎦⎥⎥⎥⎤​=∑j=1k​eθjT​xi​1​⎣⎢⎢⎢⎢⎡​eθ1T​xi​eθ2T​xi​⋮eθkT​xi​​⎦⎥⎥⎥⎥⎤​(1)

其中， θ 1 , θ 2 , … , θ k ∈ θ \theta_1,\theta_2,\ldots,\theta_k \in \theta θ1​,θ2​,…,θk​∈θ是模型的參數，乘以 1 ∑ j = 1 k e θ j T x i \frac{1}{\sum_{j=1}^{k} e^{\theta_{j}^{T} x_i}} ∑j=1k​eθjT​xi​1​是為了讓機率位于[0,1]并且機率之和為 1，softmax 回歸将輸入資料 x i x_i xi​ 歸屬于類别 j j j 的機率為

p ( y i = j ∣ x i ; θ ) = e θ j T x i ∑ l = 1 k e θ l T x i (2) p\left(y_i=j | x_i ; \theta\right)=\frac{e^{\theta_{j}^{T} x_i}}{\sum_{l=1}^{k} e^{\theta_{l}^{T} x_i}}\tag{2} p(yi​=j∣xi​;θ)=∑l=1k​eθlT​xi​eθjT​xi​​(2)

上面的式子可以用下圖形象化的解析(來自台大李宏毅《一天搞懂深度學習》)。

二、原理

2.1 梯度下降法參數求解

softmax 回歸的參數矩陣 θ \theta θ 可以記為

θ = [ θ 1 T θ 2 T ⋮ θ k T ] (3) \theta=\left[\begin{array}{c}{\theta_{1}^{T}} \\ {\theta_{2}^{T}} \\ {\vdots} \\ {\theta_{k}^{T}}\end{array}\right]\tag{3} θ=⎣⎢⎢⎢⎡​θ1T​θ2T​⋮θkT​​⎦⎥⎥⎥⎤​(3)

定義 softmax 回歸的代價函數

L ( θ ) = − 1 m [ ∑ i = 1 m ∑ j = 1 k 1 { y i = j } log ⁡ e θ j T x i ∑ l = 1 k e θ l T x i ] (4) L(\theta)=-\frac{1}{m}\left[\sum_{i=1}^{m} \sum_{j=1}^{k} 1\left\{y_i=j\right\} \log \frac{e^{\theta_{j}^{T} x_i}}{\sum_{l=1}^{k} e^{\theta_{l}^{T} x_i}}\right]\tag{4} L(θ)=−m1​[i=1∑m​j=1∑k​1{yi​=j}log∑l=1k​eθlT​xi​eθjT​xi​​](4)

其中，1{·}是示性函數，即1{值為真的表達式}=1，1{值為假的表達式}=0。跟 logistic 函數一樣，利用梯度下降法最小化代價函數，下面求解 θ \theta θ 的梯度。 L ( θ ) L(\theta) L(θ)關于 θ j \theta_{j} θj​ 的梯度求解為

∂ L ( θ ) ∂ θ j = − 1 m ∂ ∂ θ j [ ∑ i = 1 m ∑ j = 1 k 1 { y i = j } log ⁡ e θ j T x i ∑ l = 1 k e θ l T x i ] = − 1 m ∂ ∂ θ j [ ∑ i = 1 m ∑ j = 1 k 1 { y i = j } ( θ j T x i − log ⁡ ∑ l = 1 k e θ l T x i ) ] = − 1 m [ ∑ i = 1 m 1 { y i = j } ( x i − ∑ j = 1 k e θ j T x i ⋅ x i ∑ l = 1 k e θ l T x i ) ] = − 1 m [ ∑ i = 1 m x i 1 { y i = j } ( 1 − ∑ j = 1 k e θ j T x i ∑ l = 1 k e θ l T x i ) ] = − 1 m [ ∑ i = 1 m x i ( 1 { y i = j } − ∑ j = 1 k 1 { y i = j } e θ j T x i ∑ l = 1 k e θ l T x i ) ] = − 1 m [ ∑ i = 1 m x i ( 1 { y i = j } − e θ j T x i ∑ l = 1 k e θ l T x i ) ] = − 1 m [ ∑ i = 1 m x i ( 1 { y i = j } − p ( y i = j ∣ x i ; θ ) ) ] (5) \begin{aligned} \frac{\partial L(\theta)}{\partial \theta_{j}} &=-\frac{1}{m} \frac{\partial}{\partial \theta_{j}}\left[\sum_{i=1}^{m} \sum_{j=1}^{k} 1\left\{y_i=j\right\} \log \frac{e^{\theta_{j}^{T} x_i}}{\sum_{l=1}^{k} e^{\theta_{l}^{T} x_i}}\right] \\ &=-\frac{1}{m} \frac{\partial}{\partial \theta_{j}}\left[\sum_{i=1}^{m} \sum_{j=1}^{k} 1\left\{y_i=j\right\}\left(\theta_{j}^{T} x_i-\log \sum_{l=1}^{k} e^{\theta_{l}^{T} x_i}\right)\right] \\ &=-\frac{1}{m}\left[\sum_{i=1}^{m} 1\left\{y_i=j\right\}\left(x_i-\sum_{j=1}^{k} \frac{e^{\theta_{j}^{T} x_i} \cdot x_i}{\sum_{l=1}^{k} e^{\theta_{l}^{T} x_i}}\right)\right] \\ &=-\frac{1}{m}\left[\sum_{i=1}^{m} x_i1\left\{y_i=j\right\}\left(1-\sum_{j=1}^{k} \frac{e^{\theta_{j}^{T} x_i}}{\sum_{l=1}^{k} e^{\theta_{l}^{T} x_i}}\right)\right] \\ &=-\frac{1}{m}\left[\sum_{i=1}^{m} x_i\left(1\left\{y_i=j\right\}-\sum_{j=1}^{k} 1\left\{y_i=j\right\} \frac{e^{\theta_{j}^{T} x_i}}{\sum_{l=1}^{k} e^{\theta_{l}^{T} x_i}}\right)\right] \\ &=-\frac{1}{m}\left[\sum_{i=1}^{m} x_i\left(1\left\{y_i=j\right\}- \frac{e^{\theta_{j}^{T} x_i}}{\sum_{l=1}^{k} e^{\theta_{l}^{T} x_i}}\right)\right] \\ &=-\frac{1}{m}\left[\sum_{i=1}^{m} x_i\left(1\left\{y_i=j\right\}-p\left(y_i=j | x_i ; \theta\right)\right)\right] \end{aligned}\tag{5} ∂θj​∂L(θ)​​=−m1​∂θj​∂​[i=1∑m​j=1∑k​1{yi​=j}log∑l=1k​eθlT​xi​eθjT​xi​​]=−m1​∂θj​∂​[i=1∑m​j=1∑k​1{yi​=j}(θjT​xi​−logl=1∑k​eθlT​xi​)]=−m1​[i=1∑m​1{yi​=j}(xi​−j=1∑k​∑l=1k​eθlT​xi​eθjT​xi​⋅xi​​)]=−m1​[i=1∑m​xi​1{yi​=j}(1−j=1∑k​∑l=1k​eθlT​xi​eθjT​xi​​)]=−m1​[i=1∑m​xi​(1{yi​=j}−j=1∑k​1{yi​=j}∑l=1k​eθlT​xi​eθjT​xi​​)]=−m1​[i=1∑m​xi​(1{yi​=j}−∑l=1k​eθlT​xi​eθjT​xi​​)]=−m1​[i=1∑m​xi​(1{yi​=j}−p(yi​=j∣xi​;θ))]​(5)

感謝 CSDN 部落客[2]提供了另外一種求解方法，具體如下

∂ L ( θ ) ∂ θ j = − 1 m [ ∑ i = 1 m ∂ ∂ θ j ( 1 { y i = j } log ⁡ e θ j T x i ∑ l = 1 k e θ l T x i + ∑ c ≠ j k 1 { y i = c } log ⁡ e θ c T x i ∑ l = 1 k e θ l T x i ) ] = − 1 m [ ∑ i = 1 m ( 1 { y i = j } ( x i − e θ j T x i ⋅ x i ∑ l = 1 k e θ l T x i ) + ∑ c ≠ j k 1 { y i = c } ( − e θ j T x i ⋅ x i ∑ l = 1 k e θ l T x i ) ) ] = − 1 m [ ∑ i = 1 m x i ( 1 { y i = j } ( 1 − e θ j T x i ∑ l = 1 k e θ l T x i ) − ∑ c ≠ j k 1 { y i = c } e θ j T x i ∑ l = 1 k e θ l T x i ) ] = − 1 m [ ∑ i = 1 m x i ( 1 { y i = j } − 1 { y i = j } p ( y i = j ∣ x i ; θ ) − ∑ c ≠ j k 1 { y i = c } p ( y i = j ∣ x i ; θ ) ) ] = − 1 m [ ∑ i = 1 m x i ( 1 { y i = j } − ∑ j = 1 k 1 { y i = j } p ( y i = j ∣ x i ; θ ) ) ] = − 1 m [ ∑ i = 1 m x i ( 1 { y i = j } − p ( y i = j ∣ x i ; θ ) ) ] (6) \begin{aligned} \frac{\partial L(\theta)}{\partial \theta_{j}} &=-\frac{1}{m}\left[\sum_{i=1}^{m} \frac{\partial}{\partial \theta_{j}}\left(1\left\{y_i=j\right\} \log \frac{e^{\theta_{j}^{T} x_i}}{\sum_{l=1}^{k} e^{\theta_{l}^{T} x_i}}+\sum_{c \neq j}^{k} 1\left\{y_i=c\right\} \log \frac{e^{\theta_{c}^{T} x_i}}{\sum_{l=1}^{k} e^{\theta_{l}^{T} x_i}}\right)\right]\\ &=-\frac{1}{m}\left[\sum_{i=1}^{m}\left(1\left\{y_i=j\right\}\left(x_i-\frac{e^{\theta_{j}^{T} x_i} \cdot x_i}{\sum_{l=1}^{k} e^{\theta_{l}^{T} x_i}}\right)+\sum_{c \neq j}^{k} 1\left\{y_i=c\right\}\left(-\frac{e^{\theta_{j}^{T} x_i} \cdot x_i}{\sum_{l=1}^{k} e^{\theta_{l}^{T} x_i}}\right)\right)\right]\\ &=-\frac{1}{m}\left[\sum_{i=1}^{m} x_i\left(1\left\{y_i=j\right\}\left(1-\frac{e^{\theta_{j}^{T} x_i}}{\sum_{l=1}^{k} e^{\theta_{l}^{T} x_i}}\right)-\sum_{c \neq j}^{k} 1\left\{y_i=c\right\}\frac{e^{\theta_{j}^{T} x_i}}{\sum_{l=1}^{k} e^{\theta_{l}^{T} x_i}}\right)\right] \\ &=-\frac{1}{m}\left[\sum_{i=1}^{m} x_i\left(1\left\{y_i=j\right\}-1\left\{y_i=j\right\} p\left(y_i=j | x_i ; \theta\right)-\sum_{c \neq j}^{k} 1\left\{y_i=c\right\}p\left(y_i=j | x_i ; \theta\right)\right)\right] \\ &=-\frac{1}{m}\left[\sum_{i=1}^{m} x_i\left(1\left\{y_i=j\right\}-\sum_{j=1}^{k} 1\left\{y_i=j\right\} p\left(y_i=j | x_i ; \theta\right)\right)\right] \\ &=-\frac{1}{m}\left[\sum_{i=1}^{m} x_i\left(1\left\{y_i=j\right\}-p\left(y_i=j | x_i ; \theta\right)\right)\right] \end{aligned}\tag{6} ∂θj​∂L(θ)​​=−m1​⎣⎡​i=1∑m​∂θj​∂​⎝⎛​1{yi​=j}log∑l=1k​eθlT​xi​eθjT​xi​​+c​=j∑k​1{yi​=c}log∑l=1k​eθlT​xi​eθcT​xi​​⎠⎞​⎦⎤​=−m1​⎣⎡​i=1∑m​⎝⎛​1{yi​=j}(xi​−∑l=1k​eθlT​xi​eθjT​xi​⋅xi​​)+c​=j∑k​1{yi​=c}(−∑l=1k​eθlT​xi​eθjT​xi​⋅xi​​)⎠⎞​⎦⎤​=−m1​⎣⎡​i=1∑m​xi​⎝⎛​1{yi​=j}(1−∑l=1k​eθlT​xi​eθjT​xi​​)−c​=j∑k​1{yi​=c}∑l=1k​eθlT​xi​eθjT​xi​​⎠⎞​⎦⎤​=−m1​⎣⎡​i=1∑m​xi​⎝⎛​1{yi​=j}−1{yi​=j}p(yi​=j∣xi​;θ)−c​=j∑k​1{yi​=c}p(yi​=j∣xi​;θ)⎠⎞​⎦⎤​=−m1​[i=1∑m​xi​(1{yi​=j}−j=1∑k​1{yi​=j}p(yi​=j∣xi​;θ))]=−m1​[i=1∑m​xi​(1{yi​=j}−p(yi​=j∣xi​;θ))]​(6)

2.2 模型參數特點

softmax 回歸有一個不尋常的特點：它有一個“備援“的參數集。為了便于闡述這一特點，假設我們從參數向量 θ j \theta_{j} θj​ 中減去向量 ψ \psi ψ ，那麼對于機率函數，我們有

( y i = j ∣ x i ; θ ) = e ( θ j − ψ ) T x i ∑ l = 1 k e ( θ l − ψ ) T x i = e θ j T x i e − ψ T x i ∑ l = 1 k e θ l T x i e − ψ T x i = e θ j T x i ∑ l = 1 k e θ l T x i (7) \begin{aligned} \left(y_i=j | x_i ; \theta\right) &=\frac{e^{\left(\theta_{j}- \psi\right)^{T} x_i}}{\sum_{l=1}^{k} e^{\left(\theta_{l}- \psi\right)^{T} x_i}}\\ &=\frac{e^{\theta_{j}^{T} x_i }e^{-\psi^{T} x_i }}{\sum_{l=1}^{k} e^{\theta_{l}^{T} x_i}e^{-\psi^{T} x_i}}\\ &=\frac{e^{\theta_{j}^{T} x_i}}{\sum_{l=1}^{k} e^{\theta_{l}^{T} x_i}} \end{aligned} \tag{7} (yi​=j∣xi​;θ)​=∑l=1k​e(θl​−ψ)Txi​e(θj​−ψ)Txi​​=∑l=1k​eθlT​xi​e−ψTxi​eθjT​xi​e−ψTxi​​=∑l=1k​eθlT​xi​eθjT​xi​​​(7)

換句話說，從參數向量中的每個元素 θ j \theta_j θj​ 中減去 ψ \psi ψ 一點也不會影響到假設的類别預測！這表明了 softmax 回歸的參數中是有多餘的。正式地說， softmax 模型是過參數化的（ overparameterized​ 或參數備援的），這意味着對任何一個拟合資料的假設而言，多種參數取值有可能得到同樣的假設 h θ h_\theta hθ​，即從輸入 x x x 經過不同的模型參數的假設計算進而得到同樣的分類預測結果。

進一步說，若代價函數 L ( θ ) L(\theta) L(θ) 被某組模型參數 ( θ 1 , θ 2 , … , θ k ) (\theta_1, \theta_2,\ldots, \theta_k) (θ1​,θ2​,…,θk​) 最小化，那麼對任意的 ψ \psi ψ ，代價函數也可以被 ( θ 1 − ψ , θ 2 − ψ , … , θ k − ψ ) (\theta_1 - \psi, \theta_2 - \psi,\ldots, \theta_k - \psi) (θ1​−ψ,θ2​−ψ,…,θk​−ψ) 最小化。是以， L ( θ ) L(\theta) L(θ) 的最小值時的參數并不唯一。（有趣的是， L ( θ ) L(\theta) L(θ) 仍是凸的，并且在梯度下降中不會遇到局部最優的問題，但是 Hessian​ 矩陣是奇異或不可逆的，這将會導緻在牛頓法的直接實作上遇到數值問題。）

注意到，通過設定 ψ = θ k \psi = \theta_k ψ=θk​ ，總是可以用 θ k − ψ = 0 ⃗ \theta_k - \psi = \vec{0} θk​−ψ=0

代替 θ k \theta_k θk​ ，而不會對假設函數有任何影響。是以，可以去掉參數向量 θ \theta θ 中的最後一個（或該向量中任意其它任意一個）元素 θ k \theta_{k} θk​ ，而不影響假設函數的表達能力。實際上，因參數備援的特性，與其優化全部的 k ⋅ n k\cdot n k⋅n 個參數 ( θ 1 , θ 2 , … , θ k ) (\theta_1,\theta_2,\ldots,\theta_k) (θ1​,θ2​,…,θk​) （其中 θ k ∈ ℜ n \theta_k \in \Re^{n} θk​∈ℜn），也可令 θ k = 0 ⃗ \theta_k = \vec{0} θk​=0

，隻優化剩餘的 ( k − 1 ) ⋅ n (k-1) \cdot n (k−1)⋅n 個參數，算法依然能夠正常工作。

2.3 正則化

當訓練資料不夠多的時候，容易出現過拟合現象，拟合系數往往非常大👉[過拟合原因]，為此在損失函數後面加上一個正則項，即

L ( θ ) = − 1 m [ ∑ i = 1 m ∑ j = 1 k 1 { y i = j } log ⁡ e θ j T x i ∑ l = 1 k e θ l T x i ] + λ ∑ i = 1 k ∑ j = 1 n θ i j 2 (8) L(\theta)=-\frac{1}{m}\left[\sum_{i=1}^{m} \sum_{j=1}^{k} 1\left\{y_i=j\right\} \log \frac{e^{\theta_{j}^{T} x_i}}{\sum_{l=1}^{k} e^{\theta_{l}^{T} x_i}}\right]+\lambda\sum_{i=1}^{k} \sum_{j=1}^{n}\theta_{ij}^{2} \tag{8} L(θ)=−m1​[i=1∑m​j=1∑k​1{yi​=j}log∑l=1k​eθlT​xi​eθjT​xi​​]+λi=1∑k​j=1∑n​θij2​(8)

那麼新的損失函數的梯度為

∂ L ( θ ) ∂ θ j = − 1 m [ ∑ i = 1 m x i ( 1 { y i = j } − p ( y i = j ∣ x i ; θ ) ) ] + λ θ j (9) \frac{\partial L(\theta)}{\partial \theta_{j}} =-\frac{1}{m}\left[\sum_{i=1}^{m} x_i\left(1\left\{y_i=j\right\}-p\left(y_i=j | x_i ; \theta\right)\right)\right]+\lambda\theta_j \tag{9} ∂θj​∂L(θ)​=−m1​[i=1∑m​xi​(1{yi​=j}−p(yi​=j∣xi​;θ))]+λθj​(9)

⚠️注意：上式中的 θ j \theta_j θj​ 中的 θ j 0 \theta_{j0} θj0​ 不應該被懲罰，因為他是一個常數項，是以在實際使用的時候僅僅需要對 θ j 1 , θ j 2 , … , θ j n \theta_{j1},\theta_{j2},\dots,\theta_{jn} θj1​,θj2​,…,θjn​ 進行懲罰即可，這個會在後面的 python 代碼中提到😃。

2.4 softmax 與 logistic 回歸的關系

文章開頭說過，softmax 回歸是 logistic 回歸的一般形式，logistic 回歸是 softmax 回歸在 k = 2 k=2 k=2 時的特殊形式，下面通過公式推導來看下當 k = 2 k=2 k=2 時 softmax 回歸是如何退化成 logistic 回歸。

當 k = 2 k=2 k=2 時，softmax 回歸的假設函數為

h θ ( x i ) = [ p ( y i = 1 ∣ x i ; θ ) p ( y i = 2 ∣ x i ; θ ) ] = 1 e θ 1 T x i + e θ 2 T x i [ e θ 1 T x i e θ 2 T x i ] (10) h_{\theta}\left(x_i\right)=\left[\begin{array}{c}{p\left(y_i=1 | x_i ; \theta\right)} \\ {p\left(y_i=2 | x_i ; \theta\right)} \end{array}\right]=\frac{1}{e^{\theta_{1}^{T} x_i}+e^{\theta_{2}^{T} x_i}}\left[\begin{array}{c}{e^{\theta_{1}^{T} x_i}} \\ {e^{\theta_{2}^{T} x_i}} \end{array}\right]\tag{10} hθ​(xi​)=[p(yi​=1∣xi​;θ)p(yi​=2∣xi​;θ)​]=eθ1T​xi​+eθ2T​xi​1​[eθ1T​xi​eθ2T​xi​​](10)

前面說過 softmax 回歸的參數具有備援性，從參數向量 θ 1 , θ 2 \theta_1,\theta_2 θ1​,θ2​ 中減去向量 θ 1 \theta_1 θ1​完全不影響結果。現在我們令 θ ′ = θ 2 − θ 1 \theta'=\theta_2-\theta_1 θ′=θ2​−θ1​，并且兩個參數向量都減去 θ 1 \theta_1 θ1​，則有

h θ ( x i ) = 1 e 0 ⃗ T x i + e ( θ 2 − θ 1 ) T x i [ e 0 ⃗ T x i e ( θ 2 − θ 1 ) T x i ] = [ 1 1 + e ( θ 2 − θ 1 ) T x i e ( θ 2 − θ 1 ) T x i 1 + e ( θ 2 − θ 1 ) T x i ] = [ 1 1 + e ( θ 2 − θ 1 ) T x i 1 − 1 1 + e ( θ 2 − θ 1 ) T x i ] = [ 1 1 + e ( θ ′ ) T x i 1 − 1 1 + e ( θ ′ ) T x i ] (11) \begin{aligned} h_{\theta}(x_i) &=\frac{1}{e^{\vec{0}^{T} x_i}+e^{\left(\theta_{2}-\theta_{1}\right)^{T} x_i}}\left[\begin{array}{c}{e^{\vec{0}^{T} x_i}} \\ {e^{\left(\theta_{2}-\theta_{1}\right)^{T} x_i}}\end{array}\right] \\ &=\left[\begin{array}{c}{\frac{1}{1+e^{\left(\theta_{2}-\theta_{1}\right)^{T} x_i}}} \\ {\frac{e^{\left(\theta_{2}-\theta_{1}\right)^{T} x_i}}{1+e^{\left(\theta_{2}-\theta_{1}\right)^{T} x_i}}}\end{array}\right] \\ &=\left[\begin{array}{c}{\frac{1}{1+e^{\left(\theta_{2}-\theta_{1}\right)^{T} x_i}}} \\ 1-{\frac{1}{1+e^{\left(\theta_{2}-\theta_{1}\right)^{T} x_i}}}\end{array}\right] \\ &=\left[\begin{array}{c}{\frac{1}{1+e^{\left(\theta'\right)^{T} x_i}}} \\ 1-{\frac{1}{1+e^{\left(\theta'\right)^{T} x_i}}}\end{array}\right] \\ \end{aligned}\tag{11} hθ​(xi​)​=e0

Txi​+e(θ2​−θ1​)Txi​1​[e0

Txi​e(θ2​−θ1​)Txi​​]=[1+e(θ2​−θ1​)Txi​1​1+e(θ2​−θ1​)Txi​e(θ2​−θ1​)Txi​​​]=[1+e(θ2​−θ1​)Txi​1​1−1+e(θ2​−θ1​)Txi​1​​]=[1+e(θ′)Txi​1​1−1+e(θ′)Txi​1​​]​(11)

這樣就化成了 logistic 回歸。

三、實作

3.1 python 手動實作

這裡的資料使用的是 sklearn 的算法包生成的随機資料，其中，訓練資料為 3750×2 的資料，測試資料為 1250×2 的資料，生成代碼如下

def gen_dataset():
    from sklearn.datasets import load_iris
    from sklearn.model_selection import train_test_split
    from sklearn.datasets import make_blobs
    import matplotlib.pyplot as plt
    np.random.seed(13)
    X, y = make_blobs(centers=4, n_samples = 5000)
    # 繪制資料分布
    plt.figure(figsize=(6,4))
    plt.scatter(X[:,0], X[:,1],c=y)
    plt.title("Dataset")
    plt.xlabel("First feature")
    plt.ylabel("Second feature")
    plt.show()

    # 重塑目标以獲得具有 (n_samples, 1)形狀的列向量
    y = y.reshape((-1,1))
    # 分割資料集
    X_train, X_test, y_train, y_test = train_test_split(X, y)
    train_dataset = np.append(X_train,y_train, axis = 1)
    test_dataset = np.append(X_test,y_test, axis = 1)
    np.savetxt("train_dataset.txt", train_dataset, fmt="%.4f %.4f %d")
    np.savetxt("test_dataset.txt", test_dataset, fmt="%.4f %.4f %d")
           

資料分布情況如下圖所示

softmax 算法的核心部分就是求解梯度矩陣，我們設輸入資料為 X = { x 1 , x 2 , … , x m } X=\{x_1,x_2,\ldots,x_m\} X={x1​,x2​,…,xm​}，這是一個 m × n m×n m×n 的矩陣，輸出類别為 y = { y 1 , y 2 , … , y m } y=\{y_1,y_2,\ldots,y_m\} y={y1​,y2​,…,ym​}，其中 y i y_i yi​ 是一個 1 × k 1×k 1×k 的one-hot 矩陣， k k k 表示類别個數，那麼 y y y 其實是一個 m × k m×k m×k 的矩陣，輸入資料對應的機率為 P = { p 1 , p 2 , … , p m } P=\{p_1,p_2,\ldots,p_m\} P={p1​,p2​,…,pm​}， 同樣的這也是一個 m × k m×k m×k 的矩陣。那麼根據公式(9)，可以知道 θ j \theta_j θj​ 的梯度為

∂ L ( θ ) ∂ θ j = − 1 m ( y i − P i ) T X + λ θ j (12) \frac{\partial L(\theta)}{\partial \theta_{j}} =-\frac{1}{m}\left(y_i-P_i\right)^TX+\lambda\theta_j \tag{12} ∂θj​∂L(θ)​=−m1​(yi​−Pi​)TX+λθj​(12)

由此可以推導出 θ \theta θ 的參數矩陣為

∂ L ( θ ) ∂ θ = − 1 m ( y − P ) T X + λ θ (13) \frac{\partial L(\theta)}{\partial \theta} =-\frac{1}{m}\left(y-P\right)^TX+\lambda\theta \tag{13} ∂θ∂L(θ)​=−m1​(y−P)TX+λθ(13)

注意到這裡也考慮了 θ j \theta_j θj​ 的第 0 項 ，是以在寫代碼的時候需要把 θ \theta θ 的第 0 列的懲罰項減去。

softmax 回歸的代碼如下

def load_dataset(file_path):
    dataMat = []
    labelMat = []
    fr = open(file_path)
    for line in fr.readlines():
        lineArr = line.strip().split()
        dataMat.append([1.0, float(lineArr[0]), float(lineArr[1])])
        labelMat.append(int(lineArr[2]))
    return dataMat, labelMat


def train(data_arr, label_arr, n_class, iters = 1000, alpha = 0.1, lam = 0.01):
    '''
    @description: softmax 訓練函數
    @param {type} 
    @return: theta 參數
    '''    
    n_samples, n_features = data_arr.shape
    n_classes = n_class
    # 随機初始化權重矩陣
    weights = np.random.rand(n_class, n_features)
    # 定義損失結果
    all_loss = list()
    # 計算 one-hot 矩陣
    y_one_hot = one_hot(label_arr, n_samples, n_classes)
    for i in range(iters):
        # 計算 m * k 的分數矩陣
        scores = np.dot(data_arr, weights.T)
        # 計算 softmax 的值
        probs = softmax(scores)
        # 計算損失函數值
        loss = - (1.0 / n_samples) * np.sum(y_one_hot * np.log(probs))
        all_loss.append(loss)
        # 求解梯度
        dw = -(1.0 / n_samples) * np.dot((y_one_hot - probs).T, data_arr) + lam * weights
        dw[:,0] = dw[:,0] - lam * weights[:,0]
        # 更新權重矩陣
        weights  = weights - alpha * dw
    return weights, all_loss
        

def softmax(scores):
    # 計算總和
    sum_exp = np.sum(np.exp(scores), axis = 1,keepdims = True)
    softmax = np.exp(scores) / sum_exp
    return softmax


def one_hot(label_arr, n_samples, n_classes):
    one_hot = np.zeros((n_samples, n_classes))
    one_hot[np.arange(n_samples), label_arr.T] = 1
    return one_hot


def predict(test_dataset, label_arr, weights):
    scores = np.dot(test_dataset, weights.T)
    probs = softmax(scores)
    return np.argmax(probs, axis=1).reshape((-1,1))


if __name__ == "__main__":
    #gen_dataset()
    data_arr, label_arr = load_dataset('train_dataset.txt')
    data_arr = np.array(data_arr)
    label_arr = np.array(label_arr).reshape((-1,1))
    weights, all_loss = train(data_arr, label_arr, n_class = 4)

    # 計算預測的準确率
    test_data_arr, test_label_arr = load_dataset('test_dataset.txt')
    test_data_arr = np.array(test_data_arr)
    test_label_arr = np.array(test_label_arr).reshape((-1,1))
    n_test_samples = test_data_arr.shape[0]
    y_predict = predict(test_data_arr, test_label_arr, weights)
    accuray = np.sum(y_predict == test_label_arr) / n_test_samples
    print(accuray)

    # 繪制損失函數
    fig = plt.figure(figsize=(8,5))
    plt.plot(np.arange(1000), all_loss)
    plt.title("Development of loss during training")
    plt.xlabel("Number of iterations")
    plt.ylabel("Loss")
    plt.show()
           

函數輸出的測試資料準确率為

0.9952
           

程式中記錄了每個循環的損失函數，其變化曲線如下圖所示。

3.2 sklearn 算法包實作

sklearn

的實作比較簡單，與 logistic 回歸的代碼類似。

def softmax_lib():
    data_arr, label_arr = load_dataset('train_dataset.txt')
    from sklearn import linear_model
    model_softmax_regression = linear_model.LogisticRegression(solver='lbfgs',multi_class="multinomial",max_iter=10)
    model_softmax_regression.fit(data_arr, label_arr)
    test_data_arr, test_label_arr = load_dataset('test_dataset.txt')
    y_predict = model_softmax_regression.predict(test_data_arr)
    accurcy = np.sum(y_predict == test_label_arr) / len(test_data_arr)
    print(accurcy)
           

輸出結果為

0.9848
           

本文的完整代碼和資料去👉[我的 github]檢視

四、參考

[1] https://zhuanlan.zhihu.com/p/34520042

[2] https://blog.csdn.net/u012328159/article/details/72155874

[3] https://zhuanlan.zhihu.com/p/56139075

[4] https://zh.wikipedia.org/wiki/Softmax%E5%87%BD%E6%95%B0

[5] http://deeplearning.stanford.edu/tutorial/supervised/SoftmaxRegression/