我们知道机器学习算法的实质就是模型表征+优化算法+模型评估三个部分，其中优化算法的工作就是从模型表征空间里找出模型评估指标最好的模型。虽然现有的优化算法已经集成到各类机器学习平台中，但是如果想要成为称职的算法工程师，从原理上了解优化算法必不可少。
[TOC]
1. 监督学习的损失函数
监督学习中损失函数作为模型评估的指标，定义了模型和数据的匹配程度。定义损失函数$L(\cdot , \cdot)=Y\times Y\to\mathbb{R}\ge0$，如果$L(f(x_i, \theta), y_i)$越小，则说明该模型在样本点中匹配的越好。
1.1 监督学习中的损失函数有哪些？
1.1.1 分类问题
对于二分类问题，$Y={1, -1}$，我们希望 $sign f(x_i, \theta)=y_i$，
-
0-1 Loss
最自然的是0-1损失，$L_{0-1}(f, y)=\mathbb{I}_{fy\leq0}$
其中，$\mathbb{I}_{p}$是指示函数(Indicator Function)，当且仅当$p$为真时取值为1，否则取之为0。0-1损失可以直观刻画分类错误率，但是因为其非凸非光滑，算法很难对该函数进行优化。
-
Hinge Loss
Hinge Loss是0-1 Loss相对饱和的凸上界, $L_{hinge}(f,y)=max(0, 1-fy)$
当$fy\ge1$时，不做任何惩罚。当$fy=1$时，Hinge函数不可导，所以不能用梯度下降法进行优化，而是使用次梯度下降法(Subgradient Descent Method)。
-
Logistic Loss
Logstic Loss也是0-1 Loss的凸上界，$L_{logistic}(f, y) = log_2(1+e^{-fy})$
优点是Logistic Loss处处光滑，适合用梯度下降进行优化，但是问题也是该函数会对所有样本点进行惩罚，所以对异常值会敏感一些。
-
Cross Entropy
当预测值在$f\in[-1,1]$的时候，另一个常用的损失函数是交叉熵。 $L_{cross\ entropy}(f, y)=-log_2\left(\frac{1-fy}{2}\right)$ 。交叉熵也是0-1Loss的凸上界。
下图给出上述四种函数的曲线形式
1.1.2 回归问题
对于回归问题，$y=\mathbb{R}$，我们希望$f(x_i, \theta)\simeq y_i$，
-
Square Loss
最常见的是平方损失函数，$L_{square}(f, y)=(f-y)^2$
平方损失函数是光滑函数，所以可以用梯度下降进行优化。当预测值离真实值越远时，惩罚力度越大。同样这也会导致它对异常点比较敏感。
-
Absolute Loss
为了从一定程度上解决平方损失函数的问题，可以采用绝对损失函数，$L_{absolute}(f, y)=\lvert f-y\rvert$。
绝对损失函数相当于做中值回归，而平方损失函数相当于做均值回归。因而绝对损失函数对异常点相对更加鲁棒，但是绝对损失函数在$f=y$这一点不连续。
-
Huber Loss
基于上面两种损失函数的优劣，我们可以推出Huber损失函数，$L_{Huber}=\begin{cases}(f-y)^2, &\lvert f-y\rvert\leq\delta\\ 2\delta\lvert f-y\rvert-\delta^2,&\vert f-y\vert\gt\delta\end{cases}$
容易看出，通过预先设定阈值$\delta$，可以当$\lvert f-y\rvert$较小的时候为平方损失，当$\lvert f-y\rvert$较大的时候为线性损失，函数整体处处可导，并且对异常值的容忍相对较高。
下面是三种回归的损失的图例
2. 机器学习中的优化问题
To be updated
3. 经典优化算法
假设现在有一道无约束优化问题：$min_\theta L(\theta)$
其中目标函数$L(\theta)$是光滑的。问求解这个问题的优化算法有哪些，他们的使用场景分别是什么。
经典的优化算法可以分为直接法和迭代法两类。
3.1 直接法
直接法要求目标函数需要满足两个条件：
-
第一个条件是要求$L(\cdot)$是凸函数。如果$L(\cdot)$是凸函数，那么$\theta^*$是最优解的充分必要条件就是$L(\cdot)$在$\theta^*$处的梯度为0，即 $\nabla L(\theta^*)=0$
-
问了能够求出上式中的$\theta^*$，第二个条件就是上式有闭式解。
满足上面两个条件的经典例子是岭回归(Ridge Regression)，其目标函数为 $L(\theta)=|X\theta-y|_2^2+\lambda|\theta|^2$
并可以推导出最优解为 $\theta^*=(X^TX+\lambda I)^{-1}X^Ty$
3.2 迭代法
迭代法就是迭代地修正对最优解的估计。因为直接法的条件要求限制了它的应用空间，所以很多实际问题会采用迭代法。
假设当前对最优解的估计值为$\theta_t$，希望求解最优化问题 $\delta_t={arg\ min}_\delta L(\theta_t+\delta)$
来得到更好的估计值 $\theta_{t+1}=\theta_t+\delta_t$
迭代法又分为一阶法和二阶法两类。
一阶法是对函数$L(\theta_t+\delta)$做一阶泰勒展开，得到近似式 $L(\theta_t+\delta)\simeq L(\theta_t)+\nabla L(\theta_t)^T\delta$
由于该近似式仅在$\delta$较小时才比较准确，所以在求解$\delta_t$时一般要加上L2正则项. \(\begin{aligned}L(\theta_t+\delta) &=arg\ min_\delta \big(L(\theta_t)+\nabla L(\theta_t)^T\delta + \frac{1}{2\alpha}\|\delta\|_2^2\big)\\&=-\alpha\nabla L(\theta_t)\end{aligned}\)
所以一阶法迭代的公式为 $\theta_{t+1}=\theta_t - \alpha\nabla L(\theta_t)$
其中，$\alpha$是学习率。一阶法也叫梯度下降法，梯度就是目标函数的一阶信息。
二阶法就是对函数$L(\theta_t+\delta)$做二阶泰勒展开，得到近似式 $\begin{equation}L(\theta_t+\delta)\simeq L(\theta_t)+\nabla L(\theta_t)^T\delta+\frac{1}{2}\delta^2\nabla^2 L(\theta_t)^T\delta\end{equation}$
其中，$\delta^2\nabla^2 L(\theta_t)$是函数$L$在$\theta_t$处的Hessian矩阵。通过求解近似优化问题，可得 \(\begin{aligned}L(\theta_t+\delta)&={arg\ min}_\delta(L(\theta_t)+\nabla L(\theta_t)^T\delta+\frac{1}{2}\delta^2\nabla^2 L(\theta_t)^T\delta)\\&=-\alpha\nabla^2L(\theta_t)^{-1}L(\theta_t)\end{aligned}\)
二阶法又称为牛顿法，Hessian矩阵就是目标函数的二阶信息。二阶法的收敛速度往往大于一阶法，但是当据维度过高，Hessian矩阵的求逆计算复杂度很大，而且当函数非凸的时候，二阶法有可能收敛到鞍点。
4. 梯度验证
给定目标函数$min_{\theta\in\mathbb{R}^n}L(\theta)$，如何根据求目标函数值的方式来验证求目标函数梯度的功能是否正确呢？
根据梯度的定义，目标函数的梯度为 $\nabla L(\theta)=\left[\frac{\partial L(\theta)}{\partial\theta_1},\dots,\frac{\partial L(\partial\theta)}{\theta_n}\right]^T$
其中对于$i=1,2,\dots,n$，梯度的第$i$个元素的定义为 $\frac{\partial L(\theta)}{\partial\theta_i}=lim_{h\to0}\frac{L(\theta+he_i)-L(\theta-he_i)}{2h}$
其中，$e_i$是单位向量，与$\theta$维度相同，仅在第$i$个元素取1，其余元素全为0。可以令$h$等于一个非常小的值，如$10^{-7}$，则 $\frac{\partial L(\theta)}{\partial\theta_i}\simeq\frac{L(\theta+he_i)-L(\theta-he_i)}{2h}$
可以利用泰勒展开近似计算误差，令单变量函数 $\tilde L(x)=L(\theta+xe_i)$，根据泰勒展开和拉格朗日余子式，
$L(\theta+he_i)=\tilde L(h)=\tilde L(0)+\tilde L’(0)h+\frac{1}{2}\tilde L’‘(0)h_2+\frac{1}{6}\tilde L^{(3)}(p_i)h^3$，其中，$p_i\in(0,h)$。
类似地，$L(\theta-he_i)=\tilde L(h)=\tilde L(0)-\tilde L’(0)h+\frac{1}{2}\tilde L’‘(0)h^2-\frac{1}{6}\tilde L^{(3)}(q_i)h^3$。其中，$q_i\in(-h,0)$。
上面两式相减，等式两边同时除以$2h$，由于$\tilde L(0)=\frac{\partial L(\theta)}{\partial\theta_i}$，可得 $\frac{L(\theta+he_i)-L(\theta-he_i)}{2h}=\frac{\partial L(\theta)}{\partial\theta_i}+\frac{1}{12}\left(\tilde L^{(3)}(p_i)+\tilde L^{(3)}(q-i)\right)h^2$。
当$h$充分小时，$p_i,q_i$都接近0，可以近似地认为$h^3$前面的系数是常数$M$，因此可以近似地认为误差为 $\left\lvert\frac{L(\theta+he_i)-L(\theta-he_i)}{2h}-\frac{\partial L(\theta)}{\partial\theta_i}\right\rvert=Mh^2$
当$h$很小的时候，每当$h$减小为原来的$10^{-1}$，近似误差都约减小原来的$10^{-2}$。可以近似地认为近似误差是$h$的高阶无穷小。
所以，实际应用中，我们随机初始化参数$\theta$，取较小的$h$，如$10^{-7}$，对$i=1,2,\dots,n$ 依次验证下式是否成立：$\left\lvert\frac{L(\theta+he_i)-L(\theta-he_i)}{2h}-\frac{\partial L(\theta)}{\partial\theta_i}\right\rvert\lt h$。
如果对于某个$i$，等式不成立，有以下可能：
- 改下标对应的$M$过大；
- 该梯度分量计算不正确。
此时可以固定$\theta$，减小$h$为原来的$10^{-1}$，再次计算下标$i$对应的近似误差，若近似误差约减小为原来的$10^{-2}$，则对应第一种可能，应该用更小的$h$重新做一次梯度验证；否则对应第二种可能，应该检查求梯度的代码是否有错误。
5. 随机梯度下降法
5.1当数据量特别大时，经典的梯度下降法会出现哪些问题，需要做哪些改进？
机器学习中，优化算法的目标函数通常可以写为 $L(\theta)=\mathbb{E}{(x,y)\sim P{data}}L(f(x,\theta), y)$。
其中，$\theta$表示模型的参数，$x$表示输入数据，$y$表示模型的目标输出，$P_{data}$表示数据分布，$f(x,\theta)$表示模型的实际输出，$L$表示模型在数据$(x,y)$上的损失，$\mathbb{E}$表示期望。$L(\theta)$表示模型在数据集上的平均损失。最优化的目标就是再到模型平均损失最小的参数$\theta^\ast$: $\theta^\ast=arg\space minL(\theta)$
经典的梯度下降法用所有训练数据的平均损失来近似目标函数，即
\[\begin{align} &L(\theta)=\frac{1}{M}\sum^M_{i=1}L(f(x_i,\theta), y_i), \\ &\nabla L(\theta)=\frac{1}{M}\sum^M_{i=1}\nabla L(f(x_i,\theta), y_i) \end{align}\]其中$M$是训练样本总数，根据上式可得参数更新的公式为(这里用到了直接法，一阶泰勒展开)：$\theta_{t+1}=\theta_t-\alpha\nabla L(\theta)$. 所以，梯度下降法在对参数进行更新的时候需要遍历整个训练数据的所有样本。当$M$很大的时候，时间和算力的开销很大，在实际应用场景里不实用。
为了解决这个问题，随机梯度下降法(Stochastic Gradient Decent, SGD)用单个训练样本的损失来近似平局损失，即
\[\begin{align} &L(\theta;x_i,y_i)=L(f(x_i,\theta), y_i) \\ &\nabla L(\theta;x_i,y_i)=\nabla L(f(x_i,\theta), y_i) \end{align}\]因此，SGD仅使用单个训练样本就可以对模型参数进行更新。一来加快了收敛速率，二来也可以处理样本数据不断更新的场景。
实际应用中，为了降低随机梯度的方差，提高稳定性，利用高度优化的矩阵计算操作，多个训练样本数据会被同时处理，这种方法被称为小批量梯度下降法(Mini-Batch Gradient Decent)。
假设需要同时处理m个训练样本数据${(x_1,y_1),(x_2,y_2),\dots,(x_m,y_m)}$，目标函数和其梯度公式为
\[\begin{align} &L(\theta)=\frac{1}{m}\sum^m_{j=1}L(f(x_{i_j},\theta), y_{i_j}), \\ &\nabla L(\theta)=\frac{1}{m}\sum^m_{j=1}\nabla L(f(x_{i_j},\theta), y_{i_j}) \end{align}\]对于Mini-Batch Gradient Decent有三点注意事项：
-
如何选取Batch Size m？
不同场景下最优的m往往不同，所以需要通过调参来选择。但是一般m取2的幂次可以充分利用矩阵运算操作。
-
如何选择m个训练数据？
为了避免数据特定顺序给训练带来影响，所以一般情况是每次遍历前，先对数据进行随机排序然后顺序选取m个样本，或者直接利用随机函数从样本数据总体中抽取m个样本。
-
如何选择学习率$\alpha$?
为了加快收敛速率，同时增加最后模型的精确率，通常需要应用衰减式的学习速率：一开始采用较大的学习率，当误差曲线进入平台期后，适当减小学习速率以获得更加精细的调整。最优的学习速率也需要通过调参获得。
所以综上，当遇到数据量较大的问题时采用Mini-Batch Gradient Decent可以大大提高训练的速率。
6. 随机梯度下降法的加速
6.1 深度学习中，SGD偶尔会失效，无法给出满意的结果，这是什么原因？
这里就通过对比批量梯度下降法(BGD)(传统的梯度下降法)和随机梯度下降法(SGD)来进行解释。
BGD在全部训练集${x_i,y_i}^n_{i=1}$上计算准确的精度，即 $\sum_{i=1}^n\nabla_\theta f(\theta;x_i,y_i)+\nabla_\theta \phi(\theta)$
其中$f(\theta;x_i,y_i)$是每个样本$(x_i,y_i)$的损失函数，$\phi(\theta)$为正则项
SGD则采样单个样本来估计当前梯度，即 $\nabla_\theta f(\theta;x_i,y_i)+\nabla_\theta \phi(\theta)$
结合上式和下图可以可以看出，BGD通过利用整个数据集，牺牲时间和内存上的性能来换取比较稳定的参数收敛轨迹；而SGD每步仅采用单个或者mini-batch的数据，提高计算速率，减少内存开销，问题是使参数收敛曲线抖动非常厉害，而且有停留在局部最优点的可能。
进一步说，深度学习中的陷阱很多。除了局部最优点，还有更可怕的是“山谷”和鞍点这两类地形。在“山谷”中，准确的梯度估计是沿着山道向下直到谷底，但是粗糙的梯度估计会使最优估计在山壁间来回震荡，导致收敛不稳定和收敛速度慢；鞍点是因为梯度比较平稳，接近0，使得收敛渐渐停止或者朝着相反的方向移动。
6.2 为了解决SGD的痛点，后来提出了哪些改进？分别有什么特点？
SGD的本质其实就是采用迭代的方法更新参数，更新公式为 $\theta_{t+1}=\theta_t-\eta g_t$
其中，当前估计的负梯度$-g_t$表示步子的方向，学习速率$\eta$控制步幅。改进的SGD仍然基于这样的更新公式。
6.2.1 SGD+Momentum
从名字上就可以直观地理解，加上动量(Momentum)，使当梯度为0时，优化方向也可以根据“惯性”继续向前移动，从而可以冲出局部最优解，或者鞍点。参数更新的公式为 $\begin{align}
&v_t=\gamma v_{t-1}+\eta g_t \
&\theta_{t+1}=\theta_t-v_t
\end{align}$
具体来说，前进速度$v_t$由两部分组成，一个是学习速率$\eta$乘当前估计的梯度$g_t$；另一个是带衰减系数$\gamma$的前一次速率$v_{t-1}$。
下图对比了带或不带动量的SGD收敛曲线。可以明显看出，带有动量的SGD方法收敛速度更快，曲线也更平滑，但是带有动量的方法往往会先超过最优点，然后再回来。
6.2.2 AdaGrad
动量的获得是基于历史信息，除此之外，还可以通过对参数空间环境的感知，自适应地确定学习速率，使不同参数具有不同的更新速率。比如训练词嵌入模型时，有的词或者词组频繁出现，而有的极少出现。数据的稀疏性导致相应参数的梯度的稀疏性。不频繁出现的词或者词组的参数的梯度在大多数情况下为0，从而这些参数被更新的频率很低。
实际应用中，我们希望更新频率低的参数可以具有较高的更新步幅，而更新频率高的参数可以减小更新步幅。Adagrad利用历史梯度平方和来衡量不同参数的梯度的稀疏性，取值越小表明越稀疏。所以参数具体更新公式为 $\theta_{t+1,i}=\theta_i-\frac{\eta}{\sqrt{\sum_{k=0}^tg_{k,i}^2+\epsilon}}g_{t,i}$
其中$\theta_{t+1,i}$表示时间$t+1$时，参数向量$\theta_{t+1}$的第$i$个参数，$g_{k,i}$表示$k$时刻梯度向量$g_k$的第$i$个维度。另外，分母中求和的形式实现了退火的过程，这是很多优化问题中常用的策略，意味着随着时间的推移，学习速率$\frac{\eta}{\sqrt{\sum_{k=0}^tg_{k,i}^2+\epsilon}}$会越来越小，从而保证算法最终收敛。
6.2.3 Adam
Adam把惯性保持和环境感知这两个有点集于一身。一方面，Adam记录了梯度的一阶矩(First Moment)，即过往梯度和当前梯度的平均，这体现了惯性的保持；另一方面，Adam还记录了梯度的二阶矩(Second Moment)，即过往梯度平方和当前梯度平方的平均，类似AdaGrad，体现了环境感知能力，为不同参数产生不同的学习速率。一阶矩和二阶矩采用了滑动窗口内求平均值的思想进行融合。即当前梯度和近一段时间内梯度的平均值，时间久远的梯度对当前平均值的贡献呈指数衰减。具体来说，一阶矩和二阶矩采用指数衰退平均(Exponential decay average)技术，计算公式为 $\begin{align}
&m_t=\beta_1m_{t-1}+(1-\beta_1)g_t \
&v_t=\beta_2v_{i-1}+(1-\beta_2)g_t^2
\end{align}$
其中$\beta_1$，$\beta_2$是衰减系数，$m_t$是一阶矩，$v_t$是二阶矩。一阶矩相当于估计$\mathbb{E}[g_t]$：由于当前梯度$g_t$是随机采样得到的估计结果，因此更关注它在统计意义上的期望；二阶矩相当于估计$\mathbb{E}[g_t^2]$，这点和AdaGrad不同，不是$g_t^2$从开始到现在的加和，而是他的期望。它的物理意义是，
- $\lVert m_t\rVert$大且$v_t$大，当前梯度大而且稳定，表明遇到一个大坡，前进方向明确；s
- $\lVert m_t\rVert$趋于零且$v_t$大，当前梯度不稳定，表明可能遇到一个峡谷，容易引起反弹震荡；
- $\lVert m_t\rVert$趋于零且$v_t$趋于零，当前梯度趋于零，表明可能遇到局部最低点，也有可能遇到坡度极缓的平底，要注意不要陷入平原(plataeu)；
- $\lVert m_t\rVert$大且$v_t$趋于零，这种情况不可能出现。
具体来说Adam的更新公式为 $\theta_{t+1}=\theta_t-\frac{\eta\cdot\hat{m}_t}{\sqrt{\hat{v}_t+\epsilon}}$ ，其中，$\hat{m}_t=\frac{m_t}{1-\beta_1^t}$，$\hat{v}_t=\frac{v_t}{1-\beta_2^t}$。
7. L1正则化与稀疏性
首先一个很基础的问题——为什么希望模型参数具有稀疏性？
所谓稀疏性，实际上就是指模型的参数很多为0。相当于对模型进行了一次参数选择，只留下比较重要的特征。这样可以提高模型的泛化能力，降低过拟合的可能。
7.1 为什么L1正则化会可以是参数具有稀疏性？
7.1.1 从解空间形状的角度解释
现在比较直观的解释是从解空间出发。
在二维的情况下，黄色的区域分别是L2和L1正则化约束下的解空间，绿色的等高线是凸优化问题中目标函数的等高线。如下图所示，多边形解空间更容易在尖角处碰撞出稀疏解。
但是为什么正则项会带来解空间约束？为什么L1和L2对应的形状不一样？深入探究的话，这些问题还需要进一步去解答。
可以通过KKT条件给出一种解释。
实际上，“带正则项”和“带约束条件“是等价的。为了约束$w$的可能的取值空间防止过拟合，需要给这个最优化问题加上一个约束，就是$w$的L2范数的平方不大于$m$：
$\begin{array}{lcl}min\sum_{i=1}^N(y_i-w^Tx_i)^2\quad &(7.1)\\ s.t.\quad \vert w\vert^2\le m \quad &(7.2)\end{array}$
解决带约束条件的凸优化问题，写出拉格朗日函数 $\sum_{i=1}^N(y_i-w^Tx_i)^2+\lambda(|w|^2-m)$
如果$w^\ast$和$\lambda^\ast$分别是原问题和对偶问题的最优解，则根据KKT条件，他们应该满
$\begin{array}{lcl}\nabla_w\big(\sum_{i=1}^N(y_i-w^{\ast T}x_i)^2+\lambda^\ast(|w^\ast|^2-m1)\big)=0\quad &(7.3)\\ \lambda^\ast\ge0\quad &(7.4)\end{array}$
其中(7.2)是$w^\ast$作为带L2正则项的优化问题最优解的条件，而$\lambda^\ast$为L2正则项前面的正则参数。所以式(7.2)中L2正则化恰好为参数$w$定义了一个圆形的解空间，而L1正则化则定义了一个菱形的解空间。
7.1.2 从函数叠加的角度解释
另一种非常直观的方式是从函数叠加的角度来解释。这里仅考虑一维的情况。
下图中，令原始目标函数$L(w)$用棕线表示，则最小值点可以用蓝点表示，对应的$w^\ast$不为0。
首先考虑加上L2正则项，则目标函数变为 $L(w)+C\lVert w\rVert^2$，其函数在图中用黄线表示，其中最小值点用黄点表示，对应$w^\ast$不为0，但是绝对值变小了。
然后再考虑加上L1正则项，则目标函数变为 $L(w)+C\lvert w\rvert$，对应函数在图中用绿线表示，其中最小值点用红点表示，对应的$w^\ast$为0。
产生上述现象的原因是，对加入L1正则项的目标函数求导，正则项部分在原点左边的导数为$-C$，原点右边的导数为$C$。因此如果原始目标函数的导数的绝对值小于$C$，那么导数在原点处为0，w取得最小值。但是对于L2正则化，因为正则项在原点左右的导数为$2Cw$，且在原点为0，所以如果原始目标函数在原点不为0，带L2正则项的目标函数就很难在原点取得导数为0。所以L2正则项主要的作用是减少$w$的绝对值大小。
在一些在线梯度下降算法里也会用截断梯度法来产生稀疏性。
7.1.3 从贝叶斯的角度解释
从贝叶斯角度来解释，简单来说，L1正则化相当于给$w$引入拉普拉斯先验，L2正则化相当于给$w$引入高斯先验，二拉普拉斯先验使得参数为0的概率更大。
如上图所示，高斯分布在取得极值点处是平滑的，即高斯分布认为$w$在极值点取不同值的概率是相近的，所以L2正则化会偏向于减小$w$；而拉普拉斯分布的极值点是一个尖峰，即参数$w$更偏向于在极值处取0。
