CS229机器学习笔记 监督学习(Supervised Learning)（1）

目录

• 监督学习模型（Supervised Learning Model）
  • 线性回归（Liner Regression）
    • 从概率理论推导成本函数
    • 局部加权线性回归（Locally Weighted Linear Regression）
  • 逻辑回归（Logistic Regression）
  • 感知器学习算法（The Perceptron Learning Algorithm）
  • Softmax 回归（Softmax Regression）
• 寻找最优参数 $\theta$
  • 梯度下降/上升法（Gradient Descent/Gradient Ascent）
    • 批量梯度下降法（Batch Gradient Descent）
    • 随机梯度下降法（Stochastic Gradient Descent）
    • 小批次梯度下降（Mini-Batch Gradient Descent）
  • 正规方程（The Normal Equations）
  • 牛顿法（Newton's Method）
• 广义线性模型（General Liner Model）
  • 指数族（The Exponential Family）
    • 高斯分布（Gaussian）
    • 伯努利分布（Bernoulli）
    • 其它分布
  • 构建广义线性模型（Constructing GLMs）
    • 构建普通最小二乘模型
    • 构建逻辑回归模型
    • 构建Softmax回归模型
• 总结
• 参考

监督学习本质上就是，输入训练集（数据的集合），输出一个合适的函数 h（假设,hypothesis）。而这个函数 h，通过输入数据特征可以得到一个接近真实合理的评估输出。

举个例子，以下是训练集的数据是俄勒冈州波特兰市的 40 套房屋的面积、卧室数量和价格，我们希望通过输入面积来预测价格：

居住面积（平方英尺）	价格（千美元）
21042104	400400
16001600	330330
24002400	369369
14161416	232232
...	...

数据$(x,y)$：

• 数据的特征（函数输入）：$x = [ x_0,x_1,x_2,...,x_n ]^T$，共有n个特征$x_1,x_2,...,x_n$（特殊地，$x_0=1$，不属于特征），例子中的居住面积就是数据的一个特征。
• 数据的标签（函数真实输出）：$y$，例子中的价格就是数据的标签。

训练集（数据的集合）：

• 设共有m个数据，包含$(x^{(1)},y^{(1)})$...$(x^{(m)},y^{(m)})$，例子中的 m = 40

监督学习模型（Supervised Learning Model）

线性回归（Liner Regression）

回归问题：输入训练集，输出一个合适的 $h$ 函数（假设，hypothesis），并且这个函数的输出是实数。

为了解决回归问题，我们必须得定义好用什么样的 $h$ 函数。

线性回归一类模型则是假设数据集中特征与结果存在着线性关系：

$y^{(i)}=\theta^T x^{(i)}+ \epsilon ^{(i)}$

上式中$\epsilon ^{(i)}$是误差项，用于存放由于建模所忽略的变量导致的效果（比如可以是卖主的心情、货币的贬值等杂七杂八的因素）

我们将通过训练集的数据拟合成的线性函数作为假设函数 $h$ ：

$\begin{aligned}h_\theta(x) &= \theta_0 + \theta_1 \times x_1 + \theta_2 \times x_2 + ... + \theta_n \times x_n = \sum_{i=0}^n \theta _{i} x_{i} \\ &= \theta^Tx \end{aligned}$

其中，n 是数据特征的数量，$\theta$ 是 n+1 维向量，代表着 $h$ 函数的参数，$x$ 是n+1维向量，代表着 $h$ 函数的输入。

在房子的例子中，图中直线即为通过一定方法拟合成的线性函数：

定义好基本的 $h$ 函数（线性函数）形状后，我们需要得到合适的 $\theta$ 参数，尽可能使得 $h_\theta(x)$ 逼近 $y$ 。

为此，可以定义一个 $J(\theta)$ 最小二乘成本函数（cost function）：

$J(\theta) = \frac{1}{2}\sum_{i=1}^{m}(h_{\theta}(x^{(i)})-y^{(i)})^2$

当成本函数 $J(\theta)$ 越小时，自然而然会认为 $h_\theta(x)$ 越贴近 $y$。

接下来就需要找到合适的参数 $\theta$ 来尽可能让成本函数 $J(\theta)$ 最小化。当然寻找参数 $\theta$ 有很多，在下文的第二章将会介绍几个方式，例如梯度下降法。

需要用到最小二乘成本函数 $J(\theta)$ 的这个回归模型也叫做普通最小二乘回归模型(ordinary least squares regression model)。

从概率理论推导成本函数

为什么这个最小二乘成本函数是合理的，我们可以给这个假设函数做一系列的统计学假设，然后就能推出这个最小二乘成本函数。

我们可以进一步假设误差项 $\epsilon ^{(i)}$ 是独立同分布的 (IID) ，服从高斯分布（Gaussian distribution），其平均值为0，方差为 $\sigma ^2$。

那么，$\epsilon ^{(i)}$的密度函数就是：

$p(\epsilon ^{(i)} )= \frac 1{\sqrt{2\pi\sigma}} exp (- \frac {(\epsilon ^{(i)} )^2}{2\sigma^2})$

这也意味着存在下面的等量关系：

$p(y ^{(i)} |x^{(i)}; \theta)= \frac 1{\sqrt{2\pi\sigma}} exp (- \frac {(y^{(i)} -\theta^T x ^{(i)} )^2}{2\sigma^2})$

这里的记号 $p(y ^{(i)} |x^{(i)}; \theta)$ 表示的是这是一个对于给定 $x^{(i)}$ 的 $y^{(i)}$ 的分布，用 $\theta$ 进行了参数化。注意不能用 $\theta (p(y ^{(i)} |x^{(i)},\theta))$ 来当做条件，因为 $\theta$ 并不是一个随机变量。这个$y^{(i)}$的分布还可以写成 $y^{(i)} | x^{(i)}; \theta ∼ N (\theta ^T x^{(i)}, \sigma^2)$。

