现在做在线学习和CTR常常会用到逻辑回归，而传统的批量算法无法有效地处理超大规模的数据集和在线数据流，google在2010年～2013年从理论研究到实际工程化实现的FTRL算法，在处理诸如逻辑回归之类的带非光滑正则化项（例如L1范数，做模型复杂度控制和稀疏化）的凸优化问题上性能非常出色。
参考：
https://zhuanlan.zhihu.com/p/32903540
https://zhuanlan.zhihu.com/p/61724627
Online Learning的学习算法主要有两类：在线凸优化和在线Bayesian。在线凸优化方法有很多，像FOBOS算法、RDA、FTRL等；在线Bayesian 方面有比如AdPredictor 算法、基于内容的在线矩阵分解算法等。本文是Google在2013年KDD上发表的论文，描述了FTRL-Proximal算法的工程化的实践。国内外公司纷纷进行了各自的实验，比如亚马逊，Yahoo，阿里，百度等也应用了该算法在其搜索领域、推荐领域等的产品中，且都取得了明显提升。另外在各大CTR竞赛，Kaggle比赛上时常能看到它的身影。
先验知识
逻辑回归
模型表达式
$y=\sigma(f(\boldsymbol{x}))=\sigma\left(\boldsymbol{w}^{T} \boldsymbol{x}\right)=\frac{1}{1+e^{-\boldsymbol{w}^{T} \boldsymbol{x}}}$
假设我们已经训练好了一组权值$w^T$ 。只要把我们需要预测的 x 代入到上面的方程，输出的y值就是这个标签为A的概率，我们就能够判断输入数据是属于哪个类别。
模型损失函数
我们把单个样本看做一个事件，那么这个事件发生的概率就是：$P(y \mid \boldsymbol{x})=\left\{\begin{array}{r}
p, y=1 \\
1-p, y=0
\end{array}\right.$
这个函数不方便计算，它等价于：$P\left(y_{i} \mid \boldsymbol{x}_{i}\right)=p^{y_{i}}(1-p)^{1-y_{i}}$。如果我们采集到了一组数据一共N个，即$\left\{\left(\boldsymbol{x}_{1}, y_{1}\right),\left(\boldsymbol{x}_{2}, y_{2}\right),\left(\boldsymbol{x}_{3}, y_{3}\right) \ldots\left(\boldsymbol{x}_{N}, y_{N}\right)\right\}$，这个合事件发生的总概率为：$\begin{aligned}
P_{\text {总 }} &=P\left(y_{1} \mid \boldsymbol{x}_{1}\right) P\left(y_{2} \mid \boldsymbol{x}_{2}\right) P\left(y_{3} \mid \boldsymbol{x}_{3}\right) \ldots P\left(y_{N} \mid \boldsymbol{x}_{N}\right) \\
&=\prod_{n=1}^{N} p^{y_{n}}(1-p)^{1-y_{n}}
\end{aligned}$
我们通过两边取对数来把连乘变成连加的形式：
$\begin{aligned}
F(\boldsymbol{w})=\ln \left(P_{\text {总 }}\right) &=\ln \left(\prod_{n=1}^{N} p^{y_{n}}(1-p)^{1-y_{n}}\right) \\
&=\sum_{n=1}^{N} \ln \left(p^{y_{n}}(1-p)^{1-y_{n}}\right) \\
&=\sum_{n=1}^{N}\left(y_{n} \ln (p)+\left(1-y_{n}\right) \ln (1-p)\right)
\end{aligned}$
其中$p=\frac{1}{1+e^{-\boldsymbol{w}^{T} \boldsymbol{x}}}$。这个函数又叫做它的损失函数。损失函数可以理解成衡量我们当前的模型的输出结果，跟实际的输出结果之间的差距的一种函数。这里的损失函数的值等于事件发生的总概率，我们希望它越大越好。但是跟损失的含义有点儿违背，因此也可以在前面取个负号。
最大似然估计MLE
我们知道损失函数 F(w) 是正比于总概率 P_总 的，而 F(w) 又只有一个变量 w 。也就是说，通过改变 w 的值，就能得到不同的总概率值 P_总 。那么当我们选取的某个 $w^$ 刚好使得总概率 P_总 取得最大值的时候。我们就认为这个 $w^$ 就是我们要求得的 w 的值，这就是最大似然估计的思想，即：$\boldsymbol{w}^{*}=\arg \max _{w} F(\boldsymbol{w})=-\arg \min _{w} F(\boldsymbol{w})$
下面推导一下F(w)的梯度：
计算p的导数：
$\begin{aligned}
p^{\prime}=f^{\prime}(\boldsymbol{w}) &=\left(\frac{1}{1+e^{-\boldsymbol{w}^{T} \boldsymbol{x}}}\right)^{\prime} \\
&=-\frac{1}{\left(1+e^{-\boldsymbol{w}^{T} \boldsymbol{x}}\right)^{2}} \cdot\left(1+e^{-\boldsymbol{w}^{T} \boldsymbol{x}}\right)^{\prime} \\
&=-\frac{1}{\left(1+e^{-\boldsymbol{w}^{T} \boldsymbol{x}}\right)^{2}} \cdot e^{-\boldsymbol{w}^{T} \boldsymbol{x}} \cdot\left(-\boldsymbol{w}^{T} \boldsymbol{x}\right)^{\prime} \\
&=-\frac{1}{\left(1+e^{-\boldsymbol{w}^{T} \boldsymbol{x}}\right)^{2}} \cdot e^{-\boldsymbol{w}^{T} \boldsymbol{x}} \cdot(-\boldsymbol{x}) \\
&=\frac{e^{-\boldsymbol{w}^{T} \boldsymbol{x}}}{\left(1+e^{-\boldsymbol{w}^{T} \boldsymbol{x}}\right)^{2}} \cdot \boldsymbol{x} \\
&=\frac{1}{1+e^{-\boldsymbol{w}^{T} \boldsymbol{x}}} \cdot \frac{e^{-\boldsymbol{w}^{T} \boldsymbol{x}}}{1+e^{-\boldsymbol{w}^{T} \boldsymbol{x}}} \cdot \boldsymbol{x} \\
&=p(1-p) \boldsymbol{x}
\end{aligned}$计算F(w)的导数：
$\begin{aligned}
\nabla F(\boldsymbol{w}) &=\nabla\left(\sum_{n=1}^{N}\left(y_{n} \ln (p)+\left(1-y_{n}\right) \ln (1-p)\right)\right) \\
&=\sum\left(y_{n} \ln ^{\prime}(p)+\left(1-y_{n}\right) l n^{\prime}(1-p)\right) \\
&=\sum\left(\left(y_{n} \frac{1}{p} p^{\prime}\right)+\left(1-y_{n}\right) \frac{1}{1-p}(1-p)^{\prime}\right) \\
&=\sum_{N}\left(y_{n}(1-p) \boldsymbol{x}_{n}-\left(1-y_{n}\right) p \boldsymbol{x}_{n}\right) \\
&=\sum_{n=1}^{N}\left(y_{n}-p\right) \boldsymbol{x}_{n}
\end{aligned}$
$\nabla F(\boldsymbol{w})=\sum_{n=1}^{N}\left(y_{n}-\frac{1}{1+e^{-\boldsymbol{w}^{T} \boldsymbol{x}_{n}}}\right) \boldsymbol{x}_{n}$
参数求解
梯度下降法（GD）
随便初始化一个$w_0$，然后给定一个步长$\eta$，通过不断地修改$\boldsymbol{w}_{t+1}<-\boldsymbol{w}_{t}$，从而最后靠近到达取得最大值的点，即不断进行下面的迭代过程，直到达到指定次数，或者梯度等于0为止：$\boldsymbol{w}_{t+1}=\boldsymbol{w}_{t}+\eta \nabla F \quad(\boldsymbol{w})$
随机梯度下降法（SGD）
如果我们能够在每次更新过程中，加入一点点噪声扰动，可能会更加快速地逼近最优值。在SGD中，我们不直接使用$\nabla F(w)$，而是采用另一个输出为随机变量的替代函数$G(w)$：$\boldsymbol{w}_{t+1}=\boldsymbol{w}_{t}+\eta G(\boldsymbol{w})$。
当然，这个替代函数G(w)需要满足它的期望值等于$\nabla F(\boldsymbol{w})$，相当于这个函数围绕着$\nabla F(\boldsymbol{w})$的输出值随机波动。
具体地，在SGD中，我们每次只要均匀地、随机选取其中一个样本(x_i, y_i)，即把它的值乘以N，就相当于获得了梯度的无偏估计值，即$E(G(\boldsymbol{w}))=\nabla F(\boldsymbol{w})$，因此SGD的更新公式为：$\boldsymbol{w}_{t+1}=\boldsymbol{w}_{t}+\eta N\left(y_{n}-\frac{1}{1+e^{-\boldsymbol{w}^{T} \boldsymbol{x}_{n}}}\right) \boldsymbol{x}_{n}$
这样我们前面的求和就没有了，同时$\eta N$都是常数，N的值刚好可以并入$\eta$当中,因此SGD的迭代更新公式为：$\boldsymbol{w}_{t+1}=\boldsymbol{w}_{t}+\eta\left(y_{n}-\frac{1}{1+e^{-\boldsymbol{w}^{T} \boldsymbol{x}_{n}}}\right) \boldsymbol{x}_{n}$
其中$(x_i, y_i)$是对所有样本随机抽样的一个结果。
L1、L2正则化的直观解释
参考：
https://blog.csdn.net/red_stone1/article/details/80755144
https://www.zhihu.com/question/37096933
https://blog.csdn.net/xpy870663266/article/details/104573194
https://www.cnblogs.com/heguanyou/p/7688344.html
稀疏解指的是参数矩阵的大部分参数值为0，这样得到的模型相对较小。稀疏的特征会大大减少predict时的内存和复杂度。先说几个结论：
- L1正则化会产生更稀疏的解，因此基于L1正则化的学习方法相当于嵌入式的特征选择方法.
- L2正则化计算更加方便，只需要计算向量内积，L1范数的计算效率特别是遇到非稀疏向量时非常低
- L1正则化相当于给权重设置拉普拉斯先验，L2正则化相当于给权重设置高斯先验
- L1倾向于学得稀疏的权重矩阵，L2倾向于学得更小更分散的权重
L2 正则化直观解释
L2 正则化公式：$L=E_{i n}+\lambda \sum_{j} w_{j}^{2}$，这个式子是怎么得来的呢？我们知道，正则化的目的是限制参数过多或者过大，避免模型更加复杂。例如，使用多项式模型，如果使用 10 阶多项式，模型可能过于复杂，容易发生过拟合。所以，为了防止过拟合，我们可以将其高阶部分的权重 w 限制为 0，这样，就相当于从高阶的形式转换为低阶。
为了达到这一目的，最直观的方法就是限制 w 的个数，但是这类条件属于 NP-hard 问题，求解非常困难。所以，一般的做法是寻找更宽松的限定条件：$\sum_{j} w_{j}^{2} \leq C$。因此，我们的目标转换为：最小化训练样本误差 Ein，但是要遵循 w 平方和小于 C 的条件，用一张图表示如下：
蓝色椭圆区域是最小化 Ein 区域，红色圆圈是 w 的限定条件区域。在没有限定条件的情况下，一般使用梯度下降算法，在蓝色椭圆区域内会一直沿着 w 梯度的反方向前进，直到找到全局最优值 wlin。例如空间中有一点 w（图中紫色点），此时 w 会沿着 -∇Ein 的方向移动，如图中蓝色箭头所示。但是，由于存在限定条件，w 不能离开红色圆形区域，最多只能位于圆上边缘位置，沿着切线方向。w 的方向如图中红色箭头所示。
那么w 最终会在什么位置取得最优解呢？w 是沿着圆的切线方向运动，如上图绿色箭头所示。运动方向与 w 的方向（红色箭头方向）垂直。运动过程中，根据向量知识，只要 -∇Ein 与运行方向有夹角，不垂直，则表明 -∇Ein 仍会在 w 切线方向上产生分量，那么 w 就会继续运动，寻找下一步最优解。只有当 -∇Ein 与 w 的切线方向垂直时，-∇Ein在 w 的切线方向才没有分量，这时候 w 才会停止更新，到达最接近 wlin 的位置，且同时满足限定条件。
-∇Ein 与 w 的切线方向垂直，即 -∇Ein 与 w 的方向平行。如上图所示，蓝色箭头和红色箭头互相平行。这样，根据平行关系得到：$\nabla E_{i n}+\lambda w=0$，这样，我们就把优化目标和限定条件整合在一个式子中了。也就是说只要在优化 Ein 的过程中满足上式，就能实现正则化目标。
根据最优化算法的思想：梯度为 0 的时候，函数取得最优值。已知 ∇Ein 是 Ein 的梯度，观察上式，λw 是否也能看成是某个表达式的梯度呢？当然可以！λw 可以看成是 1/2λw*w 的梯度：$\frac{\partial}{\partial w}\left(\frac{1}{2} \lambda w^{2}\right)=\lambda w$。这样，我们根据平行关系求得的公式，构造一个新的损失函数：$E_{a u g}=E_{i n}+\frac{\lambda}{2} w^{2}$
L1 正则化直观解释
L1 正则化公式也很简单，直接在原来的损失函数基础上加上权重参数的绝对值：$L=E_{i n}+\lambda \sum_{j}\left|w_{j}\right|$，我仍然用一张图来说明如何在 L1 正则化下，对 Ein 进行最小化的优化：
L1 与 L2 解的稀疏性
直观视角
以二维情况讨论，上图左边是 L2 正则化，右边是 L1 正则化。从另一个方面来看，满足正则化条件，实际上是求解蓝色区域与黄色区域的交点，即同时满足限定条件和 Ein 最小化。对于 L2 来说，限定区域是圆，这样，得到的解 w1 或 w2 为 0 的概率很小，很大概率是非零的。
对于 L1 来说，限定区域是正方形，方形与蓝色区域相交的交点是顶点的概率很大，这从视觉和常识上来看是很容易理解的。也就是说，方形的凸点会更接近 Ein 最优解对应的 w 位置，而凸点处必有 w1 或 w2 为 0。这样，得到的解 w1 或 w2 为零的概率就很大了。所以，L1 正则化的解具有稀疏性。
扩展到高维，同样的道理，L2 的限定区域是平滑的，与中心点等距；而 L1 的限定区域是包含凸点的，尖锐的。这些凸点更接近 Ein 的最优解位置，而在这些凸点上，很多 wj 为 0。
例子视角
上面说的稍微有点抽象，我们以下图中一个具体的Loss函数为例：
则最优的 x 在绿点处，x 非零。现在施加 L2 regularization，新的损失函数$\left(L+C x^{2}\right)$如图中蓝线所示：
最优的 x 在黄点处，x 的绝对值减小了，但依然非零。而如果施加 L1 regularization，则新的损失函数$(L+C|x|)$如图中粉线所示：
最优的 x 就变成了 0。这里利用的就是绝对值函数的尖峰。两种 regularization 能不能把最优的 x 变成 0，取决于原先的费用函数在 0 点处的导数。如果本来导数不为 0，那么施加 L2 regularization 后导数依然不为 0，最优的 x 也不会变成 0。而施加 L1 regularization 时，只要保证施加L1后x=0处左右两边导数异号，且施加L1后x=0处左右两边导数分别是f’(0)-C和f’(0)+C，即C>|f’(0)|，那么就会形成极小值点。
概率视角
为f(w)加入正则化项，相当于为f(w)的参数加入先验，即要求w满足某以分布。L1正则化项相当于为w加入拉普拉斯先验，L2正则化项相当于加入高斯分布的先验，两种分布的概率函数如下图：
很明显可以观察出，在紫色区域中，P_G(w) < P_L(w)，说明Gauss分布中，值大的w更少。即L2与L1相比，值大的w更少，因此L2比L1更平滑。
在红色区域中，P_L(w) < P_G(w)，结合图来看，Gauss分布中值很小的w和值为0的w概率很近。而Laplace分布中，值很小的w小于值为0的w。这说明laplace要w更多为0，而Gauss分布要w很小但不一定为0。因此L1比L2更稀疏。
梯度视角
在L2正则下对参数求导为：
$\begin{aligned}
C=& C_{0}+\frac{\lambda}{2 n} \sum_{w} w^{2} \\
\frac{\partial C}{\partial w} &=\frac{\partial C_{0}}{\partial w}+\frac{\lambda}{n} w \\
\frac{\partial C}{\partial b} &=\frac{\partial C_{0}}{\partial b} \\
w & \rightarrow w-\eta \frac{\partial C_{0}}{\partial w}-\frac{\eta \lambda}{n} \\
&=\left(1-\frac{\eta \lambda}{n}\right) w-\eta \frac{\partial C_{0}}{\partial w} .
\end{aligned}$
由上式可以看出加入L2正则项后梯度更新的变化，即在每步执行通常的梯度更新之前先缩放权重向量。在L1正则下对参数求导为：
- $C=C_{0}+\frac{\lambda}{n} \sum_{w}|w|$
- $\frac{\partial C}{\partial w}=\frac{\partial C_{0}}{\partial w}+\frac{\lambda}{n} \operatorname{sgn}(w)$
- $w \rightarrow w^{\prime}=w-\frac{\eta \lambda}{n} \operatorname{sgn}(w)-\eta \frac{\partial C_{0}}{\partial w}$
可以看到L1正则化的效果与L2正则化的效果很不一样，不再是线性地缩放每个$w_i$，而是减去了一项与权重w i同号的常数，因此当$w_i$ >0时，更新权重使得其减小，当$w_i< 0$时，更新权重使得其增大。从这点来说也可以看出L1正则化更有可能使得权重为0，而L2正则化虽也使得权重减小，但缩放操作使得其仍保持同号。
可以想象用梯度下降时，当w小于1时，L2正则项的惩罚效果越来越小，L1正则项惩罚效果依然很大，L1可以惩罚到0，而L2很难。
正则化参数 λ
正则化是结构风险最小化的一种策略实现，能够有效降低过拟合。损失函数实际上包含了两个方面：一个是训练样本误差。一个是正则化项。其中，参数 λ 起到了权衡的作用。
以 L2 为例，若 λ 很小，对应上文中的 C 值就很大。这时候，圆形区域很大，能够让 w 更接近 Ein 最优解的位置。若 λ 近似为 0，相当于圆形区域覆盖了最优解位置，这时候，正则化失效，容易造成过拟合。相反，若 λ 很大，对应上文中的 C 值就很小。这时候，圆形区域很小，w 离 Ein 最优解的位置较远。w 被限制在一个很小的区域内变化，w 普遍较小且接近 0，起到了正则化的效果。但是，λ 过大容易造成欠拟合。欠拟合和过拟合是两种对立的状态。
【注】：有关Laplace先验与L1正则化和Gauss先验与L2正则化可参考https://www.cnblogs.com/heguanyou/p/7688344.html
在线学习和稀疏性简介
稀疏性在机器学习中是很看重的事情，尤其我们做工程应用，稀疏的特征会大大减少predict时的内存和复杂度。这一点其实很容易理解，说白了，即便加入L1范数，因为是浮点运算，训练出的w向量也很难出现绝对的零。到这里，大家可能会想说，那还不容易，当计算出的w对应维度的值很小时，我们就强制置为零不就稀疏了么。对的，其实不少人就是这么做的，后面的Truncated Gradient和FOBOS都是类似思想的应用。GD求解的模型虽然精度相对较高，但具有训练太费时、不易得到稀疏解和对不可微点迭代效果欠佳等缺点；SGD则存在模型解的精度低、收敛速度慢和很难得到稀疏解的缺点。虽然学术界提出了很多加快收敛或者提高模型精度的方法（如添加momentum项、添加nesterov项、Adagrad算法、Adadelta算法、GSA算法等），但这些方法在提高模型解的稀疏性方面效果有限，而FTRL在这方面则更加有效。
Online Gradient Descent(OGD)
在线梯度下降法类似于随机梯度下降SGD，一次只处理一个样本，但不同的是不要求序列的样本满足独立同分布的性质，其更新公式为：$w_{t+1}=w_{t}-\eta_{t} \mathbf{g}_{t}$，其中$\eta_{t}=\frac{1}{\sqrt{t}}$是全局非增的。OGD方法在许多问题上均有不错的表现，但是获得的解不够稀疏。在OGD目标函数中添加L1正则化项，可以带来一定的稀疏性，使得模型的参数接近零。但由于参数的更新是两个浮点数进行相加，很难有机会产生精确零解。
Truncated Gradient
Truncated Gradient具体算法如下：
从上面的公式可以看出，$\lambda, \theta$决定了 W 的稀疏程度，如果两者都很大，那么稀疏性就会越强。特别的，当$\lambda=\theta$时，此时只需要控制一个参数就可以控制稀疏性。
FOBOS
在该算法中，权重的更新分成两个步骤：
- $W^{(t+0.5)}=W^{(t)}-\eta^{(t)} G^{(t)}$
- $W^{(t+1)}=\operatorname{argmin}_{W}\left\{\frac{1}{2}\left|W-W^{(t+0.5)}\right|_{2}^{2}+\eta^{(t+0.5)} \Psi(W)\right\}$
第一个步骤实际上是一个标准的梯度下降（SGD），第二个步骤是对第一个步骤的结果进行局部调整。写成一个公式那就是：$W^{(t+1)}=\operatorname{argmin}_{W}\left\{\frac{1}{2}\left|W-W^{(t)}+\eta^{(t)} G^{(t)}\right|_{2}^{2}+\eta^{(t+0.5)} \Psi(W)\right\}$
假设$F(W)=\frac{1}{2}\left|W-W^{(t)}+\eta^{(t)} G^{(t)}\right|_{2}^{2}+\eta^{(t+0.5)} \Psi(W)$，如果$W^{t+1}$存在一个最优解，那么可以推断0向量一定属于F(W)的次梯度集合：$0 \in \partial F(W)=W-W^{(t)}+\eta^{(t)} G^{(t)}+\eta^{(t+0.5)} \partial \Psi(W)$。因为$W^{(t+1)}=\operatorname{argmin}_{W} F(W)$，那么可以得到权重更新的另外一种形式：$W^{(t+1)}=W^{(t)}-\eta^{(t)} G^{(t)}-\eta^{(t+0.5)} \partial \Psi\left(W^{(t+1)}\right)$
从上面的公式可以看出，更新后的$W^{t+1}$不仅和$W^{(t)}$有关，还和自己$\Psi\left(W^{(t+1)}\right)$有关。这也许就是”前向后向切分”这个名称的由来。
L1-FOBOS
FOBOS的L1正则化形式为：
$\begin{array}{l}
W^{(t+1)}=\operatorname{argmin}_{W}\left\{\frac{1}{2}\left|W-W^{(t+0.5)}\right|_{2}^{2}+\eta^{(t+0.5)} \Psi(W)\right\} \\
=\operatorname{argmin}_{W} \sum_{i=1}^{N}\left(\frac{1}{2}\left(w_{i}-v_{i}\right)^{2}+\tilde{\lambda}\left|w_{i}\right|\right)
\end{array}$
其中$\Psi(W)$是L1范数，中间向量是$W^{(t+0.5)}=\left(v_{1}, \ldots, v_{N}\right) \in \mathbb{R}^{N}$，并且参数$\tilde{\lambda}=\eta^{(t+0.5)} \lambda$。可以对特征权重W的每一个维度进行单独求解：$w_{i}^{(t+1)}=\operatorname{argmin}_{w_{i}}\left(\frac{1}{2}\left(w_{i}-v_{i}\right)^{2}+\tilde{\lambda}\left|w_{i}\right|\right) \text { for all } 1 \leq i \leq N$
通过数学计算可证明，在L1正则化下，FOBOS 的特征权重的各个维度的更新公式是：
$\begin{array}{l}
w_{i}^{(t+1)}=\operatorname{sgn}\left(v_{i}\right) \max \left(0,\left|v_{i}\right|-\tilde{\lambda}\right) \\
=\operatorname{sgn}\left(w_{i}^{(t)}-\eta^{(t)} g_{i}^{(t)}\right) \max \left\{0,\left|w_{i}^{(t)}-\eta^{(t)} g_{i}^{(t)}\right|-\eta^{(t+0.5)} \lambda\right\} \text { for all } 1 \leq i \leq N
\end{array}$
其中$g_{i}^{(t)}$是梯度$G^{(t)}$)在第i个维度的分量。
从公式上看，如果$\left|w_{i}^{(t)}-\eta^{(t)} g_{i}^{(t)}\right| \leq \eta^{(t+0.5)} \lambda$，则有$w_{i}^{(t+1)}=0$。意思就是如果这次训练产生梯度的变化不足以令权重值发生足够大的变化时，就认为在这次训练中该维度不够重要，应该强制其权重是0。如果$\left|w_{i}^{(t)}-\eta^{(t)} g_{i}^{(t)}\right| \geq \eta^{(t+0.5)} \lambda$，则有$w_{i}^{(t+1)}=\left(w_{i}^{(t)}-\eta^{(t)} g_{i}^{(t)}\right)-\eta^{(t+0.5)} \lambda \cdot \operatorname{sgn}\left(w_{i}^{(t)}-\eta^{(t)} g_{i}^{(t)}\right)$。
那么L1-FOBOS和Truncated Gradient有什么关系呢？令$\theta=+\infty, k=1, \lambda_{T G}^{(t)}=\eta^{(t+0.5)} \lambda$，可以得到L1-FOBOS和Truncated Gradient完全一致，换句话说L1-FOBOS是Truncated Gradient在一些特定条件下的形式。
Regularized Dual Averaging (RDA)
RDA又叫做正则对偶平均算法，特征权重的更新策略是：$W^{(t+1)}=\operatorname{argmin}_{W}\left\{\frac{1}{t} \sum_{r=1}^{t} G^{(r)} \cdot W+\Psi(W)+\frac{\beta^{(t)}}{t} h(W)\right\}$
其中$G^{(t)} \cdot W$指的是向量$G^{(t)}$和W的内积，$\frac{1}{t} \sum_{r=1}^{t} G^{(r)} \cdot W$包括了之前所有梯度的平均值，$\Psi(W)$是正则项，h(W)是一个严格凸函数，$\left\{\beta^{(t)} \mid t \geq 1\right\}$是一个非负递增序列。
L1-RDA
令$\Psi(W)=|W|_{1}$，$h(W)=\frac{1}{2}|W|_{2}^{2}$（满足凸函数要求），$\beta^{(t)}=\gamma \sqrt{t} \text { with } \gamma>0$（满足非负递增序列要求），则RDA的L1正则化形式为：
$W^{(t+1)}=\operatorname{argmin}_{W}\left\{\frac{1}{t} \sum_{r=1}^{t} G^{(r)} \cdot W+\lambda|W|_{1}+\frac{\gamma}{2 \sqrt{t}}|W|_{2}^{2}\right\}$
此时可以采取同样的策略，把各个维度拆分成N个独立的问题来处理，即：
$w_{i+1}^{(t+1)}=\operatorname{argmin}_{w_{i}}\left\{\bar{g}_{i}^{(t)}+\lambda\left|w_{i}\right|+\frac{\gamma}{2 \sqrt{t}} w_{i}^{2}\right\}$
其中$\lambda>0, \gamma>0, \bar{g}_{i}^{(t)}=\frac{1}{t} \sum_{r=1}^{t} g_{i}^{(r)}$。
通过数学推导可以证明：
- $w_{i}^{(t+1)}=0, \text { if }\left|\bar{g}_{i}^{(t)}\right|<\lambda$
- $w_{i}^{(t+1)}=-\frac{\sqrt{t}}{\gamma}\left(\bar{g}_{i}^{(t)}-\lambda \cdot \operatorname{sgn}\left(\bar{g}_{i}^{(t)}\right)\right), \text { otherwise }$
意思就是：当某个维度的累加平均值$\left|\bar{g}_{i}^{(t)}\right|$小于$\lambda$时，该维度的权重被设置成零，此时就可以产生特征权重的稀疏性。
RDA有如下特点：
- 非梯度下降类方法，属于更加通用的一个primal-dual algorithmic schema的一个应用
- 克服了SGD类方法所欠缺的exploiting problem structure，especially for problems with explicit regularization。
- 能够更好地在精度和稀疏性之间做trade-off
L1-RDA和L1-FOBOS
L1-RDA和L1-FOBOS有什么区别和联系吗？从截断方式来看，在 RDA 的算法中，只要梯度的累加平均值小于参数$\lambda$就直接进行截断，说明 RDA 更容易产生稀疏性；同时，RDA 中截断的条件是考虑梯度的累加平均值，可以避免因为某些维度训练不足而导致截断的问题，这一点与 TG，FOBOS 不一样。通过调节参数$\lambda$可以在精度和稀疏性上进行权衡。
L1-RDA和L1-FOBOS的形式上的统一。在 L1-FOBOS中，假设$\eta^{(t+0.5)}=\eta^{(t)}=\Theta\left(t^{-0.5}\right)$，可以得到$W^{(t+1)}=\operatorname{argmin}_{W}\left\{\frac{1}{2}\left|W-W^{(t)}+\eta^{(t)} G^{(t)}\right|_{2}^{2}+\eta^{(t)} \lambda|W|_{1}\right\}$。通过计算可以把L1-FOBOS写成：$W^{(t+1)}=\operatorname{argmin}_{W}\left\{G^{(t)} \cdot W+\lambda|W|_{1}+\frac{1}{2 \eta^{(t)}}\left|W-W^{(t)}\right|_{2}^{2}\right\}$。同理L1-RDA可以写成$W^{(t+1)}=\operatorname{argmin}_{W}\left\{G^{(1: t)} \cdot W+t \lambda|W|_{1}+\frac{1}{2 \eta^{(t)}}|W-0|_{2}^{2}\right\}$（其中$G^{(1: t)}=\sum_{s=1}^{t} G^{(s)}$）。
如果令$\sigma^{(s)}=\frac{1}{\eta^{(s)}}-\frac{1}{\eta^{(s-1)}}$、$\sigma^{(1: t)}=\frac{1}{\eta^{(t)}}$、$\sigma^{(1: t)}=\sum_{s=1}^{t} \sigma^{(s)}$，则L1-FOBOS可以写成$W^{(t+1)}=\operatorname{argmin}_{W}\left\{G^{(t)} \cdot W+\lambda|W|_{1}+\frac{1}{2} \sigma^{(1: t)}\left|W-W^{(t)}\right|_{2}^{2}\right\}$，L1-RDA可以写成$W^{(t+1)}=\operatorname{argmin}_{W}\left\{G^{(1: t)} \cdot W+t \lambda|W|_{1}+\frac{1}{2} \sigma^{(1: t)}|W-0|_{2}^{2}\right\}$。
可以看出：
- L1-FOBOS考虑了梯度和当前模的L1 正则项，L1-RDA考虑的是累积梯度。
- L1-FOBOS的第三项限制W不能够距离已经迭代过的$W^{(t)}$值太远，L1-RDA的第三项限制的是W不能够距离0太远。
Follow The (Proximally) Regularized Leader (FTRL-Proximal)
参考：
https://zhuanlan.zhihu.com/p/25315553
https://zhuanlan.zhihu.com/p/32903540
理论及实验均证明，L1-FOBOS这类基于梯度下降的算法有比较高的精度，但L1-RDA却能在损失一定精度的情况下产生更好的稀疏性。把这二者的优点结合成一个算法，这就是FTRL算法的来源,如下图：
FTRL 综合考虑了 FOBOS 和 RDA 对于梯度和正则项的优势和不足，其权重更新公式是：$W^{(t+1)}=\operatorname{argmin}_{W}\left\{G^{(1: t)} \cdot W+\lambda_{1}|W|_{1}+\frac{\lambda_{2}}{2}|W|_{2}^{2}+\frac{1}{2} \sum_{s=1}^{t} \sigma^{(s)}\left|W-W^{(s)}\right|_{2}^{2}\right\}$
上面的公式出现了L2范数，不过这一项的引入不会影响 FTRL 的稀疏性，只是使得求解结果更加“平滑”。通过数学计算并且放弃常数项可以得到上面的优化问题相当于求使得下面式子的最小的参数 W：$\left(G^{(1: t)}-\sum_{s=1}^{t} \sigma^{(s)} W^{(s)}\right) \cdot W+\lambda_{1}|W|_{1}+\frac{1}{2}\left(\lambda_{2}+\sum_{s=1}^{t} \sigma^{(s)}\right)|W|_{2}^{2}$
假设$Z^{(t)}=G^{(1: t)}-\sum_{s=1}^{l} \sigma^{(s)} W^{(s)}=Z^{(t-1)}+G^{(t)}-\sigma^{(t)} W^{(t)}$，则上式等价于$W^{(t+1)}=\operatorname{argmin}_{W}\left\{Z^{(t)} \cdot W+\lambda_{1}|W|_{1}+\frac{1}{2}\left(\lambda_{2}+\sum_{s=1}^{t} \sigma^{(s)}\right)|W|_{2}^{2}\right\}$，写成分量的形式就是：$\operatorname{argmin}_{w_{i}}\left\{z_{i}^{(t)} w_{i}+\lambda_{1}\left|w_{i}\right|+\frac{1}{2}\left(\lambda_{2}+\sum_{s=1}^{t} \sigma^{(s)}\right) w_{i}^{2}\right\}$
通过数学推导可以证明：
- $w_{i}^{(t+1)}=0, \text { if }\left|z_{i}^{(t)}\right|<\lambda_{1}$
- $w_{i}^{(t+1)}=-\left(\lambda_{2}+\sum_{s=1}^{t} \sigma^{(s)}\right)^{-1} \cdot\left(z_{i}^{(t)}-\lambda_{1} \cdot \operatorname{sgn}\left(z_{i}^{(t)}\right)\right), \text { otherwise }$
上式也可以说明，引入 L2 正则化并没有对 FTRL 的稀疏性产生影响。
Per-Coordinate Learning Rates
在 SGD 的算法里面使用的是一个全局的学习率$\eta^{(t)}=O\left(t^{-0.5}\right)$，意味着学习率是一个正数并且逐渐递减，对每一个维度都是一样的。而在 FTRL 算法里面，每个维度的学习率是不一样的。如果特征 A 比特征 B变化快，那么在维度 A 上面的学习率应该比维度 B 上面的学习率下降得更快。在 FTRL 中，维度i的学习率的定义如下：$\eta_{i}^{(t)}=\alpha /\left(\beta+\sqrt{\sum_{s=1}^{t}\left(g_{i}^{(s)}\right)^{2}}\right)$
按照之前的定义$\sigma^{(1: t)}=1 / \eta^{(t)}$，所以$\sum_{s=1}^{t} \sigma^{(s)}=\eta_{i}^{(t)}=\alpha /\left(\beta+\sqrt{\sum_{s=1}^{t}\left(g_{i}^{(s)}\right)^{2}}\right)$
因此，FTRL的整体算法如下：
所谓的per-coordinate，意思就是FTRL是对w每一维分开训练更新的，每一维使用的是不同的学习速率。与w所有特征维度使用统一的学习速率相比，这种方法考虑了训练样本本身在不同特征上分布的不均匀性，如果包含w某一个维度特征的训练样本很少，每一个样本都很珍贵，那么该特征维度对应的训练速率可以独自保持比较大的值，每来一个包含该特征的样本，就可以在该样本的梯度上前进一大步，而不需要与其他特征维度的前进步调强行保持一致。
Logistic Regression的FTRL形式
Logistic Regression 的 FTRL算法如下：
不同算法的统一描述
不同的方法在统一的描述形式下，区别点仅在第二项和第三项的处理方式：
- 第一项：梯度或累积梯度；
- 第二项：L1正则化项的处理；
- 第三项：这个累积加和限定了新的迭代结果x不要离已迭代过的解太远（也即FTRL-Proximal中proximal的含义），或者离0太远（central），这一项其实也是low regret的需求。
节省内存的其他方法
下面介绍一些google在FTRL工程实践上的一些tricks。
Feature Inclusion
海量的特征中(数亿)往往包含一些特征，它们的可能在海量的样例中只有一个或极少个非零值，这样的特征对于模型基本无效，所以在数据预处理时要尽量去掉这些无效特征。但是在线学习中，要重新读取和重写海量数据是非常耗费资源的，而且如果真的在预处理阶段就删除了这些特征，那么在学习的整个过程中就不再用到这些特征。
针对这个问题，有以下常见思路：1）使用L1正则项来实现训练中的稀疏性，这种方法不用刻意统计特征出现的频率，将频率低的删除，因为随着训练的进行，信息量较少的特征就会被自动删除，但是经过实验发现这个方法相比FTRL-Proximal会降低模型预测的准确率。2）制造hash冲突，但是作用也不大。
本文探索了一种新思路：probabilistic feature inclusion，意思就是当新特征首次出现的时候就概率性的包含在模型中，这相当于在第一次读入数据的时候就进行了预处理，不必再重新读写。这种读入数据时预处理的方式可以在线设置，本文测试了两种方式，下图为两种方式的效果比较，显然BLOOM的方法在RAM和AucLos之间有更好的平衡：
- Poisson Inclusion：假设特征是随机到达的，也就是说应该符合Poisson分布。以p的概率随机加入模型并进行训练。
- Bloom Filter Inclusion：使用布隆过滤器来进行计数，当特征出现次数超过n时才加入模型。但布隆过滤器并非精确的，而是false positive的，因为产生冲突时，可能会加入还没有满n次的特征来进行训练。实际中我们并不要求完全精确的筛选。在实验中，Bloom Filter表现出更高的RAM节省及更少的AucLoss，更推荐。
Encoding Values with Fewer Bits
参考：
定点数与浮点数：
https://blog.csdn.net/dong_zhihong/article/details/8255690
https://blog.csdn.net/niaolianjiulin/article/details/82764511
https://zhuanlan.zhihu.com/p/63897066
我们平常使用的都是32或64位浮点数编码，但是在实际场景中，系数值普遍落在(-2,+2)之间，普通的编码方式就显得非常奢侈了。这里提出了q2.13编码方式：
1 bit表达正负号，2 bits表达整数位，13 bits表达小数点右边的值，总计是16 bits。用更少的bit来表达就意味着信息损失。但是算法中要运用到对每一小步的累计值，不做处理就会对模型产生影响。具体解决策略为：
加入一个随机的regret项:R，R符合标准高斯分布，值介于0和1之间。对系数进行round，保证其整体离散误差为0：$w_{i, \text { rounded }}=2^{-13}\left\lfloor 2^{13} w_{i}+R\right\rfloor$。实验表明，使用了q2.13相较于64位浮点编码可以节省75%的参数存储RAM。
Training Many Similar Models
我们在做特征筛选或者设置优化时，基于已有的训练数据，控制变量，生成多个不同的变体来进行对比实验是一个不错的方法。但这样做比较耗费资源。之前有一个比较省资源的方法是基于一个先验模型，使用其他模型来预测残差。但是这样不能进行特征移除和特征替换的效果。这里Google提出了一个方法，由于每个模型变体的相似度是很高的，一部分特征是变体之间共享的，一部分特征是变体独有的。如果采用一张hashtable来存储所有变体的模型参数，就可以摊销各个模型都进行独立存储key的内存消耗。而且也同时会降低看网络带宽、CPU、磁盘的消耗。
A Single Value Structure
这个是指不同模型变体之间共享一个特征参数，而不是每一个模型变体存储单独的一个。用一个位字段来标识这个特征被哪些变体使用。
更新的方式如下：每一个模型变体去进行prediction并计算loss，然后对于每一个特征给出一个迭代后的参数值。对所有的更新后的参数值进行平均并进行更新存储。
Computing Learning Rates with Counts
前面提到学习率的计算per-coordinate learning rates，要对梯度的平方要进行累计。这里给出了一个梯度平方累计近似解决方案：
$\begin{aligned}
\sum g_{t, i}^{2} &=\sum_{\text {positive events }}\left(1-p_{t}\right)^{2}+\sum_{\text {negative events }} p_{t}^{2} \\
& \approx P\left(1-\frac{P}{N+P}\right)^{2}+N\left(\frac{P}{N+P}\right)^{2} \\
&=\frac{P N}{N+P}
\end{aligned}$
Subsampling Training Data
实际CTR问题中，点击率远远低于50%，一般在1%左右，那就意味着负样本数量远超正样本数量。做数据重采样不仅能平衡正负样本比例，还可以在不影响整体准确度的情况下降低训练集大小。具体的筛选规则描述如下：
- query对应至少一个ad点击的全采样；
- query不对应任何ad点击的按r比例采样；
然后对样本进行权重修正：
$\omega_{t}=\left\{\begin{array}{ll}
1 & \text { event } t \text { is in a clicked query } \\
\frac{1}{r} & \text { event } t \text { is in a query with no clicks. }
\end{array}\right.$