[CS229学习笔记] 4.牛顿法，指数族分布与广义线性模型，softmax

本文对应的是吴恩达老师的CS229机器学习的第四课。这节课先介绍了牛顿法，然后给出了指数族的定义，并从指数族出发，介绍了广义线性模型并以此解释最小二乘、逻辑回归、softmax等模型的来源。

---

牛顿法(Newton’s Method)

牛顿法是另一种求解曲线零点的方法，其具体做法如果所示：从某一起始点 ${\theta }^{(0)}$ 开始，找到其对应的函数值 $f({\theta }^{(0)})$，然后作切线，以切线与横坐标轴的交点作为更新值，多次迭代直至收敛。

从上图中我们可以看出，参数 ${\theta }^{(0)}$ 更新的值可以由当前函数值与斜率共同得出，即：

$$\Delta =\frac{f({\theta }^{(0)})}{f^{\prime }({\theta }^{(0)})}$$

对应的 ${\theta }^{(0)}$ 更新公式就是：

$${\theta }^{(t+1)}={\theta }^{(t)}-\frac{f({\theta }^{(t)})}{f^{\prime }({\theta }^{(t)})}$$

那么对于我们上节课中提到的似然函数 $l(\theta )$，我们想找到其最大值，即其导数零点 $l^{\prime }(\theta )=0$，只需将牛顿法代入函数 $l^{\prime }(\theta )$ 中进行迭代求解即可。更新公式为：

$${\theta }^{(t+1)}={\theta }^{(t)}-\frac{l^{\prime }({\theta }^{(t)})}{l^{\prime \prime }({\theta }^{(t)})}$$

其对应的矩阵形式为（证明从略）：

$${\theta }^{(t+1)}={\theta }^{(t)}-H^{-1}{\nabla }_{\theta }l$$

其中 $H=\frac{{\partial }^{2}l}{\partial {\theta }_i\partial {\theta }_j}$ 是海森矩阵(Hessian matrix)。

牛顿法的优点在于其收敛速度很快，为二次收敛速度；缺点是海森矩阵的求逆操作非常慢（求逆操作时间复杂度为$o(n^3)$）。

---

指数族分布(exponential family of distributions)与广义线性模型(generalized linear model)

接下来，课程给出了一个非常重要的概念：指数族分布。指数族分布是一系列非常常见且常用的连续及离散概率分布的统称。其包含了高斯分布（正态分布）、伯努利分布、二项分布、多项分布、指数分布、泊松分布、伽马分布、贝塔分布、狄利克雷分布、威希特分布等等。其形式如下：

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

其中 $\eta$ 被称为分布的自然参数；$T(y)$ 被称为充分统计量。

同时，我们先给出广义线性模型的三个假定条件

1. 首先，我们假定 $y$ 的分布服从某一指数族分布，即 $y|x;\theta \backsim \text{ExponentialFamily}(\eta )$
2. 给定数据点 $x$，我们的目标是预测 $T(y)$ 在给定 $x$ 下的数学期望，即 $h_{\theta }(x)=\mathbb{E}[y|x]$
3. 对于广义线性模型，我们假定自然参数 $\eta$ 与输入 $x$ 线性相关，即 $\eta ={\theta }^{T}x$

接下来，我们就可以利用指数族分布及广义线性模型的定义来作一些公式推导。

首先我们回到上节课提到的逻辑函数。当我们看到逻辑函数时，不免要问：为什么逻辑函数是这样的一个指数形式而不是其他的形式？接下来，我们就来证明两个常见的分布：伯努利分布、高斯分布属于指数分布族，并从指数族分布出发，分别推出两个之前课程中提到过的广义线性模型：逻辑函数、最小均方差。