上面仅仅是单个数据样本 $x^{(i)}$ 和 $y^{(i)}$的假设，倘若所有数据样本都给出同样的独立性假设，然后给定一个 $X$ 为设计矩阵（design matrix），包含了全部 $x^{(i)}$，那么 $\vec{y}$ 的分布要当做 $\theta$ 的函数的时候，就称它为 似然函数（likelihood function) ：

$\begin{aligned} L(\theta) &=L(\theta;X,\vec{y})=p(\vec{y}|X;\theta)\ \end{aligned}$

把上面的等式改写成下面的形式：

$\begin{aligned} L(\theta) &=\prod ^m_{i=1}p(\vec{y}|X;\theta)\\ &=\prod^m _{i=1} \frac 1{\sqrt {2\pi \sigma }} exp(- \frac {(y^{(i)}-\theta^T x^{(i)})^2}{2\sigma^2})\ \end{aligned}$

现在有了 $x^{(i)}$ 和 $y^{(i)}$ 之间关系的概率模型后，就可以用 最大似然法（maximum likelihood）来选择对参数 $\theta$ 的最佳猜测：即需要选择一个 $\theta$ 能够让函数$L(\theta)$取到最大值。

为了让求导运算简化，我们通过找对数函数 $l(\theta)$ 的最大值来间接找到 $L(\theta)$ 的最大值：

$\begin{aligned} l(\theta) &= \log L(\theta)\ \\ &=\log \prod ^m _{i=1} \frac 1{\sqrt {2\pi \sigma }} exp(- \frac {(y^{(i)}-\theta^T x^{(i)})^2}{2\sigma^2})\ \\ &= \sum^m_{i=1} \log \frac 1{\sqrt {2\pi \sigma }} exp(- \frac {(y^{(i)}-\theta^T x^{(i)})^2}{2\sigma^2}) \\ &= m \log \frac 1{\sqrt {2\pi \sigma}}- \frac 1{\sigma^2}\cdot \frac 12 \sum^m_{i=1} (y^{(i)}-\theta^Tx^{(i)})^2\ \end{aligned}$

但是仔细观察就会发现，对 $l(\theta)$ 取得最大值也就意味着下面这个子式取到最小值：

$\frac 12 \sum^m _{i=1} (y^{(i)}-\theta^Tx^{(i)})^2$

这个子式正是刚刚定义的 $J(\theta)$ 成本函数。

在对数据进行概率假设的基础上，最小二乘回归得到的 $\theta$ 和最大似然法估计的 $\theta$ 是一致的。对于验证最小二乘法是否为一个良好并且合理的过程来说，这些概率假设并不是必须的，此外可能（也确实）有其他的自然假设能够用来评判最小二乘方法。

局部加权线性回归（Locally Weighted Linear Regression）

左下角的图显示了使用一次多项式 $y = \theta_0 + \theta_1x$ 来对数据集进行拟合，但这个数据集的趋势显然不是一条严格的直线，因此用一次多项式进行的拟合效果不太好。若增加一个二次项，用二次多项式 $y = \theta_0 + \theta_1x +\theta_2x^2$来拟合，那么拟合的效果更加好。

通常来说，我们对特征补充得越多，效果就越好。不过增加太多特征也可能是不合理的：最右边的图就是使用了五次多项式 $y = \sum^5_{j=0} \theta_jx^j$来进行拟合，虽然这个拟合曲线完美地通过了所有当前数据集中的数据，但这个曲线不是一个合理的预测工具（完全符合已知数据，但很可能不符合预测数据）。

换句话说，最左边的图像就是一个 欠拟合（under fitting） 的例子，比如明显能看出拟合的模型漏掉了数据集中的结构信息；而最右边的图像就是一个 过拟合（over fitting） 的例子。因此一个学习算法要保证能良好运行，特征的选择是非常重要的。

在原始版本的线性回归算法中，要对一个查询点 $x$ 进行预测 $h_\theta(x)$ 的输出，步骤如下：

1. 使用参数 $\theta$ 进行拟合，让数据集中的值与拟合算出的值的差值平方$(y^{(i)} − \theta^T x^{(i)} )^2$最小（最小二乘法的思想）
2. 输出$\theta^T x$

局部加权线性回归 方法假设有足够多的训练数据，对不太重要的特征进行一些筛选，步骤如下：

1. 使用参数 $\theta$ 进行拟合，让加权距离 $w^{(i)}(y^{(i)} − \theta^T x^{(i)} )^2$ 最小
2. 输出$\theta^T x$

对于权值的选取可以使用下面这个比较标准的公式：

$w^{(i)} = exp(− \frac {(x^{(i)}-x)^T(x^{(i)}-x)}{2\tau^2})$

随着点$x^{(i)}$ 到查询点 x 的距离降低，训练样本的权值的也在降低，带宽参数 $\tau$ 控制了这个降低的速度。

直观地说，式子想表达的意思是如果 $|x^{(i)} − x|$ 非常小，那么权值 $w^{(i)}$ 就接近 1；反之如果 $|x^{(i)} − x|$ 非常大，那么权值 $w^{(i)}$ 就变小。换句话说，查询点 $x$ 附近的训练样本有更高得多的权值。

局部加权线性回归算法 是一个 非参数（non-parametric） 算法：使用假设 $h$ 预测时，需要依赖训练集数据，也就是说我们需要一直保留着一定规模的训练集（整个训练集）。

而之前无权重的线性回归算法就是一种 参数（parametric） 算法： 使用假设 $h$ 预测时，只需要保存固定的有限个数的参数 $\theta$ 。因为通过训练集进行学习之后，就只需要保存参数 $\theta$，就不需要再保留训练数据样本了。

逻辑回归（Logistic Regression）

分类问题 ：输入训练集，输出一个合适的 h 函数（假设，hypothesis），并且这个函数的输出是离散值。

分类问题的第一个问题就是二值化分类问题，即函数输出 $y$ 只有两个取值：0或者1

例如，假如要建立一个垃圾邮件筛选器，那么就可以用$x^{(i)}$表示一个邮件中的若干特征，然后如果这个邮件是垃圾邮件，y 就设为1，否则 y 为 0。

