广义线性模型

 指数分布族

前面学习了线性回归和logistic回归。我们知道对于$P(y|x;\theta)$

若y属于实数，满足高斯分布，得到基于最小二乘法的线性回归,若y取{0,1}，满足伯努利分布，得到Logistic回归。

这两个算法，其实都是广义线性模型的特例。

1) 伯努利分布

http://www.cnblogs.com/mikewolf2002/p/7667944.html

2) 高斯分布

http://www.cnblogs.com/mikewolf2002/p/7669977.html

指数分布族的定义：

     若一类概率分布可以写成如下形式，那么它就属于指数分布族：

$$P(y;\eta)=b(y)exp(\eta^{T}T(y)-a(\eta))$$

     $\eta$称作分布的自然参数，通常是一个实数。$T(y)$称作充分统计量，通常$T(y)=y$，实际上是一个概率分布的充分统计量（统计学知识）。对于给定的$a，b，T$三个函数，上式定义了一个以$\eta$为参数的概率分布集合，即改变$\eta$可以得到不同的概率分布。

https://en.wikipedia.org/w/index.php?title=Exponential_family&oldid=648989632

证明伯努利分布是指数分布族：

伯努利分布式对0,1问题进行建模，它的分布律如下：

$$Ber(\phi ),p(y=1; \varphi ) = \varphi; p(y=0; \varphi ) = 1- \varphi$$

伯努利分布是对0，1问题进行建模的，对于$Bernoulli（\varphi）,y\in \{0,1\}$有$p(y=1;\varphi)=\varphi;p(y=0;\varphi)=1−\varphi$，下面将其推导成指数分布族形式：

$$\begin{align*} L(\theta) &=P(y;\varphi)=\varphi^y(1-\varphi)^{1-y}=exp(log(\varphi^y(1-\varphi)^{1-y})) \\  &=exp(ylog(\varphi)+(1-y)log(1-\varphi)) \\ &= exp(ylog(\frac{\varphi}{1-\varphi})+log(1-\varphi))\end{align*}$$

    将其与指数族分布形式对比，可以看出：

$$b(y) = 1$$

$$T(y) = y$$

$$\eta=log\frac{\phi}{1-\phi} => \phi=\frac{1}{1+e^{-\eta}}$$

$$a(\eta)=-log(1-\phi)=log(1+e^{\eta})$$

    表明伯努利分布也是指数分布族的一种。从上述式子可以看到，$\eta$的形式与logistic函数（sigmoid）一致，这是因为 logistic模型对问题的前置概率估计其实就是伯努利分布。

高斯分布

    下面对高斯分布进行推导，推导公式如下（为了方便计算，我们将方差$\sigma^2$设置为1）：

$$N(\mu,1)=\frac{1}{\sqrt{2\pi}}exp^{-\frac{1}{2}(y-\mu)^2}=\frac{1}{\sqrt{2\pi}}exp^{-\frac{1}{2}y^2+y\mu-\frac{1}{2}\mu^2}=\frac{1}{\sqrt{2\pi}}exp^{-\frac{1}{2}y^2}exp^{y\mu-\frac{1}{2}\mu^2}$$

将上式与指数族分布形式比对，可知：

$$b(y)=\frac{1}{\sqrt{2\pi}}exp^{-\frac{1}{2}y^2}$$

$$T(y)=y$$

$$\eta=\mu$$

$$a(\eta)=\frac{\mu^2}{2}$$

两个典型的指数分布族，伯努利和高斯分布。其实大多数概率分布都可以表示成指数分布族形式，如下所示：

• 伯努利分布（Bernoulli）：对 0、1 问题进行建模；
• 多项式分布（Multinomial）：对 K 个离散结果的事件建模；
• 泊松分布（Poisson）：对计数过程进行建模，比如网站访问量的计数问题，放射性衰变的数目，商店顾客数量等问题；
• 伽马分布（gamma）与指数分布（exponential）：对有间隔的正数进行建模，比如公交车的到站时间问题；
• β 分布：对小数建模；
• Dirichlet 分布：对概率分布进建模；
• Wishart 分布：协方差矩阵的分布；

广义线性模型GLM