1. 伯努利分布
   $$\text{Bernoulli}(\varphi )=P(y=1;\varphi )=\varphi$$
   $$\text{Bernoulli}(\varphi )=P(y=0;\varphi )=1-\varphi$$
   结合两式可得：$$\begin{array}{rl}P(y;\varphi )& ={\varphi }^y(1-\varphi {)}^{(1-y)}\\ & =\exp \left(\log \left({\varphi }^y(1-\varphi {)}^{(1-y)}\right)\right)\\ & =\exp \left(y\log \varphi +(1-y)\log (1-\varphi )\right)\\ & =\exp \left(\log \frac{\varphi }{1-\varphi }y+\log (1-\varphi )\right)\end{array}$$
   这就是指数族分布的形式，其中 $\eta =\log \frac{\varphi }{1-\varphi },a(\eta )=-\log (1-\varphi ),b(y)=1,T(y)=y$，注意此时 $\eta$ 为标量，因此其转置等于本身。我们考虑我们的预测函数 $h_{\theta }(x)$，在给定某概率分布的条件以及权重 $\theta$ 的情况下，我们的预测值就应该是当前参数给出的数学期望，并利用得到的 $\eta =\log \frac{\varphi }{1-\varphi }$，也就是 $\varphi =\frac{1}{1+e^{-\eta }}$，我们可以得到：
   $$\begin{array}{rl}{h}_{\theta }(\mathbf{x})=\mathbb{E}[y|\mathbf{x};\theta ]& =1\cdot P(y=1|\mathbf{x};\theta )+0\cdot P(y=0|\mathbf{x};\theta )\\ & =P(y=1|\mathbf{x};\theta )\\ & =\varphi \\ & =\frac{1}{1+e^{-\eta }}\end{array}$$
   于是我们就推导出了逻辑回归公式：$h_{\theta }(\mathbf{x})=\frac{1}{1+e^{-{\theta }^T\mathbf{x}}}$。

2. 高斯分布
   类似地，我们可以写出高斯分布的公式（简单起见，我们设其方差为1）：
   $$\begin{array}{rl}N(\mu ,{\sigma }^2=1)& =\frac{1}{\sqrt{2\pi }}\exp \left(-\frac{1}{2}(y-\mu {)}^2\right)\\ & =\underset{b(y)}{\underbrace{\frac{1}{\sqrt{2\pi }}\exp \left(-\frac{1}{2}{y}^2\right)}}\exp \left(\underset{\eta T(y)}{\underbrace{\mu y}}-\underset{a(\eta )}{\underbrace{\frac{1}{2}{\mu }^2}}\right)\end{array}$$
   同样的，我们可以用高斯分布来求预测值
   $$\begin{array}{rl}{h}_{\theta }(\mathbf{x})=\mathbb{E}[y|\mathbf{x};\theta ]& =\mu \\ & =\eta \\ & ={\theta }^T\mathbf{x}\end{array}$$
   即线性回归公式。

根据上面的推导我们可以发现：线性回归与逻辑回归从概率角度的出发点其实分别是高斯分布与二项分布。高斯分布与二项分布符合我们对连续、离散参数的规律的估计，因此逻辑回归公式的形式并不是完全的经验公式，而是有理可依的。

另外，课程还给出两个数学上的定义：

1. 正则响应函数(canonical response function)：$g(\eta )=\mathbb{E}[y;\eta ]=\frac{1}{1+e^{-\eta }}$

2. 正则关联函数(canonical link function)：$g^{-1}$，即正则响应函数的逆函数

显然，伯努利分布的正则响应函数是逻辑函数。

---

softmax的原理及推导

接下来，利用一个稍复杂的指数族分布，我们来推导一个经常在神经网络中被用来做多分类任务的广义线性模型：softmax。

首先，我们以多分类问题出发，假设我们要对多个类进行区分，那么对于多类问题，多项分布是一个很好的概率分布模型。同时，对于充分统计量，我们将其写为各个方向的基向量的形式，即

$$T(1)=\left[\begin{array}{c}1\\ 0\\ ⋮\\ 0\end{array}\right],T(2)=\left[\begin{array}{c}0\\ 1\\ ⋮\\ 0\end{array}\right],\cdots ,T(k-1)=\left[\begin{array}{c}0\\ 0\\ ⋮\\ 1\end{array}\right],T(k)=\left[\begin{array}{c}0\\ 0\\ ⋮\\ 0\end{array}\right]$$