我们也可以还按照之前的线性回归的模型来根据给定的 $x$ 来预测 $y$，只要忽略掉 $y$ 是一个散列值就可以了。然而这样构建的方法运行效率会非常低，效果很差。而且从直观上来看，$h_\theta(x)$ 的值如果大于1或者小于0就都没有意义了，因为咱们已经实现都确定了$y \in \{0,1\}$，就是说 $y$ 必然应当是 0 和 1 这两个值当中的一个。

因此，逻辑回归模型使用了另一种假设函数 $h_\theta (x)$ 的形式，来解决这个问题：

$h_\theta(x) = g(\theta^T x) = \frac 1{1+e^{-\theta^Tx}}$

$g(z)= \frac 1 {1+e^{-z}}$

这个 $g(z)$ 函数叫做逻辑函数 （logistic function）或双弯曲S型函数（sigmoid function）。

当$z\to +\infty$的时候 $g(z)$ 趋向于1，而当$z\to -\infty$时 $g(z)$ 趋向于0：

$g$ 函数当然也可以使用其他从0到1之间光滑递增的函数，但是后面（下文广义线性模型处）我们会知道选择这个逻辑函数是很合理的。

之前已经看到了在一系列统计学假设的前提下，最小二乘法回归可以通过最大似然估计来推出。

那么接下来就给逻辑回归模型做一系列的统计学假设，首先假设：

$P(y=1|x;\theta)=h_{\theta}(x)$

$P(y=0|x;\theta)=1- h_{\theta}(x)$

更简洁的写法是：

$p(y|x;\theta)=(h_\theta (x))^y(1- h_\theta (x))^{1-y}$

假设 $m$ 个训练样本都是各自独立生成的，那么 似然函数 就可以写成：
$\begin{aligned} L(\theta) = p(\vec{y}| X; \theta)\ = \prod^m_{i=1} p(y^{(i)}| x^{(i)}; \theta)\ = \prod^m_{i=1} (h_\theta (x^{(i)}))^{y^{(i)}}(1-h_\theta (x^{(i)}))^{1-y^{(i)}} \ \end{aligned}$

跟之前一样，取个对数更容易计算似然函数的最大值：
$\begin{aligned} l(\theta) =\log L(\theta) \ = \sum^m_{i=1} y^{(i)} \log h(x^{(i)})+(1-y^{(i)})\log (1-h(x^{(i)})) \end{aligned}$

之后就需要通过 最大似然法 找到合适的参数 $\theta$ 来尽可能让 $l(\theta)$ 最大化，实际上就相当于让似然函数 $L(\theta)$ 最大化。

感知器学习算法（The Perceptron Learning Algorithm）

假如对逻辑回归方法修改一下，“强迫”它输出的值要么是 0 要么是 1。

要实现这个目的，只需要将逻辑回归里的函数 $g$ 的定义修改一下，改成一个阈值函数（threshold function）：

$g(z)= \begin{cases}1,&{if \ z \geq 0 } \\ 0,&{if \ z < 0 }\end{cases}$

像之前一样令 $h_\theta(x) = g(\theta^T x)$，但用刚刚上面的阈值函数作为 $g$ 的定义，然后如果我们用了下面的更新规则：

$\theta_j := \theta_j +\alpha(y^{(i)}-h_\theta (x^{(i)}))x_j^{(i)}$

这样我们就得到了感知器学习算法。在 1960 年代，这个“感知器（perceptron）”被认为是对大脑中单个神经元工作方法的一个粗略建模。

注意：虽然感知器学习算法可能看上去表面上跟我们之前讲的其他算法挺相似，但实际上这是一个和逻辑回归以及最小二乘线性回归等算法在种类上都完全不同的算法；尤其重要的是，很难对感知器的预测赋予有意义的概率解释，也很难作为一种最大似然估计算法来推出感知器学习算法。

Softmax 回归（Softmax Regression）

分类问题的第二个问题是，函数输出 $y$ 的取值可以是 $k$ 个离散整数中的一个，也就是$y ∈ \{1, 2, ..., k\}$。

例如，我们这次要进行的分类就比把邮件分成垃圾邮件和正常邮件两类这种二值化分类要更加复杂一些，比如可能是要分成三类，例如垃圾邮件、个人邮件、工作相关邮件。这样响应变量依然还是离散的，但取值就不只有两个了。

Softmax 回归就是解决这类问题的模型，但是由于其模型推导复杂（学完广义线性模型后才推导出），就不多说明。

结果上，Softmax 回归的 似然函数 的对数将是：

$\begin{aligned} l(\theta) = \sum^m_{i=1} \prod ^k_{l=1}(\frac{e^{\theta_l^Tx^{(i)}}}{\sum^k_{j=1} e^{\theta_j^T x^{(i)}}})^{1\{y^{(i)}=1\}}\ \end{aligned}$

之后就需要通过 最大似然法 找到合适的参数 $\theta$ 来尽可能让 $l(\theta)$ 最大化。

寻找最优参数 $\theta$

在对某个问题建立好了学习模型之后，我们就需要对似然函数 $l(\theta)$ 进行最大化并找到对应的参数 $\theta$ ，也就是最大似然估计。

下面几个就是寻找参数 $\theta$ 的算法。

梯度下降/上升法（Gradient Descent/Gradient Ascent）

梯度下降/上升法（Gradient Descent/Gradient Ascent）的核心思想是：进行若干次迭代，每次迭代求出某个函数的梯度，然后让 $\theta$ 沿着梯度下降/上升的方向迈出一步并更新 $\theta$ 。

对于线性回归模型的最小二乘成本函数 $J(\theta)$ 来说，这需要的是梯度下降法。

用于 $J(\theta)$ 的梯度下降法就属于最小均方算法（LMS Algorithm）：