选定了一个指数分布族后，怎样来推导出一个GLM呢？

假设：

（1）$y|x; \theta~ExponentialFamily（\eta）$，即假设试图预测的变量$y$在给定$x$，以$\theta$作为参数的条件概率，属于以$\eta$作为自然参数的指数分布族，例：若要统计网站点击量$y$，用泊松分布建模。

（2） 给定$x$，目标是求出以$x$为条件的$T(y)$的期望$E[T(y)|x]$，即让学习算法输出$h(x) = E[T(y)|x]$

（3）$\eta=\theta^Tx$,即自然参数和输入特征$x$之间线性相关，关系由θ决定。仅当$\eta$是实数时才有意义,若$\eta$是一个向量，则$\eta_i=\theta_i^Tx$

推导伯努利分布的GLM：

$y|x; \theta~ExponentialFamily（\eta）$，伯努利分布属于指数分布族,对给定的$x,\theta$，学习算法进行一次预测的输出：

$$h_\theta(x)=E[y|x;\theta]=\phi=\frac{1}{1+e^{-\eta}}=\frac{1}{(1+e^{-\theta^Tx})}$$

得到logistic回归算法。

正则响应函数：$g(\eta)=E[y;\eta]$,将自然参数$\eta$和$y$的期望联系起来

正则关联函数:$g^{-1}$

推导多项式分布的GLM：

多项式分布是在$k$个可能取值上的分布，即$y\in {1,..,k}$,如将收到的邮件分成$k$类，诊断某病人为$k$种病中的一种等问题。

（1）将多项式分布写成指数分布族的形式：

设多项式分布的参数：$\phi_1,\phi_2,…,\phi_k$ ，且$\sum\limits_{i=1}^k\phi_i=1$ ，$\phi_i$表示第$i$个取值的概率分布，最后一个参数可以由前$k-1$个推导出，所以只将前$k-1$个$\phi_1,\phi_2,…,\phi_{k-1}$视为参数。

多项式分布是少数几个$T(y)!=y$的分布，$T(1) \sim T(k)$都定义成一个$k-1$维的列向量，表示为：

$$T(1)=\begin{bmatrix} 1\\ 0\\ 0\\ \vdots\\ 0\end{bmatrix}T(2)=\begin{bmatrix} 0\\ 1\\ 0\\ \vdots\\ 0\end{bmatrix}T(3)=\begin{bmatrix} 0\\ 0\\ 1\\ \vdots\\ 0\end{bmatrix}T(k-1)=\begin{bmatrix} 0\\ 0\\ 0\\ \vdots\\ 1\end{bmatrix}T(k)=\begin{bmatrix} 0\\ 0\\ 0\\ \vdots\\ 0\end{bmatrix}$$

这样定义$T(y)$是为了将多项式分布写成指数分布族形式。

*定义符号：指示函数，$1\{.\}$

1{True} = 1, 1{False} = 0，即大括号内命题为真，值为1,否则为0。

例：1{2=3} = 0, 1{1+1=2} = 1

用$T(y)_i$表示$T(y)$的第$i$个元素，则$T(y)_i=1\{y=i\}$

根据参数$\phi$的意义（$\phi_i$表示第$i$个取值的概率分布），可推出：

$$\begin{align*} p(y;\phi) &=\phi_1^{1\{y=1\}}\phi_2^{1\{y=2\}}...\phi_k^{1\{y=k\}} \\  &=\phi_1^{1\{y=1\}}\phi_2^{1\{y=2\}}...\phi_k^{1-\sum\limits_{i=1}^{k-1}1\{y=i\}} \\ &=\phi_1^{(T(y))_1}\phi_2^{(T(y))_2}...\phi_k^{1-\sum\limits_{i=1}^{k-1}(T(y))_i} \\ &=>exp(T(y))_1log(\phi_1)+ exp(T(y))_2log(\phi_2)+...+(1-\sum\limits_{i=1}^{k-1}(T(y))_i)log(\phi_k) \\ &=exp(T(y))_1log(\phi_1/\phi_k)+ exp(T(y))_2log(\phi_2/\phi_k)+...+exp(T(y))_{k-1}log(\phi_{k-1}/\phi_k)+log(\phi_k)) \\ &=b(y)exp(\eta^{T}T(y)-a(\eta))\end{align*}$$