这里，我们介绍一个新的函数：指示函数(indicator function)。我们定义指数函数的形式为 1 { ⋅ } \mathbb{1}\{\cdot\} 1{⋅}，用来表示括号内的语句是否为真，例如 1 { 2 = 3 } = 0 , 1 { 2 + 1 = 3 } = 1 \mathbb{1}\{2=3\}=0, \mathbb{1}\{2+1=3\}=1 1{2=3}=0,1{2+1=3}=1。于是，多项分布的充分统计量可以写为：

T ( y ) i = 1 { y = i } for i in  { 1 , 2 , ⋯   , k − 1 } T(y)_i=\mathbb{1}\{y=i\} \qquad\text{for i in } \{1,2,\cdots,k-1\} T(y)i​=1{y=i}for i in {1,2,⋯,k−1}

这里简单说明一下，$T(y)$ 是一个向量，且除了第 $y$ 个值为 $1$ 外其他值都为 $0$，因此第 $i$ 个分量可以用指数函数 1 { y = i } \mathbb{1}\{y=i\} 1{y=i} 来描述。假设预测结果可以被分为 $k$ 类，每一类的概率 ${\varphi }_i$，即$P(y=i)={\varphi }_i$。由于所有类的概率之和总为 $1$，因此第 $k$ 类总可以用前 $1-k$ 类的概率来表示，即 ${\varphi }_k=1-{\varphi }_1-\cdots -{\varphi }_{k-1}$。把每一类的概率写成同一个函数，则可得多项分布概率公式为：