$\theta \leftarrow \theta - \alpha \cdot \nabla_{\theta} J(\theta) = \theta - \alpha [\frac{\partial J(\theta)}{\partial \theta_{0}},\frac{\partial J(\theta)}{\partial \theta_{1}},...,\frac{\partial J(\theta)}{\partial \theta_{n}}]^T$

也就是说，对于 $\theta$ 的每个分量 $i$ 都可以有：

$\theta_{i} \leftarrow \theta_{i} - \alpha \frac{\partial J(\theta)}{\partial \theta_{i}} = \theta_{i} - \alpha \sum_{j=1}^{m}(h_\theta(x^{(j)})-y^{(j)})x^{(j)}_i$

其中，$\alpha$ 为学习率。如果学习率比较小，那么将会缓慢收敛，需要更多次迭代；如果学习率比较大，虽然可以更快速移动，但 $J(\theta)$ 有可能在之后每次迭代中没有单调减少，也就是说可能不会收敛。

更新 $\theta$ 后便可以得到更小的成本函数 $J(\theta)$，也就是说，更新 $\theta$ 后的 $h_\theta(x)$ 越贴近 $y$ 。

如图，$\theta$ 为二维向量，梯度下降法进行一次迭代相当于从起始点（红色半山腰的十字点）每次沿下降方向走一步，就向流水一样，最终经过若干次迭代后，流向了 $J(\theta)$ 的低洼处（逼近局部最低点）：

类似地，求函数最大值可以梯度上升法（Gradient Ascent），在前面的几个模型为了最大化似然函数 $l(\theta)$ ，就可使用梯度上升法。梯度上升法里，更新方程中用的是加号而不是减号，因为我们现在是在找一个函数的最大值，而不是找最小值了：

$\theta := \theta +\alpha \nabla _\theta l(\theta)$

例如，用于逻辑回归的随机梯度上升法如下：

$g$ 函数的导数：$\begin{aligned} {g}'(z) = \frac d{dz} \times \frac 1{1+e^{-z}}\ = \frac 1{(1+e^{-z})^2} \times (e^{-z}))\ = \frac 1{(1+e^{-z})^2} \times (1- \frac 1{(1+e^{-z})})\ & = g(z)(1-g(z))\ \end{aligned}$

那么关于$\theta_j$ 的似然函数 $l(\theta)$ 的导数为：

$\begin{aligned} \frac{\partial}{\partial \theta_{j}} l(\theta) &=\left(y \frac{1}{g\left(\theta^{T} x\right)}-(1-y) \frac{1}{1-g\left(\theta^{T} x\right)}\right) \frac{\partial}{\partial \theta_{j}} g\left(\theta^{T} x\right) \\ &=\left(y \frac{1}{g\left(\theta^{T} x\right)}-(1-y) \frac{1}{1-g\left(\theta^{T} x\right)}\right) g\left(\theta^{T} x\right)\left(1-g\left(\theta^{T} x\right)\right) \frac{\partial}{\partial \theta_{j}} \theta^{T} x \\ &=\left(y\left(1-g\left(\theta^{T} x\right)\right)-(1-y) g\left(\theta^{T} x\right)\right) x_{j} \\ &=\left(y-h_{\theta}(x)\right) x_{j} \end{aligned}$

然后就得到了随机梯度上升规则：

$\theta_j := \theta_j + \alpha (y^{(i)}-h_\theta (x^{(i)}))x_j^{(i)}$

实际上，梯度下降/上升法又可以细分成几种迭代方式。

为了减少篇幅，下面都只介绍适用于 $J(\theta)$ 的下降方式，但实际上只要稍微改变，就可应用于上升方式。

批量梯度下降法（Batch Gradient Descent）

批量梯度下降法的一次迭代需要应用训练集所有数据样本，每次迭代都将得到几乎最好的梯度方向，因此所需迭代次数较少，但是缺点是当训练集的规模$m$很大时，一次迭代的开销将是巨大的。

批量梯度下降公式（一次迭代）：

$\theta_{i} \leftarrow \theta_{i} - \alpha \sum_{j=1}^{m}(h_\theta(x^{(j)})-y^{(j)})x^{(j)}_i$

随机梯度下降法（Stochastic Gradient Descent）

随机梯度下降，也叫增量梯度下降法（Incremental Gradient Descent）。

一次迭代只需要应用训练集的一个数据样本，所需迭代次数较多（一般与训练集的数据样本数一致）。每次迭代的效果没有考虑整体数据样本，因此每次迭代的梯度方向看起来有些随机（但整体还是朝正确方向的），随机梯度下降的迭代效果几乎总是稍差于批量梯度下降的一次迭代，但是这种方法的迭代开销极小，适用于大规模训练集。

随机梯度下降运行若干次迭代后也有可能无法收敛，这时候${\theta}$会一直在$J_{(\theta)}$最小值附近震荡；一个可选的解决方式是，随着迭代次数的增加，我们可以逐渐让学习速率$\alpha$逐渐趋变小，这样也能保证我们最后得到的参数会收敛到最小值，而不是在某个范围内一直在震荡。不过通常情况下，在某个范围附近震荡的值大多数其实也足够逼近了最小值了。

随机梯度下降公式（一次迭代，取第 j 个数据样本）：

$\theta_i \leftarrow \theta_i -\alpha(h_{\theta}(x^{(j)})-y^{(j)})x_i^{(j)}$

小批次梯度下降（Mini-Batch Gradient Descent）

实际上很多时候，批量梯度下降法和随机梯度下降法都不能满足庞大训练集的需求（前者迭代开销过大，后者所需迭代次数过多），而最佳方式是这两个方法的折中，即小批次梯度下降法。

所谓小批次梯度下降法，就是将数据拆分成很多小批次，每个批次取具有相同数量（$b$个）的数据样本，在每次迭代时取一小批次作为一次批量梯度算法的全部样本进行计算。这种方式的好处是，迭代开销并不会很大（远小于批量梯度算法），但是一次迭代已经有足够好的收敛效果（比随机梯度算法好得多，因为考虑多了一些数据样本）