可得：

$$\eta=\begin{bmatrix} log(\phi_1/\phi_k)\\ log(\phi_2/\phi_k)\\ \vdots\\ log(\phi_{k-1}/\phi_k) \end{bmatrix}$$

$$a(\eta)=-log(\phi_k)$$

$$b(y)=1$$

证明多项式分布是指数分布族。

再用$\eta$表示$\phi$：

$$\eta_i=log\frac{\phi_i}{\phi_k} \\ e^{\eta_i}=\frac{\phi_i}{\phi_k} \\ \phi_ke^{\eta_i}=\phi_i \\ \phi_k\sum\limits_{i=1}^ke^{\eta_i}=\sum\limits_{i=1}^k\phi_i=1 \\ \phi_i=\frac{e^{\eta_i}}{\sum\limits_{j=1}^ke^{\eta_j}}$$

（2）根据上述假设（3）中自然参数和输入$x$的线性关系，可求得：

$$\begin{align*} p(y=i|x;\theta) &=\phi_i \\    &=\frac{e^{\eta_i}}{\sum\limits_{j=1}^ke^{\eta_j}} \\ &=\frac{e^{\theta_i^Tx}}{\sum\limits_{j=1}^ke^{\theta_j^Tx}} \\ \end{align*}$$

（3）根据上述假设（2）中的输出$h(x) = E[T(y)|x]$，可求得：

$$\begin{align*} h_\theta(x) &=E[T(y)|x;\theta] \\    &=\begin{bmatrix} 1\{y=1\}|x;\theta\\ 1\{y=2\}|x;\theta\\ \vdots\\ 1\{y=k-1\}|x;\theta \end{bmatrix} \\ &= \begin{bmatrix} \phi_1\\ \phi_2\\ \vdots\\ \phi_{k-1} \end{bmatrix}\\ &=\begin{bmatrix} \frac{e^{\theta_1^Tx}}{\sum\limits_{j=1}^ke^{\theta_j^Tx}} \\ \frac{e^{\theta_2^Tx}}{\sum\limits_{j=1}^ke^{\theta_j^Tx}} \\ \vdots\\ \frac{e^{\theta_{k-1}^Tx}}{\sum\limits_{j=1}^ke^{\theta_j^Tx}} \end{bmatrix} \end{align*}$$

称这种回归算法为softmax回归，是logistic回归的推广。

Softmax回归的训练方法和logistic相似，通过极大似然估计找到参数$\theta$,其对数似然性为：

$$\begin{align*} l(\theta) &=\sum\limits_{i=1}^{m}logp(y^{(i)}|x^{(i)};\theta); \\    &=\sum\limits_{i=1}^{m}log\prod\limits_{l=1}^k\bigg(\frac{e^{\theta_l^Tx^{(i)}}}{\sum\limits_{j=1}^ke^{\theta_j^Tx^{(i)}}}\bigg)^{1\{y^{(i)}=l\}} = \sum\limits_{i=1}^m\bigg(\sum\limits_{l=1}^k\bigg( 1\{y^{(i)}=l\}\log\bigg( \frac{e^{\theta_l^Tx^{(i)}}}{\sum\limits_{j=1}^ke^{\theta_j^Tx^{(i)}} }\bigg)\bigg)\bigg)\end{align*}$$

其中$\theta$是一个$k\times(n+1)$的矩阵，代表这$k$个类的所有训练参数，每个类的参数是一个$n+1$维的向量。

跟Logistic回归一样，softmax也可以用梯度下降法或者牛顿迭代法求解，但在softmax回归中直接用上述对数似然函数求导是不能更新参数的，因为它存在冗余的参数，通常用牛顿方法中的Hessian矩阵也不可逆，是一个非凸函数，那么可以通过添加一个权重衰减项来修改代价函数，使得代价函数是凸函数，并且得到的Hessian矩阵可逆

再通过梯度上升或牛顿方法找对数似然性的最大值，即可求出参数$\theta$。

详细的公式推导见：http://ufldl.stanford.edu/wiki/index.php/Softmax%E5%9B%9E%E5%BD%92