P ( y ; ϕ ) = ϕ 1 1 { y = 1 } ϕ 2 1 { y = 2 } ⋯ ϕ k 1 { y = k } = ϕ 1 T ( y ) 1 ϕ 2 T ( y ) 2 ⋯ ϕ k T ( y ) k = exp ⁡ ( log ⁡ ( ϕ 1 T ( y ) 1 ) ) exp ⁡ ( log ⁡ ( ϕ 2 T ( y ) 2 ) ) ⋯ exp ⁡ ( log ⁡ ( ϕ k T ( y ) k ) ) = exp ⁡ ( T ( y ) 1 log ⁡ ( ϕ 1 ) ) exp ⁡ ( T ( y ) 2 log ⁡ ( ϕ 2 ) ) ⋯ exp ⁡ ( T ( y ) k log ⁡ ( ϕ k ) ) = exp ⁡ ( T ( y ) 1 log ⁡ ( ϕ 1 ) + T ( y ) 2 log ⁡ ( ϕ 2 ) + ⋯ + T ( y ) k − 1 log ⁡ ( ϕ k − 1 ) + T ( y ) k log ⁡ ( ϕ k ) ) = exp ⁡ ( T ( y ) 1 log ⁡ ( ϕ 1 ) + T ( y ) 2 log ⁡ ( ϕ 2 ) + ⋯ + T ( y ) k − 1 log ⁡ ( ϕ k − 1 ) + ( 1 − T ( y ) 1 − T ( y ) 2 − ⋯ − T ( y ) k − 1 ) log ⁡ ( ϕ k ) ) = exp ⁡ ( T ( y ) 1 log ⁡ ( ϕ 1 ϕ k ) + T ( y ) 2 log ⁡ ( ϕ 2 ϕ k ) + ⋯ + T ( y ) k − 1 log ⁡ ( ϕ k − 1 ϕ k ) ⏟ η T T ( y ) + 1 ⋅ log ⁡ ( ϕ k ) ⏟ − a ( η ) ) \begin{aligned}P(y;\phi)&=\phi_1^{\mathbb{1}\{y=1\}}\phi_2^{\mathbb{1}\{y=2\}}\cdots\phi_k^{\mathbb{1}\{y=k\}}\\ &=\phi_1^{T(y)_1}\phi_2^{T(y)_2}\cdots\phi_k^{T(y)_k}\\ &=\exp\left(\log\left(\phi_1^{T(y)_1}\right)\right)\exp\left(\log\left(\phi_2^{T(y)_2}\right)\right)\cdots\exp\left(\log\left(\phi_k^{T(y)_k}\right)\right)\\ &=\exp\left(T(y)_1\log\left(\phi_1\right)\right)\exp\left(T(y)_2\log\left(\phi_2\right)\right)\cdots\exp\left(T(y)_k\log\left(\phi_k\right)\right)\\ &=\exp\left(T(y)_1\log\left(\phi_1\right)+T(y)_2\log\left(\phi_2\right)+\cdots+T(y)_{k-1}\log\left(\phi_{k-1}\right)+T(y)_k\log\left(\phi_k\right)\right)\\ &=\exp\left(T(y)_1\log\left(\phi_1\right)+T(y)_2\log\left(\phi_2\right)+\cdots+T(y)_{k-1}\log\left(\phi_{k-1}\right)+\left(1-T(y)_1-T(y)_2-\cdots-T(y)_{k-1}\right)\log\left(\phi_k\right)\right)\\ &=\exp\left(\underbrace{T(y)_1\log\left(\frac{\phi_1}{\phi_k}\right)+T(y)_2\log\left(\frac{\phi_2}{\phi_k}\right)+\cdots+T(y)_{k-1}\log\left(\frac{\phi_{k-1}}{\phi_k}\right)}_{\eta^TT(y)}+\underbrace{1\cdot\log\left(\phi_k\right)}_{-a(\eta)}\right)\end{aligned} P(y;ϕ)​=ϕ11{y=1}​ϕ21{y=2}​⋯ϕk1{y=k}​=ϕ1T(y)1​​ϕ2T(y)2​​⋯ϕkT(y)k​​=exp(log(ϕ1T(y)1​​))exp(log(ϕ2T(y)2​​))⋯exp(log(ϕkT(y)k​​))=exp(T(y)1​log(ϕ1​))exp(T(y)2​log(ϕ2​))⋯exp(T(y)k​log(ϕk​))=exp(T(y)1​log(ϕ1​)+T(y)2​log(ϕ2​)+⋯+T(y)k−1​log(ϕk−1​)+T(y)k​log(ϕk​))=exp(T(y)1​log(ϕ1​)+T(y)2​log(ϕ2​)+⋯+T(y)k−1​log(ϕk−1​)+(1−T(y)1​−T(y)2​−⋯−T(y)k−1​)log(ϕk​))=exp⎝⎜⎜⎜⎛​ηTT(y) T(y)1​log(ϕk​ϕ1​​)+T(y)2​log(ϕk​ϕ2​​)+⋯+T(y)k−1​log(ϕk​ϕk−1​​)​​+−a(η) 1⋅log(ϕk​)​​⎠⎟⎟⎟⎞​​

于是，我们便完成了将多项分布写为指数族分布的形式，同时，注意到我们此时的 $\eta$ 可以写为矩阵形式，即：

$$\eta =\left[\begin{array}{c}\log \frac{{\varphi }_1}{{\varphi }_k}\\ \log \frac{{\varphi }_2}{{\varphi }_k}\\ ⋮\\ \log \frac{{\varphi }_{k-1}}{{\varphi }_k}\end{array}\right]$$

每个方向上有 ${\eta }_i=\log \frac{{\varphi }_i}{{\varphi }_k}$。对其作指数变换，我们可以得到 ${\varphi }_i=e^{{\eta }_i}{\varphi }_k$。由于 ${\varphi }_k$ 可以用 ${\varphi }_1,\cdots ,{\varphi }_{1-k}$ 来表示，因此我们想到用求和的方式求出 ${\varphi }_k$。于是我们对各个 ${\varphi }_i$ 求和，可得 ${\sum }_i^k{\varphi }_i=\left({\sum }_i^k{e}^{{\eta }_i}\right){\varphi }_k$，根据概率的定义，${\sum }_i^k{\varphi }_i=1$，因此 ${\varphi }_k=\frac{1}{{\sum }_i^k{e}^{{\eta }_i}}$。于是我们可以得到多项分布正则响应函数为：

$${\varphi }_i=\frac{e^{{\eta }_i}}{{\sum }_{j=1}^k{e}^{{\eta }_j}}$$