小批次梯度下降公式（一次迭代，取第 j 个到 第 j+b-1 个数据样本）：

$\theta_{i} \leftarrow \theta_{i} - \alpha \sum_{k=0}^{b-1}(h_\theta(x^{(j+k)})-y^{(j+k)})x^{(j+k)}_i$

正规方程（The Normal Equations）

寻找 $\theta$ 的一种方法，不用迭代计算，而是直接利用公式计算即可得到向量$\theta$。

推导步骤：

1. 得到矩阵$X_{m×(n+1)}$（m个数据，每个数据有n个特征）：$X = X_{m×(n+1)}=[x_1,...,x_m]^T$和向量$y$（m个标签）：$y$
2. 令矩阵导数$\nabla_{\theta} J(\theta) = 0$，即相当于求驻点：
   $\nabla_{\theta} J(\theta) = \nabla_{\theta}\frac{1}{2}(X\theta-y)^T(X\theta-y) = X^TX\theta-X^Ty = 0$
3. 得到向量$\theta = (X^TX)^{-1}X^Ty$

正规方程公式：

$\theta = (X^TX)^{-1}X^Ty$

假如$X^TX$不可逆，这通常意味着多余的特征（线性相关的特征），这时候需要考虑删减掉数据的某个或某些特征。

与梯度下降法相较：

• 优点：无需选取学习速率 α ，无需进行迭代
• 缺点：需要求$X^TX$的逆矩阵，当$x$输入特征数量 n 较多时复杂度相当高 $O(n^3)$，一般只适用于特征量n较小的情况（小于10000）。

牛顿法（Newton's Method）

寻找 $\theta$ 的一种方法，它的思想史基于牛顿法（牛顿迭代法解决的是求一个方程零点的问题，即找出一个 $\theta$ ，来尽可能接近使 $f(\theta_{correct}) = 0$ 的 $\theta_{correct}$ ）

牛顿法对$\theta$进行如下的更新，来让 $\theta$ 不断逼近 $f(\theta_{correct}) = 0$ 的 $\theta_{correct}$ ：

