假設你想要建立一個模型，根據某特征\(x\)，例如商品促銷活動，近期廣告，天氣等來預測給定時間內顧客到達商場的數量\(y\),我們知道泊松分布可以很好的描述這個問題。那么我們怎樣來建立這個問題的回歸模型呢？幸運的是泊松分布是指數族分布，所以我們可以使用廣義線性回歸模型（GLM），本文將介紹廣義線性模型來解決這個問題。
更一般的，在考慮回歸和分類問題，我們需要考慮在特征\(x\)下\(y\)的值，為了導出\(GLM\)模型，我們將會給出3個假設：
-
\(y|x;\theta \sim ExponentialFamily(η)\),給出定\(\theta\),\(y|x\)服從指數族分布，並以\(\eta\)為參數
-
給定\(x\)，我們的目標是預測\(T(y)\)的期望值，在大多數例子里，我們有\(T(y)=y\),這就意味着我們學習的輸出\(h(x)=E[y|x]\)。例如在邏輯回歸中，我們有\(h_\theta(x)=p(y=0|x) \cdot 0+p(y=1|x) \cdot 1=E[y|x;\theta]\).
-
參數\(\eta\)與輸入\(x\)是線性關系\(\eta = \theta^Tx\)（如果\(\eta\)是一個向量，則\(\eta_i=\theta^Tx\)）.
上面第三條不像一個假設，更像一個約定，可以認為是“設計的假設”。這三個假設能讓我們推出\(GLM\)模型，具這個模型有許多不錯的特性，例如易於學習等。我們很快會發現，邏輯回歸和最小二乘模型都可以作為\(GLM\)推導出來。
一、指數分布族介紹
指數分布族是指可以表示為指數形式的概率分布。指數分布的形式如下：
\[p(y;\eta)=b(y)\exp\{\eta^TT(y)-a(\eta)\} \]
其中\(\eta\)是自然參數（natrue parameter）,\(T(y)\)是充分統計量，一般情況下\(T(y)=y\),當\(a,b,T\)確定時，上式就定義了一個以\(\eta\)為參數的函數族。下面討論將伯努力分布和高斯分布化為指數分布形式。
伯努力分布是對\(0,1\)問題進行建模的，設\(y \sim Bernoulli(\phi)\),即
\[p(y=1;\phi)=\phi \quad\quad\\ p(y=0;\phi)=1-\phi \]
我們可以得到
\[p(y;\phi)=\phi^y(1-\phi)^{1-y}=\exp\{y\;ln\phi+(1-y)ln(1-\phi)\}\\ =exp\{ y\ln(\frac{\phi}{1-\phi}) +ln(1-\phi)\} \]
其中
\[T(y)=y \quad\quad\quad\quad\quad\quad\quad\quad\quad\quad\\ \eta =ln(\frac{\phi}{1-\phi}) \Longrightarrow \phi=\frac{1}{1+e^{-\eta}}\\ a=-ln(1-\phi) =ln(1+e^{\eta})\quad\quad \]
這說明伯努力分布是指數分布族的一種，\(\phi\)的形式與邏輯回歸中的logitisc函數一樣，因為邏輯回歸對問題的潛質概率分布其實就是伯努力分布。
下面對高斯分布進行推導,設\(y \sim N(\mu,\sigma^2)\)，為了方便我們令\(\sigma=1\)
\[p(y;\mu)=\frac{1}{\sqrt{2\pi}}\exp\{-\frac{1}{2}(y-\mu)^2\} \\ =\frac{1}{\sqrt{2\pi}}\exp\{ \mu y-\frac{y^2+\mu^2}{2} \} \\ =\frac{1}{\sqrt{2\pi}}\exp\{-\frac{y^2}{2}\}\exp\{\mu y - \frac{\mu^2}{2}\} \]
於是得到
\[T(y)=\frac{1}{\sqrt{2\pi}}\exp\{-\frac{y^2}{2}\} \\ T(y)=y \quad\quad \quad\quad \quad\quad \\ \eta = \mu\quad\quad \quad\quad \quad\quad \quad \\ a(\mu)=\frac{\mu^2}{2} \quad\quad \quad\quad \\ \]
二、一般最小二乘模型
為了顯示一般最小二乘模型是廣義線性回歸模型的特例，我們設目標變量\(y\)是連續變量，我們用高斯分布來給\(y\)建模，設在\(x\)給定條件下\(y\)服從高斯分布\(N(\mu,\sigma^2)\).所以這里的指數族分布就是高斯分布（高斯分布屬於指數族分布）。正如前面看到的，將高斯分布寫成指數族分布得到\(\mu = \eta\),所以
\[h_\theta(x)=E[y|x;\theta] \quad\quad\quad(1) \\ =\mu \quad\quad \quad(2) \\ =\eta \quad\quad\quad(3) \\ =\theta^Tx \quad\quad (4) \]
注意，等式(1)是根據假設2得到的，由於\(y|x;\theta \sim N(\mu,\sigma^2)\),所以\(E[y|x;\theta]=\mu\)，於是得到等式（2），等式(3)是根據假設1得到的（前面推導表明，高斯分布表示為指數族時有\(\mu=\eta\)）,等式(4)是根據假設3得到的。
三、邏輯回歸模型
現在我們來考慮邏輯回歸模型，邏輯回歸模型考慮二分類問題\(y\in \{0,1\}\)，由於\(y\)是二值函數，我們考慮用伯努力分布來建模，設\(y|x\)服從參數為\(\phi\)的伯努力分布。將伯努力分布寫成指數族分布的形式得到
\[\phi = \frac{1}{1+e^{-\eta}} \]
注意到\(y|x;\theta \sim B(\phi)\),於是\(E[y|x;\theta]=\phi\).因此按照普通最小二乘法推導相似，我們有
\[h_\theta(x)=E[y|x;\theta]=\phi=\frac{1}{1+e^{-\eta}}=\frac{1}{1+e^{-\theta^Tx}} \]
因次，我們得到了假設函數\(h_\theta(x)=1/(1+e^{-\theta^Tx})\)，之前我們一直在疑惑我們為什么會想到邏輯回歸模型的邏輯函數的呢？這里就給出了答案：我們假設\(y|x\)下服從伯努力分布，它是廣義線性模型的定義得到的。
四、SoftMax回歸
我們現在來看廣義線性回歸的另一個例子。考慮一個多分類問題，其中\(y \in \{1,2,…,k\}\).例如我們不再將郵件分為垃圾郵件和非垃圾郵件2類，而是分為垃圾郵件、個人郵件、和工作郵件3個類別。我們用廣義線性回歸模型來為這個問題建模，首先我們需要將多項式分布表示成指數族分布的形式。
為表示\(k\)中可能結果的多項式分布，我們需要\(k\)個參數變量\(\phi_1,\phi_2,…,\phi_k\)來表示每個結果發生的概率。需要注意的是，這\(k\)個參數並不是完全獨立的，實際上只有\(k-1\)個自由變元（\(\phi_k=1-\sum^{k-1}_{i=1}\phi_i\)），所以我們用\(k-1\)個變元來表示多項式分布
\[\phi_i=p(y=i;\phi) \quad i=1,2,...,k-1 \\ p(y=k;\phi)=1-\sum_{i=1}^{k-1}\phi_i \quad \quad \quad \quad \quad \]
為了將多項式分布表示為指數族分布，我們定義\(T(y) \in R^{k-1}\).
\[T(1) = \left [ \begin{matrix} 1 \\ 0 \\ 0 \\ \vdots \\ 0 \end{matrix} \right ] \quad T(2) = \left [ \begin{matrix} 0 \\ 1 \\ 0 \\ \vdots \\ 0 \end{matrix} \right ] \quad T(3) = \left [ \begin{matrix} 0 \\ 0 \\ 1 \\ \vdots \\ 0 \end{matrix} \right ] \quad \cdots \quad T(k-1) = \left [ \begin{matrix} 0 \\ 0 \\ 0 \\ \vdots \\ 1 \end{matrix} \right ] \quad T(k) = \left [ \begin{matrix} 0 \\ 0 \\ 0 \\ \vdots \\ 0 \end{matrix} \right ] \quad \]
和之前不太一樣，這里\(T(y) \not = y\),並且\(T(y)\)是一個\(k-1\)維的向量，我們將會用\((T(y))_i\)表示\(T(y)\)的第\(i\)個分量值。
下面引入一個非常有用的指標函數。指標函數\(I\{\}\),當括號里面是一個真值表達式，為真返回1,為假返回0.例如
\(I\{2=3\}=0,I\{2<3\}=1\).我們可以將\(T(y)\)和\(y\)的關系寫成指標函數的形式\((T(y))_i=I\{y=i\}\)，更一般的有
\[E[((T(y))_i]=p(y=i)=\phi_i) \]
下面我們來推導將多項式分布表示為指數分布的形式
\[p(y;\phi)=\phi_1^{I\{y=1\}}\phi_2^{I\{y=2\}}\cdots\phi_k^{I\{y=k\}} \\ =\phi_1^{I\{y=1\}}\phi_2^{I\{y=2\}}\cdots\phi_k^{I\{1-\sum_{i=1}^{k-1}I\{y=i\}\}} \\ =\phi_1^{(T(y))_1}\phi_2^{(T(y))_2}\cdots\phi_k^{I\{1-\sum_{i=1}^{k-1}{(T(y))_i}\}} \\ =exp\{(T(y))_1ln(\phi_1)+(T(y))_2ln(\phi_2)+ \\ \cdots +(1-\sum_{i=1}^{k-1}{(T(y))_i})ln(\phi_k)\} \\ =exp\{(T(y))_1ln(\phi_1/\phi_k)+(T(y))_2ln(\phi_2/\phi_k)+ \\ \cdots +(T(y))_{k-1}ln(\phi_{k-1}/\phi_k)+ln(\phi_k)\} \\ =b(y)exp\{\eta^TT(y)-a(\eta)\} \]
其中
\[\eta = \left [ \begin{matrix} ln(\phi_1/\phi_k)\\ ln(\phi_2/\phi_k)\\ ln(\phi_3/\phi_k) \\ \vdots \\ ln(\phi_{k-1}/\phi_k) \end{matrix} \right ] \quad a(\eta)=-ln(\phi_k) \quad b(y)=1 \]
上面將多項式分布化為指數族分布的形式。"link function"表示為
\[\eta_i=ln(\phi_i/\phi_k) \quad i=1,2,...,k-1 \]
為了統一，我們定義\(\eta_k=ln(\phi_k/\phi_k)=0\),接下來求的響應函數為
\[\phi_i=e^{\eta_i}\phi_k\quad i=1,2,...,k \]
從而
\[\sum_{j=1}^k\phi_j=\phi_k\sum_{j=1}^ke^{\eta_j}=1 \]
於是我們得到\(\phi_k=1/\sum_{j=1}^ke^{\eta_j}\),同時得到
\[\phi_i=\frac{e^{\eta_i}}{\sum_{j=1}^ke^{\eta_j}} \quad i=1,2,...,k-1 \]
將\(\eta_i\)映射到\(\phi_i\)的函數稱為Softmax函數。再使用假設3，\(\eta_i\)是\(x\)的線性函數，即\(\eta_i=\theta^T_ix\)（其中\(i=1,2,..,k-1\)並且\(\theta_i \in R^{n+1}\)）。為了符號方便，我們定義\(\theta_k=0\),這樣\(\eta_k=\theta_k^Tx=0\),因此我們的模型為
\[p(y=i|x;\theta)=\phi_i =\frac{e^{\eta_i}}{\sum_{j=1}^ke^{\eta_j}} \\ =\frac{e^{\theta_i^Tx}}{\sum_{j=1}^ke^{\theta^T_jx}} \]
上面模型可以應用到\(k\)分類問題，稱為softmax回歸，它是邏輯回歸的推廣形式。我們的模型將會輸出
\[h_\theta(x)=E[T(y)|x;\theta]=E([I\{y=1\},I\{y=2\},...,I\{y=k-1\}]^T|x;\theta) \\ =[\phi_1,\phi_2,...,\phi_{k-1}]^T \quad \quad \quad \\ = \left[ \begin{matrix} \frac{e^{\theta_1^Tx}}{\sum_{j=1}^ke^{\theta^T_jx}} \\ \frac{e^{\theta_2^Tx}}{\sum_{j=1}^ke^{\theta^T_jx}} \\ \vdots \\ \frac{e^{\theta_{k-1}^Tx}}{\sum_{j=1}^ke^{\theta^T_jx}} \end{matrix} \right] \quad\quad\quad\quad\quad\quad \]
我們的模型將輸出概率\(p(y=i|x;\theta)\)（\(i=1,2,…,k-1\)），並且\(p(y=k|x;\theta)=1-\sum_{i=1}^{k-1}p(y=i|x;\theta)\)計算得到。
下面來討論下參數\(\theta\)的估計，和最小二乘與邏輯回歸類似，我們假設有\(m\)個樣本\(\{(x^{(i)},y^{(i)})|i=1,2,…,m\}\)用來學習參數\(\theta_i\)（\(i=1,2,..,k-1\)）.直接寫出對數似然函數為
\[l(\theta)=\sum_{i=1}^{m}ln\;p(y^{(i)}|x^{(i)};\theta) \quad\quad\;\;\quad\quad\;\;\quad\quad\quad\;\;\;\\ =\sum_{i=1}^m\ln\;\prod_{l=1}^k(\frac{e^{\theta_l^Tx^{(i)}}}{\sum_{j=1}^ke^{\theta^T_jx^{(i)}}} )^{I\{y^{(i)}=l\}}\quad\quad\; \\ =\sum_{i=1}^m\sum_{l=1}^k I\{y^{(i)}=l\}\;ln\;(\frac{e^{\theta_l^Tx^{(i)}}}{\sum_{j=1}^ke^{\theta^T_jx^{(i)}}} ) \\ =\sum_{i=1}^m\sum_{l=1}^k I\{y^{(i)}=l\}\;ln\;(h_\theta(x))_l \quad\quad\;\; \]
我們現在可以根據梯度上升或者牛頓法來求解參數\(\theta\)。