同样的，我们假定 $\eta$ 与 $\mathbf{x}$ 线性相关，则有：

$${\varphi }_i=\frac{e^{{\theta }_i^T\mathbf{x}}}{{\sum }_{j=1}^k{e}^{{\theta }_j^T\mathbf{x}}}$$

这就是我们常用的softmax函数。下面我们同样写出softmax对应的预测函数：

h θ ( x ) = E [ T ( y ) ∣ x ; θ ] = E [ 1 { y = 1 } 1 { y = 2 } ⋮ 1 { y = k − 1 } ∣ x ; θ ] = [ ϕ 1 ϕ 2 ⋮ ϕ k − 1 ] = [ exp ⁡ ( θ 1 T x ) ∑ j = 1 k exp ⁡ ( θ j T x ) exp ⁡ ( θ 2 T x ) ∑ j = 1 k exp ⁡ ( θ j T x ) ⋮ exp ⁡ ( θ k − 1 T x ) ∑ j = 1 k exp ⁡ ( θ j T x ) ] \begin{aligned}h_\theta(\bold{x})&=\mathbb{E}\left[T(y)|\bold{x};\theta\right]=\mathbb{E}\left[\left.\begin{matrix}\mathbb{1}\{y=1\}\\\mathbb{1}\{y=2\}\\\vdots\\\mathbb{1}\{y=k-1\}\end{matrix}\right|\bold{x};\theta\right]=\begin{bmatrix}\phi_1\\\phi_2\\\vdots\\\phi_{k-1}\end{bmatrix}\\ &=\begin{bmatrix}\frac{\exp{\left(\theta_{1}^T\bold{x}\right)}}{\sum_{j=1}^k\exp{\left(\theta_j^T\bold{x}\right)}}\\\frac{\exp{\left(\theta_{2}^T\bold{x}\right)}}{\sum_{j=1}^k\exp{\left(\theta_j^T\bold{x}\right)}}\\\vdots\\\frac{\exp{\left(\theta_{k-1}^T\bold{x}\right)}}{\sum_{j=1}^k\exp{\left(\theta_j^T\bold{x}\right)}}\end{bmatrix}\end{aligned} hθ​(x)​=E[T(y)∣x;θ]=E⎣⎢⎢⎢⎡​1{y=1}1{y=2}⋮1{y=k−1}​∣∣∣∣∣∣∣∣∣​x;θ⎦⎥⎥⎥⎤​=⎣⎢⎢⎢⎡​ϕ1​ϕ2​⋮ϕk−1​​⎦⎥⎥⎥⎤​=⎣⎢⎢⎢⎢⎢⎢⎢⎡​∑j=1k​exp(θjT​x)exp(θ1T​x)​∑j=1k​exp(θjT​x)exp(θ2T​x)​⋮∑j=1k​exp(θjT​x)exp(θk−1T​x)​​⎦⎥⎥⎥⎥⎥⎥⎥⎤​​

然后其对数似然函数可以写为：

l ( θ ) = ∑ i = 1 m log ⁡ P ( y ( i ) ∣ x ( i ) ; θ ) = ∑ i = 1 m log ⁡ ∏ l = 1 k ( e θ l T x ( i ) ∑ j = 1 k e θ j T x ( i ) ) 1 { y ( i ) = 1 } \begin{aligned}l(\theta)&=\sum_{i=1}^m\log P(y^{(i)}|x^{(i)};\theta)\\ &=\sum_{i=1}^m\log\prod_{l=1}^k\left(\frac{e^{\theta_l^T\bold{x}^{(i)}}}{\sum_{j=1}^ke^{\theta_j^T\bold{x}^{(i)}}}\right)^{\mathbb{1}\{y^{(i)}=1\}}\end{aligned} l(θ)​=i=1∑m​logP(y(i)∣x(i);θ)=i=1∑m​logl=1∏k​(∑j=1k​eθjT​x(i)eθlT​x(i)​)1{y(i)=1}​

同样地，利用最大似然函数求上式的导数，或者利用梯度下降法，即可得到softmax对应的最优 $\theta$ 值。