$\theta \leftarrow \theta - \frac {f(\theta)}{f'(\theta)}$

我们可以把它理解成用一个线性函数来对函数 $f$ 进行逼近，这条直线是 $f$ 的切线，而猜测值是$\theta$，解的方法就是找到线性方程等于零的点，把这一个零点作为$\theta$设置给下一次猜测，然后依次类推。

而似然函数 $l(\theta)$ 的最大值的点应该对应着是它的导数 $l′(\theta) = 0$ 的点，所以通过令$f(\theta) = l′(\theta)$，就可以同样用牛顿法来找到 $l(\theta)$ 的最大值，然后得到下面的更新规则：

$\theta \leftarrow \theta - \frac {l'(\theta)}{l''(\theta)}$

不过，$\theta$ 往往是一个向量，所以要对牛顿法扩展到多维情况，也叫牛顿-拉普森法（Newton-Raphson method）：

$\theta \leftarrow \theta - H^{-1}\nabla_\theta l(\theta)$

其中， $H$（Hessian矩阵）是一个 $(n + 1)\times (n + 1)$ 矩阵，其详细定义是：

$H_{ij}= \frac {\partial^2 l(\theta)}{\partial \theta_i \partial \theta_j}$

牛顿法通常都能比（批量）梯度下降法收敛得更快，是二阶收敛（如果第k次的误差是0.01，则第k+1次迭代误差将会变为0.0001），达到最小值所需要的迭代次数少很多。然而，牛顿法中的单次迭代开销往往要比梯度下降法付出更多的性能开销，因为要查找和转换一个 $n\times n$的 Hessian 矩阵（不过只要 $n$ 不大，牛顿法通常就还是更快一些），算法复杂度为 $O(n^2)$。

广义线性模型（General Liner Model）

到目前为止，我们看过了回归的案例，也看了分类的案例：

• 在回归的案例中，概率模型使用了高斯分布模型
• 在二值分类的案例中，概率模型使用了伯努利模型
• 在多值分类的案例中，概率模型使用了多项分布模型

其实这三种方法都是一个更广泛使用的模型的特例，这种更广泛使用的模型就叫做广义线性模型。广义线性模型还可以推出其他模型以应用到其他的分类和回归问题上。

指数族（The Exponential Family）

在学习 GLMs 之前，我们要先定义一下指数族分布（exponential family distributions）。

如果一个分布能用下面的方式来写出来，我们就说这类分布属于指数族：
$p(y;\eta) = b(y) exp(\eta^TT(y)-a(\eta))$

$η$ 叫做此分布的自然参数（natural parameter，也叫典范参数 canonical parameter）；$T(y)$叫做充分统计量（sufficient statistic），我们目前用的这些分布中通常$T (y) = y$；而$a(\eta)$是一个对数分割函数（log partition function）；$e^{−a(\eta)}$这个量本质上扮演了归一化常数（normalization constant）的角色，也就是确保 $p(y; \eta)$ 的总和或者积分等于1。

使用不同的 $T$, $a$ 和 $b$ 函数，就相当于定义了一个用 $η$ 进行参数化的分布族（family，或者叫集 set）。通过改变 $\eta$，我们就能得到这个分布族中的不同分布，伯努利分布和高斯分布就都属于指数分布族。

高斯分布（Gaussian）

在推导普通最小二乘线性回归的时候，$\sigma%^2$ 的值对我们最终选择的 $\theta$ 和$h_\theta(x)$都没有影响。所以我们可以给$\sigma^2$取一个任意值。

为了简化推导过程，就令$\sigma^2 = 1$。然后就有了下面的等式：

$\begin{aligned} p(y;\mu) &= \frac 1{\sqrt{2\pi}} exp (- \frac 12 (y-\mu)^2) \\ &= \frac 1{\sqrt{2\pi}} exp (- \frac 12 y^2) \times exp (\mu y -\frac 12 \mu^2) \ \end{aligned}$

如果我们把$\sigma^2$留作一个变量，高斯分布就也可以表达成指数分布的形式，其中$\eta ∈ R^2$就是一个二维向量，同时依赖$\mu$和$\sigma$。然而，对于广义线性模型GLMs方面的用途，$\sigma^2$参数就也可以看成是对指数分布族的更泛化的定义：$p(y; \eta, \tau ) = b(a, \tau ) exp((\eta^T T (y) − a(\eta))/c(\tau))$。这里面的$\tau$叫做分散度参数（dispersion parameter），对于高斯分布，$c(\tau) = \sigma^2$；不过上文中我们已经进行了简化，所以针对我们要考虑的各种案例，就不需要再进行更加泛化的定义了。

这样，我们就可以看出来高斯分布是属于指数分布族的，可以写成下面这样：

$\begin{aligned} \eta &= \mu \ T(y) \\ &= y \ a(\eta) \\ &= \mu ^2 /2\ \\ &= \eta ^2 /2\ b(y) \\ &= (1/ \sqrt {2\pi })exp(-y^2/2) \end{aligned}$

伯努利分布（Bernoulli）

伯努利分布的均值是 $\phi$，也写作 $Bernoulli(\phi)$，确定的分布是 $y ∈ \{0, 1\}$，因此有

$p(y = 1; \phi) = \phi$

$p(y = 0;\phi) = 1−\phi$

我们这样来写伯努利分布：

$\begin{aligned} p(y;\phi) &= \phi ^y(1-\phi)^{1-y}\ \\ &= exp(y \log \phi + (1-y)\log(1-\phi))\ \\ &= exp( (log (\frac {\phi}{1-\phi}))y+\log (1-\phi) )\ \end{aligned}$

这样，就能得到自然参数（natural parameter）$\eta$ ，即 $\eta = log (\frac \phi {1 − \phi})$。

如果翻转这个定义，用 $η$ 来解 $\phi$ 就会得到 $\phi = \frac{1}{(1 + e^{−\eta} )}$，这正好就是之前见到过的逻辑函数（logistic function）。

把逻辑回归作为一种广义线性模型（GLM）的时候还会遇到这个情况。

$\begin{aligned} T(y) &= y \ a( \eta) \\ &= - \log (1- \phi) \\ &= \log {(1+ e^ \eta)}\ b(y) \\ &=1 \end{aligned}$

上面这组式子就表明了伯努利分布可以写成指数族的形式，对应一组 $T$，$a$ 和 $b$ 函数。

其它分布

指数分布族里面还有很多其他的分布：

• 多项式分布（multinomial），用于对多值分类变量进行建模
• 泊松分布（Poisson），用于对计数类数据进行建模
• $\gamma$和指数分布（the gamma and the exponential），这个用于对连续的、非负的随机变量进行建模，例如时间间隔
• $\beta$和狄利克雷分布（the beta and the Dirichlet），这个是用于概率的分布

在接下来的广义线性模型（一个通用建模方案），其中的 $y$ (给定 $x$ 和 $\theta$) 可以是上面这些分布中的任意一种。

构建广义线性模型（Constructing GLMs）

假如你要构建一个模型，来估计在给定的某个小时内来到你商店的顾客人数（或者是你的网站的页面访问次数），基于某些确定的特征 $x$，例如商店的促销、最近的广告、天气、今天周几啊等等。而泊松分布（Poisson distribution）通常能适合用来对访客数目进行建模，而且也属于指数分布族的一个分布。知道了这个分布之后，怎么来建立一个模型来解决这个具体问题呢？

我们可以针对这类问题构建一个 广义线性模型（Generalized Linear Model，缩写为 GLM）。

进一步泛化，设想一个分类或者回归问题，要预测一些随机变量 $y$ 的值，作为 $x$ 的一个函数。要导出适用于这个问题的广义线性模型，就要对我们的模型、给定 $x$ 下 $y$ 的条件分布来做出以下三个假设：

(1) $y | x; \theta ∼ Exponential Family(\eta)$，即给定 $x$ 和 $\theta$, $y$ 的分布属于指数分布族，是一个参数为 $η$ 的指数分布。

(2) 给定 $x$，目的是要预测对应这个给定 $x$ 的 $T(y)$ 的期望值。由于这里用到的绝大部分情况都是 $T(y) = y$，这就意味着预测值 $h_\theta(x)$ 要满足 $h_\theta(x) = E[y|x]$。

$E[y|x$] 的意思：对给定 $x$ 时的 $y$ 值的期望。

(3) 自然参数 $η$ 和输入值 $x$ 是线性相关的，$\eta = \theta^T x$，或者如果 $\eta$ 是有值的向量，则有 $\eta_i = \theta_i^T x$。

实际上，第3个假设不太好证明，因此可以看作是一个我们在设计广义线性模型时候的一种 设计选择，而不是一个假设。

也就是说，广义线性模型 GLM 就是通过这三个假设或设计，根据关于 $y$ 的分布的不同类型的建模需求，来推导出一个合适的学习算法。这样的算法通常会有很多优点，例如学习效率高、效果更好。

构建普通最小二乘模型

普通最小二乘线性回归实际上是广义线性模型中的一种实例。

背景设置：目标变量 $y$（在广义线性模型的术语也叫做响应变量response variable）是连续的，然后我们将给定 $x$ 的 $y$ 的分布以 高斯分布 $N(\mu, \tau^2)$来建模，其中 $\mu$ 可以是依赖 $x$ 的一个函数。

在把高斯分布写成指数分布族的分布的时候，会有：

$\mu = \eta$

所以就能得到下面的等式：

$\begin{aligned} h_\theta(x) &= E[y|x;\theta] \\ &= \mu \\ &= \eta \\ &= \theta^Tx\ \end{aligned}$

第一行的等式是基于假设(2)；第二个等式是基于定理当$y|x; \theta ∼ N (\mu, \sigma ^2)$，则 y 的期望就是$\mu$；第三个等式是基于假设(1)，以及之前我们此前将高斯分布写成指数族分布的时候推导出来的性质$\mu = \eta$；最后一个等式就是基于假设(3)。

所以，假设函数 $h$ 的形式为：

$\begin{aligned} h_\theta(x) = \theta^Tx\ \end{aligned}$

构建逻辑回归模型

逻辑回归也是广义线性模型中的一种实例。

背景设置：要解决的问题是二值化分类问题，也就是 $y \in \{0, 1\}$，那么很自然就选择 伯努利分布 $Bernoulli(\phi)$ 来对给定 $x$ 的 $y$ 的分布进行建模了。

把伯努利分布写成一种指数族分布时，会有：

$\phi = \frac{1}{(1 + e^{−\eta})}$

另外，如果有 $y|x; \theta ∼ Bernoulli(\phi)$，那么 $E [y|x; \theta] = \phi$。

所以就跟刚刚推导普通最小二乘法的过程类似，有以下等式：

$\begin{aligned} h_\theta(x) &= E[y|x;\theta] \\ &= \phi \\ &= \frac{1}{(1 + e^{−\eta})} \\ &= \frac{1}{(1 + e^{−\theta^T x})}\ \end{aligned}$

所以，假设函数 $h$ 的形式为：

$h_\theta(x) = \frac{1}{(1 + e^{−\theta^T x})}$。

逻辑回归为什么会使用 $\frac{1}{(1 + e^{−z})}$ 这个逻辑函数，上面就是一种解答：一旦我们假设以 $x$ 为条件的 $y$ 的分布是伯努利分布，那么根据广义线性模型和指数分布族的定义，就会得出这个式子。

再介绍一些术语：这里给出分布均值的函数 $g$ 是一个关于自然参数的函数，$g(\eta) = E[T(y); \eta]$，这个函数也叫做规范响应函数（canonical response function），它的反函数$g^{−1}$叫做规范链接函数（canonical link function）。注：不过很多教科书用$g$表示链接函数，而用反函数$g^{−1}$来表示响应函数；但是这里为了方便说明，用的是反过来的，这也是继承了早期的机器学习中的用法。

对于高斯分布来说，它的规范响应函数正好就是识别函数（identify function）；而对于伯努利分布来说，它的规范响应函数则是逻辑函数（logistic function）。

构建Softmax回归模型

Softmax回归也是广义线性模型的一种实例。

背景设置：要解决的问题是多值分类问题，因此使用 多项式分布 来对给定 $x$ 的 $y$ 的分布进行建模。

要实现这一目的，首先还是要把多项式分布也用指数族分布来进行描述。

要对一个可能有 k 个不同输出值的多项式进行参数化，就可以用 k 个参数 $\phi_1,...,\phi_ k$ 来对应各自输出值的概率。不过这么多参数可能太多了，形式上也太麻烦，他们也未必都是互相独立的

比如对于任意一个$\phi_i$中的值来说，只要知道其他的$k-1$个值，就能知道这最后一个了，因为总和等于1，也就是$\sum^k_{i=1} \phi_i = 1$

所以就可以去掉一个参数，只用 $k-1$ 个参数：$\phi_1,...,\phi_ {k-1}$ 来对多项式进行参数化：

$\phi_i = p (y = i; \phi)$

$p (y = k; \phi) = 1 −\sum ^{k−1}{i=1}\phi i$

为了表述起来方便，我们还要设 $\phi_k = 1 − \sum_{i=1}^{k−1} \phi_i$

但一定要注意，$\phi_k$ 并不是一个参数，而是完全由其他的$k-1$个参数来确定的。

要把一个多项式表达成为指数组分布，还要按照下面的方式定义一个 $T (y) ∈ R^{k−1}$:

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

这次和之前的样例都不一样了，就是不再有$T(y) = y$；然后，$T(y)$现在是一个$k – 1$维的向量，而不是一个实数了。向量$T(y)$中的第 i 个元素写成$(T(y))_i$。

现在介绍一种非常有用的记号。指示函数（indicator function）$1\{\cdot\}$，如果参数为真，则等于1；反之则等于0。

$1\{True\} = 1, 1\{False\} = 0$，例如，1{2 = 3} = 0, 而1{3 = 5 − 2} = 1。

所以我们可以把 $T (y)$ 和 $y$ 的关系写成 $(T(y))_i = 1\{y = i\}$。

在此基础上，就有了 $E[(T(y))i] = P (y = i) = \phi_i$。

现在一切就绪，可以把多项式写成指数族分布：

$\begin{aligned} \mathrm{p}(\mathrm{y} ; \phi) &=\phi_{1}^{1\{y=1\}} \phi_{2}^{1\{y=2\}} \ldots \phi_{\mathrm{k}}^{1\{\mathrm{y}=\mathrm{k}\}} \\ &=\phi_{1}^{1\{\mathrm{y}=1\}} \phi_{2}^{1\{\mathrm{y}=2\}} \ldots \phi_{\mathrm{k}}^{1-\sum_{\mathrm{k}}^{\mathrm{k}-1} 1\{\mathrm{y}=\mathrm{i}\}} \\ &=\phi_{1}^{(\mathrm{T}(\mathrm{y}))}{ } \phi_{2}^{(\mathrm{T}(\mathrm{y}))_{2}} \ldots \phi_{\mathrm{k}}^{1-\sum_{\mathrm{i}=1}^{\mathrm{k}-1}(\mathrm{~T}(\mathrm{y}))_{\mathrm{i}}} \\ &=\exp \left((\mathrm{T}(\mathrm{y}))_{1} \log \left(\phi_{1}\right)+(\mathrm{T}(\mathrm{y}))_{2} \log \left(\phi_{2}\right)+\cdots+\left(1-\sum_{\mathrm{i}=1}^{\mathrm{k}-1}(\mathrm{~T}(\mathrm{y}))_{\mathrm{i}}\right) \log \left(\phi_{\mathrm{k}}\right)\right) \\ &=\exp \left((\mathrm{T}(\mathrm{y}))_{1} \log \left(\frac{\phi_{1}}{\phi_{\mathrm{k}}}\right)+(\mathrm{T}(\mathrm{y}))_{2} \log \left(\frac{\phi_{2}}{\phi_{\mathrm{k}}}\right)+\cdots+(\mathrm{T}(\mathrm{y}))_{\mathrm{k}-1} \log \left(\frac{\phi_{\mathrm{k}-1}}{\phi_{\mathrm{k}}}\right)+\log \left(\phi_{\mathrm{k}}\right)\right) \\ &=\mathrm{b}(\mathrm{y}) \exp \left(\eta^{\mathrm{T}} \mathrm{T}(\mathrm{y})-\mathrm{a}(\eta)\right) \end{aligned}$

其中：

$\eta=\left[\begin{array}{c}\log \left(\phi_{1} / \phi_{\mathrm{k}}\right) \\ \log \left(\phi_{2} / \phi_{\mathrm{k}}\right) \\ \vdots \\ \log \left(\phi_{\mathrm{k}-1} / \phi_{\mathrm{k}}\right)\end{array}\right]$
$\mathrm{a}(\eta)=-\log \left(\phi_{\mathrm{k}}\right)$
$\mathrm{b}(\mathrm{y})=1$

这样咱们就把多项式方程作为一个指数族分布来写了出来。

与$i (for : i = 1, ..., k)$对应的链接函数为：

$\eta_i =\log \frac {\phi_i}{\phi_k}$

为了方便起见，我们再定义$\eta_k = \log (\phi_k/\phi_k) = 0$。对链接函数取反函数然后推导出响应函数，就得到了下面的等式：

$\begin{aligned} e^{\eta_i} &= \frac {\phi_i}{\phi_k}\ \\ \phi_k e^{\eta_i} &= \phi_i \ \\ \phi_k \sum^k_{i=1} e^{\eta_i} &= \sum^k_{i=1}\phi_i= 1\ \end{aligned}$

这就说明了$\phi_k = \frac 1 {\sum^k_{i=1} e^{\eta_i}}$，然后可以把这个关系代入回到上面第二个等式，这样就得到了响应函数：

$\phi_i = \frac { e^{\eta_i} }{ \sum^k_{j=1} e^{\eta_j}}$

上面这个函数从 $η$ 映射到了 $\phi$，称为 Softmax 函数。

要完成我们的建模，还要用到前文提到的假设(3)，也就是$\eta_i$是一个 $x$ 的线性函数。所以就有了$\eta_= \theta_i^Tx (for:i = 1, ..., k − 1)$，其中的$\theta_1, ..., \theta_{k−1} ∈ R^{n+1}$就是我们建模的参数。

为了表述方便，我们这里还是定义 $\theta_k = 0$，这样就有 $\eta_k = \theta_k^T x = 0$，跟前文提到的相符。

因此，我们的模型假设了给定 $x$ 的 $y$ 的条件分布为：

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

也就是说，我们的假设函数会对每一个$i = 1,...,k$，给出$p (y = y|x; \theta)$概率的估计值。

$\begin{aligned} \mathrm{h}_{\theta}(\mathrm{x}) &=\mathrm{E}[\mathrm{T}(\mathrm{y}) \mid \mathrm{x} ; \theta] \\ &=\mathrm{E}\left[\left| \begin{array}{c}1(\mathrm{y}=1) \\ 1(\mathrm{y}=2) \\ \vdots \\ 1(\mathrm{y}=\mathrm{k}-1)\end{array} \right| \mathrm{x} ; \theta\right] \\ &=\mathrm{E}\left[\begin{array}{c}\phi_{1} \\ \phi_{2} \\ \vdots \\ \phi_{\mathrm{k}-1}\end{array}\right] \\ &=\mathrm{E}\left[\begin{array}{c}\frac{\exp \left(\theta_{1}^{\mathrm{T}} \mathrm{x}\right)}{\sum_{\mathrm{j}=1}^{\mathrm{k}} \exp \left(\theta_{\mathrm{j}}^{\mathrm{T}} \mathrm{x}\right)} \\ \frac{\exp \left(\theta_{2}^{\mathrm{T}} \mathrm{x}\right)}{\sum_{\mathrm{j}=1}^{\mathrm{k}} \exp \left(\theta_{\mathrm{j}}^T \mathrm{x}\right)} \\ \vdots \\ \frac{\exp \left(\theta^T_{k-1} \mathrm{x}\right)}{\sum_{\mathrm{j}=1}^{\mathrm{k}} \exp \left(\theta_{\mathrm{j}}^{\mathrm{T}} \mathrm{x}\right)}\end{array}\right] \end{aligned}$

虽然在前面假设的这个$h_\theta(x)$只有$k-1$维，但 $p (y = y|x; \theta)$ 可以通过用 $1− \sum^{k-1}_{i=1}\phi_i$ 来得到。

最后，与之前对普通最小二乘线性回归和逻辑回归的原始推导类似，要拟合这个模型的参数 $\theta$ ，可以先写出其似然函数的对数：

$\begin{aligned} l(\theta) =\sum^m_{i=1} \log p(y^{(i)}|x^{(i)};\theta)\ = \sum^m_{i=1} \prod ^k_{l=1}(\frac{e^{\theta_l^Tx^{(i)}}}{\sum^k_{j=1} e^{\theta_j^T x^{(i)}}})^{1\{y^{(i)}=l\}}\ \end{aligned}$

总结

广义线性模型解决监督学习问题（对于一个属于指数分布族的变量进行预测或者分类）：

1. 根据这个变量的特征找出这个变量的分布模型 $p(y|x)$
2. 利用分布模型得到机器学习的假设模型 $h_\theta (x) = E(y|x)$
3. 得到似然函数 $L(\theta)$
4. 通过极大化似然函数得到机器学习模型中的最佳参数 $\theta$，问题解决

初学ML，笔记仍在整理中...

参考

• [1] Stanford-CS229 by 吴恩达
• [2]《Python机器学习》