8.int-ml

下载地址:https://github.com/daiwk/collections/blob/master/pdfs/int-ml.pdf

本文参考自李航的《统计学习方法》、周志华的《机器学习》、Hulu的《百面机器学习》等。

概述

统计学习三要素

模型

监督学习中,模型是要学习的条件概率分布或决策函数

模型的假设空间

假设空间是所有可能的条件概率分布或决策函数

定义1

可以定义为决策函数的集合

F={fY=f(X)}\mathcal{F}=\{f|Y=f(X)\}
  • XXYY是定义在X\mathcal{X}Y\mathcal{Y}上的变量

  • F\mathcal{F}是一个参数向量决定的函数族

F={fY=fθ(X),θRn}\mathcal{F}=\{f|Y=f_\theta(X),\theta \in R^n\}

参数向量θ\theta取值于n维欧式空间RnR^n,称为参数空间

定义2

也可以定义为条件概率的集合

F={PP(YX)}\mathcal{F}=\{P|P(Y|X)\}
  • XXYY是定义在X\mathcal{X}Y\mathcal{Y}上的随机变量

  • F\mathcal{F}是一个参数向量决定的条件概率分布族

F={PPθ(YX),θRn}\mathcal{F}=\{P|P_\theta(Y|X),\theta \in R^n\}

策略

损失函数与风险函数

损失函数(loss function)或代价函数(cost function): 度量预测值f(X)f(X)与真实值YY的误差程度,记为L(Y,f(X))L(Y,f(X)),是个非负实值函数。损失函数越小,模型越好。

  • 0-1损失函数:

L(Y,f(X))={0Yf(X)1Y=f(X)L(Y,f(X))=\left\{\begin{matrix} 0 & Y\neq f(X) \\ 1 & Y = f(X) \end{matrix}\right.
  • 平方损失函数:

L(Y,f(X))=(Yf(X))2L(Y,f(X))=(Y-f(X))^2
  • 绝对损失函数:

L(Y,f(x))=Yf(X)L(Y,f(x))=|Y-f(X)|
  • 对数损失函数(logarithmic loss function)/对数似然损失函数(log-likelihood loss function):

L(Y,P(YX))=logP(YX)L(Y,P(Y|X))=-logP(Y|X)

风险函数(risk function)或期望损失(expected loss)XXYY服从联合分布P(X,Y)P(X,Y),理论上模型f(X)f(X)关于联合分布P(X,Y)P(X,Y)的平均意义下的损失:

Rexp(f)=EP[L(Y,f(X))]=X×YL(y,f(x))P(x,y)dxdyR_{exp}(f)=E_P[L(Y,f(X))]=\int _{\mathcal{X}\times \mathcal{Y}}L(y,f(x))P(x,y)dxdy

学习的目标:选择期望风险最小的模型。但联合分布P(X,Y)P(X,Y)是未知的,所以无法直接计算Rexp(f)R_{exp}(f)。所以监督学习是病态问题(ill-formed problem):一方面需要联合分布,另一方面联合分布是未知的。

给定训练集:

T={(x1,y1),...(xN,yN)}T=\{(x_1,y_1),...(x_N,y_N)\}

经验风险(expirical risk)/经验损失(expirical loss):模型f(X)f(X)关于训练集的平均损失

Remp(f)=1Ni=1NL(yi,f(xi))R_{emp}(f)=\frac{1}{N}\sum ^N_{i=1}{L(y_i,f(x_i))}

根据大数定律,当样本容量NN趋向无穷时,经验风险RempR_{emp}趋于期望风险Rexp(f)R_{exp}(f)

经验风险最小化与结构风险最小化

经验风险最小化(empirical risk minimization, ERM)经验风险最小的模型就是最优模型。所以需要求解的最优化问题是:

minfFRerm=minfF1NL(yi,f(xi))\min_{f\in \mathcal{F}}R_{erm}=\min_{f\in \mathcal{F}}\frac{1}{N}L(y_i,f(x_i))

当满足以下两个条件时,经验风险最小化就等价于极大似然估计(maximum likelihood estimation):

  • 模型是条件概率分布

  • 损失函数是对数损失函数

当样本量足够大时,ERM能有很好的效果,但样本量不够多时,为了防止过拟合,需要用下面的方法。

结构风险最小化(structual risk minimization, SRM):结构风险=经验风险+表示模型复杂度的正则化项(regularizer)或罚项(penalty term)。结构风险定义如下:

Rsrm(f)=1Ni=1NL(yi,f(xi))+λJ(f)R_{srm}(f)=\frac{1}{N}\sum^N_{i=1}L(y_i,f(x_i))+\lambda J(f)

其中,J(f)J(f)是模型的复杂度,模型越复杂,J(f)J(f)越大。λ0\lambda\ge 0是用于权衡经验风险和模型复杂度的系数。

当满足以下3个条件时,结构化风险最小化等价于)贝叶斯估计中的最大后验概率估计(maximum posterior probability estimation, MAP):

  • 模型是条件概率分布

  • 损失函数是对数损失函数,对应后验估计中的似然函数

  • 模型复杂度由模型的先验概率表示

似然函数先验概率乘积即对应贝叶斯最大后验估计的形式

参考https://www.zhihu.com/question/23536142

所以结构风险最小化就是求解优化问题:

minfF1Ni=1NL(yi,f(xi))+λJ(f)\min_{f\in\mathcal{F}}\frac{1}{N}\sum ^N_{i=1}L(y_i,f(x_i))+\lambda J(f)

算法

算法指的是学习模型的具体方法,即使用什么计算方法求解最优模型。

因为统计学习问题归结为最优化问题,所以统计学习的算法就是求解最优化问题的算法

  • 如果有显式的解析解,此最优化问题就比较简单

  • 如果没有,需要用数值计算方法求解,需要考虑如何保证找到全局最优解,并使求解过程高效

模型评估与模型选择

训练误差与测试误差

假设学习到的模型是Y=f^(X)Y=\hat{f}(X),训练误差是模型Y=f^(X)Y=\hat{f}(X)关于训练数据集的平均损失(NN是样本容量):

Remp(f^)=1Ni=1NL(yi,f^(xi))R_{emp}(\hat{f})=\frac{1}{N}\sum ^N_{i=1}L(y_i,\hat{f}(x_i))

测试误差是模型Y=f^(X)Y=\hat{f}(X)关于测试数据集的平均损失(NN'是测试样本容量):

etest(f^)=1Ni=1NL(yi,f^(xi))e_{test}(\hat{f})=\frac{1}{N'}\sum ^{N'}_{i=1}L(y_i,\hat{f}(x_i))
  • 训练误差的大小,对判断给定的问题是不是一个容易学习的问题是有意义的

  • 测试误差反映了学习方法对未知数据的预测能力,即泛化能力(generalization ability)

过拟合与模型选择

当模型复杂度增加时,训练误差会逐渐减小并趋向于0;测试误差会先减小,达到最小值后又会增大。

当模型复杂度过大时,就会出现过拟合。所以需要在学习时防止过拟合,选择复杂度适当的模型。

正则化与交叉验证

正则化

模型选择的典型方法是正则化(regularization)。正则化是结构风险最小化策略的实现,即在经验上加一个正则化项(regularizer)或罚项(penalty term)。

正则化项一般是模型复杂度单调递增函数,模型越复杂,正则化值越大。比如,正则化项可以是模型参数向量的范数。

minfF1Ni=1NL(yi,f(xi))+λJ(f)\min_{f\in \mathcal{F}}\frac{1}{N}\sum^N_{i=1}L(y_i,f(x_i))+\lambda J(f)

其中λ0\lambda \ge 0

正则化符合奥卡姆剃刀(Occam's razor)原理:在所有可能选择的模型中,能够很好地解释已知数据并且十分简单才是最好的模型,也就是应该选择的模型。

从贝叶斯估计的角度来看,正则化项对应于模型的先验概率。可以假设复杂的模型有较小的先验概率,简的模型有较大的先验概率。

交叉验证

泛化能力

泛化误差

泛化误差上界

生成模型与判别模型

分类问题

标注问题

回归问题

一个常问的问题:平方损失在做分类问题的损失函数时,有什么缺陷?

一方面,直观上看,平方损失函数对每一个输出结果都十分看重,而交叉熵只看重正确分类的结果。例如三分类问题,如果预测的是(a,b,c)(a,b,c),而实际结果是(1,0,0)(1,0,0),那么:

mse=(a1)2+b2+c2mse=(a-1)^2+b^2+c^2
crossentropy=1×loga0×logb0×logc=logacrossentropy=-1\times \log a-0\times \log b -0\times \log c=-\log a

所以,交叉熵损失的梯度只和正确分类有关。而平方损失的梯度和错误分类有关,除了让正确分类尽可能变大,还会让错误分类都变得更平均。实际中后面这个调整在分类问题中是不必要的,而回归问题上这就很重要。

另一方面,从理论角度分析。两个损失的源头不同。平方损失函数假设最终结果都服从高斯分布,而高斯分布的随机变量实际是一个连续变量,而非离散变量。如果假设结果 变量服从均值tt,方差σ\sigma的高斯分布,那么利用最大似然法可以优化其负对数似然,最终公式变为:

maxiN[12log(2πσ2)(tiy)22σ2]\max\sum^N_i[-\frac{1}{2}\log(2\pi \sigma^2)-\frac{(t_i-y)^2}{2\sigma ^2}]

除去与yy无关的项,剩下的就是平方损失函数。

感知机

k近邻法

朴素贝叶斯法

参考https://www.cnblogs.com/jiangxinyang/p/9378535.html

贝叶斯公式

p(θX)=p(Xθ)p(θ)p(X)p(\theta|X)=\frac{p(X|\theta)p(\theta)}{p(X)}

又可以写为:

posterior=likelihoodpriorevidenceposterior=\frac{likelihood * prior}{evidence}

其中:

  • posteriorposterior:通过样本XX得到参数θ\theta的概率,即后验概率

  • likelihoodlikelihood:通过参数θ\theta得到样本XX的概率,即似然函数。

  • priorprior:参数θ\theta的先验概率

  • evidenceevidencep(X)=p(Xθ)p(θ)dθp(X)=\int p(X|\theta)p(\theta)d\theta。样本XX发生的概率,是各种θ\theta条件下发生的概率的积分。

极大似然估计(MLE)

极大似然估计的核心思想是:认为当前发生的事件概率最大的事件。因此就可以给定的数据集,使得该数据集发生的概率最大来求得模型中的参数。

似然函数如下:

p(Xθ)=i=1np(xiθ)p(X|\theta)=\prod ^n_{i=1}p(x_i|\theta)

为便于计算,对似然函数两边取对数,得到对数似然函数:

LL(θ)=logp(Xθ)=logi=1np(xiθ)=i=1nlogp(xiθ)\begin{aligned} LL(\theta)&=\log p(X|\theta) \\ &= \log \prod ^n_{i=1}p(x_i|\theta) \\ &= \sum ^n_{i=1}\log p(x_i|\theta) \\ \end{aligned}

极大似然估计只关注当前的样本,也就是只关注当前发生的事情,不考虑事情的先验情况。由于计算简单,而且不需要关注先验知识,因此在机器学习中的应用非常广,最常见的就是逻辑回归

最大后验估计(MAP)

最大后验估计中引入了先验概率(先验分布属于贝叶斯学派引入的,像L1,L2正则化就是对参数引入了拉普拉斯先验分布和高斯先验分布),而且最大后验估计要求的是p(θX)p(\theta|X)

f(x)=argmaxθp(θX)=argmaxθp(Xθ)p(θ)p(X)=argmaxθp(Xθ)p(θ)\begin{aligned} f(x)&= \arg\max _{\theta}p(\theta|X) \\ &= \arg\max _{\theta}\frac{p(X|\theta)p(\theta)}{p(X)} \\ &= \arg\max _{\theta} p(X|\theta)p(\theta) \\ \end{aligned}

其中因为分母p(X)p(X)θ\theta无关,所以可以去掉,同样地,取log:

f(x)=argmaxθlogp(Xθ)p(θ)=argmaxθ{i=1nlogp(xiθ)+logp(θ)}\begin{aligned} f(x)&= \arg\max _{\theta} \log p(X|\theta)p(\theta) \\ &= \arg\max _{\theta} \{\sum ^n_{i=1}\log p(x_i|\theta)+\log p(\theta)\}\\ \end{aligned}

最大后验估计不只是关注当前的样本的情况,还关注已经发生过的先验知识

最大后验估计和最大似然估计的区别:

最大后验估计允许我们把先验知识加入到估计模型中,这在样本很少的时候是很有用的(因此朴素贝叶斯在较少的样本下就能有很好的表现),因为样本很少的时候我们的观测结果很可能出现偏差,此时先验知识会把估计的结果“拉”向先验,实际的预估结果将会在先验结果的两侧形成一个顶峰。通过调节先验分布的参数,比如beta分布的α\alphaβ\beta,我们还可以调节把估计的结果“拉”向先验的幅度,α\alphaβ\beta越大,这个顶峰越尖锐。这样的参数,我们叫做预估模型的“超参数”。

贝叶斯估计

贝叶斯估计和极大后验估计有点相似,都是以最大化后验概率为目的。区别在于:

  • 极大似然估计和极大后验估计都是只返回了的预估值

  • 极大后验估计在计算后验概率的时候,把分母p(X)p(X)忽略了,在进行贝叶斯估计的时候则不能忽略

  • 贝叶斯估计要计算整个后验概率的概率分布

决策树

logistic回归与最大熵模型

支持向量机

提升方法

gbdt

参考https://www.cnblogs.com/pinard/p/6140514.html

梯度提升树(Gradient Boosting Decison Tree, 以下简称GBDT)有很多简称,有GBT(Gradient Boosting Tree), GTB(Gradient Tree Boosting ), GBRT(Gradient Boosting Regression Tree), MART(Multiple Additive Regression Tree)等等。

gbdt概述

在Adaboost中,我们是利用前一轮迭代弱学习器的误差率来更新训练集的权重,这样一轮轮的迭代下去。

GBDT也是迭代,但是弱学习器限定了只能使用CART回归树模型,同时迭代思路和Adaboost也有所不同。

假设我们前一轮迭代得到的强学习器是ft1(x)f_{t-1}(x),损失函数是L(y,ft1(x))L(y, f_{t-1}(x)),本轮迭代的目标是找到一个CART回归树模型的弱学习器ht(x)h_t(x),使得本轮的损失函数L(y,ft(x))=L(y,ft1(x)+ht(x))L(y, f_{t}(x)) =L(y, f_{t-1}(x)+ h_t(x))最小。也就是说,要找到一个决策树,使得样本的损失L(yft1(x),ht(x))L(y-f_{t-1}(x), h_t(x))最小。

gbdt的负梯度拟合

tt轮的第ii个样本的损失函数的负梯度:

rti=[L(yi,f(xi)))f(xi)]f(x)=ft1(x)r_{ti} = -\bigg[\frac{\partial L(y_i, f(x_i)))}{\partial f(x_i)}\bigg]_{f(x) = f_{t-1}(x)}

利用

EM算法及其推广

隐马尔可夫模型

条件随机场

附录

矩阵

优化

拉格朗日乘子法

拉格朗日乘子法(Lagrange multipliers)是一种寻找多元函数在一组约束下的极值的方法。通过引入拉格朗日乘子,将dd个变量和kk个约束条件的最优化问题转化为具有d+kd+k个变量的无约束优化问题求解。

等式约束

假设xxdd维向量,要寻找xx的某个取值xx^*,使目标函数f(x)f(x)最小且同时满足g(x)=0g(x)=0的约束。

从几何角度看,目标是在由方程g(x)=0g(x)=0确定的d1d-1维曲面上,寻找能使目标函数f(x)f(x)最小化的点。

  1. 对于约束曲面g(x)=0g(x)=0上的任意点xx,该点的梯度ablag(x)abla g(x)正交于约束曲面

  2. 最优点xx^*,目标函数f(x)f(x)在该点的梯度ablaf(x)abla f(x^*)正交于约束曲面

对于第1条,梯度本身就与曲面的切向量垂直,是曲面的法向量,并且指向数值更高的等值线

证明:

参考http://xuxzmail.blog.163.com/blog/static/251319162010328103227654/

假设z=f(x,y)z=f(x,y)的等值线:Γ:f(x,y)=c\Gamma :f(x,y)=c,两边求微分:

df(x,y)=dcfxdx+fydy=0\begin{matrix} df(x,y)=dc & \\ \frac{\partial f}{\partial x} dx + \frac{\partial f}{\partial y} dy =0 & \\ & \end{matrix}

看成两个向量的内积:

fxdx+fydy={fx,fy}{dx,dy}=0\frac{\partial f}{\partial x} dx + \frac{\partial f}{\partial y} dy = \{\frac{\partial f}{\partial x}, \frac{\partial f}{\partial y}\}\cdot \{dx,dy\}=0

而内积ab=abcosθa\cdot b=|a||b|cos\theta为0说明夹角是90度,而{fx,fy}\{\frac{\partial f}{\partial x}, \frac{\partial f}{\partial y}\}是梯度向量,{dx,dy}\{dx,dy\}是等值线的切向量,所以梯度向量和切向量是垂直的。

对于第2条,可以用反证法,如下图,蓝色是g(x)=0g(x)=0,橙色是f(x)f(x)的等值线(图里假设f(x)=x2+y2f(x)=x^2+y^2),交点的ablaf(x)abla f(x^*)的梯度和g(x)g(x)的切面不垂直,那么,可以找到更小的等值线,使夹角更接近90度,也就是说,这个点不是真正的最优点xx^*

所以,在最优点xx^*处,梯度ablag(x)abla g(x)ablaf(x)abla f(x)的方向必然相同或相反,也就是存在λ0\lambda \neq 0,使得:

ablaf(x)+λg(x)=0abla f(x^*)+\lambda \nabla g(x^*)=0

其中λ\lambda是拉格朗日乘子,定义拉格朗日函数

L(x,λ)=f(x)+λg(x)L(x,\lambda)=f(x)+\lambda g(x)

其中L(x,λ)L(x,\lambda)对x的偏导ablaxL(x,λ)abla _xL(x,\lambda)置0,就得到:

ablaf(x)+λg(x)=0abla f(x)+\lambda \nabla g(x)=0

L(x,λ)L(x,\lambda)λ\lambda的偏导ablaλL(x,λ)abla _{\lambda}L(x,\lambda)置0,就得到

g(x)=0g(x)=0

所以,原约束问题可以转化为对L(x,λ)L(x,\lambda)的无约束优化问题。

不等式约束

考虑不等式约束g(x)0g(x)\le 0,最优点或者在边界g(x)=0g(x)=0上,或者在区域g(x)<0g(x)<0中。

  • 对于g(x)<0g(x)<0

相当于使f(x)f(x)取得最小值的点落在可行域内,所以约束条件相当于没有用,所以,直接对f(x)f(x)求极小值即可。因为L(x,λ)=f(x)+λg(x)L(x,\lambda)=f(x)+\lambda g(x),所以

ablaxL(x,λ)=f(x)+λg(x)abla _xL(x,\lambda)=\nabla f(x)+\lambda \nabla g(x)

因为g(x)<0g(x)<0,想要只让ablaf(x)=0abla f(x)=0,那么令 λ=0\lambda = 0 即可。

  • 对于g(x)=0g(x)=0

这就变成了等式约束,且此时ablaf(x)abla f(x^*)ablag(x)abla g(x^*)反向相反。因为在g(x)=0g(x)=0越往里值是越小的,而梯度是指向等值线高的方向,所以梯度是指向外的。而f(x)f(x)的可行域又在g(x)g(x)的里面和边界上,我们要找的是f(x)f(x)的最小值,所以f(x)f(x)的梯度是指向内部的。

ablaf(x)+λg(x)=0abla f(x)+\lambda \nabla g(x)=0,两个又是反向的,所以λ>0\lambda>0

结合g(x)0g(x)\le 0g(x)=0g(x)=0两种情况的结论,就得到了KKT(Karush-Kuhn-Tucker)条件

g(x)=0,λ>0g(x)<0,λ=0}{g(x)0λ0λg(x)=0\left.\begin{matrix} g(x)=0, \lambda >0 \\ g(x)<0, \lambda = 0 \end{matrix}\right\}\Rightarrow \left\{\begin{matrix} g(x)\le 0 \\ \lambda \ge 0\\ \lambda g(x)=0 \end{matrix}\right.

其中λg(x)=0\lambda g(x)=0是因为λ\lambdag(x)g(x)至少一个是0,而且不能都不是0。

以上三个条件有各自的名字:

  • Primal feasibility(原始可行性):g(x)0g(x)\le 0

  • Dual feasibility(对偶可行性):λ0\lambda \ge 0

  • Complementary slackness:λg(x)=0\lambda g(x)=0

带等式和不等式约束的拉格朗日乘子法

对于多个约束的情形,mm个等式约束和nn个不等式约束,可行域DRd\mathbb{D}\subset \mathbb{R}^d非空的优化问题:

minxf(x)\min_{x}f(x)
s.thi(x)=0,i=1,...,mgj(x)0,j=1,...,n\begin{matrix} s.t & h_i(x)=0,i=1,...,m \\ & g_j(x)\le 0,j=1,...,n \end{matrix}

引入拉格朗日乘子λ=(λ1,λ2,...,λm)T\lambda=(\lambda_1,\lambda_2,...,\lambda_m)^Tμ=(μ1,μ2,...,μn)T\mu=(\mu_1,\mu_2,...,\mu_n)^T,相应的拉格朗日函数为:

L(x,λ,μ)=f(x)+i=1mλihi(x)+j=1nμjgj(x)L(x,\lambda, \mu)=f(x)+\sum^m_{i=1}\lambda_ih_i(x)+\sum^n_{j=1}\mu_jg_j(x)

由不等式约束引入的KKT条件(j=1,2,...n)(j=1,2,...n)

{gj(x)0μj0μjgj(x)=0\left\{\begin{matrix} g_j(x)\le 0\\ \mu_j\ge 0\\ \mu_jg_j(x)=0 \end{matrix}\right.

梯度下降

《统计学习方法》的视角

假设f(x)f(x)有一阶连续偏导,对于无约束的最优化问题而言:minxRnf(x)\min _{x\in R^n}f(x)

f(x)f(x)x(k)x^{(k)}附近的一阶泰勒展开如下,其中gk=g(x(k))=f(x(k))g_k=g(x^{(k)})=\nabla f(x^{(k)})f(x)f(x)x(k)x^{(k)}的梯度:

f(x)=f(x(k))+gkT(xx(k))f(x)=f(x^{(k)})+g^T_k(x-x^{(k)})

所以对于x=x(k+1)x=x^{(k+1)}

f(x(k+1))=f(x(k))+gkT(x(k+1)x(k))f(x^{(k+1)})=f(x^{(k)})+g^T_k(x^{(k+1)}-x^{(k)})

x(k+1)=x(k)+λkpkx^{(k+1)}=x^{(k)}+\lambda _kp_kpkp_k是搜索方向,λk\lambda _k是步长,代入上式,有

f(x(k+1))=f(x(k))+gkT(x(k)+λkpkx(k))=f(x(k))+gkTλkpk\begin{aligned} f(x^{(k+1)})&=f(x^{(k)})+g^T_k(x^{(k)}+\lambda _kp_k-x^{(k)}) \\ &=f(x^{(k)})+g^T_k\lambda _kp_k \\ \end{aligned}

为了让每次迭代的函数值变小,可以取pk=f(x(k))p_k=-\nabla f(x^{(k)})

λk\lambda_k看成是可变化的,所以需要搜索λk\lambda _k使得

f(x(k)+λkpk)=minλ0f(x(k)+λpk)f(x^{(k)}+\lambda_kp_k)=\min_{\lambda \ge0}f(x^{(k)}+\lambda p_k)

梯度下降法:

输入:目标函数f(x)f(x),梯度g(x)=f(x)g(x)=\nabla f(x),精度要求ε\varepsilon

输出:f(x)f(x)的极小点xx^*

  1. 取初始值x(0)Rnx^{(0)}\in R^n,置k=0k=0

  2. 计算f(x(k))f(x^{(k)})

  3. 计算梯度gk=g(x(k))g_k=g(x^{(k)}),当gk<ε\left \| g_k \right \|<\varepsilon,则停止计算,得到近似解x=x(k)x^*=x^{(k)};否则,令pk=g(x(k))p_k=-g(x^{(k)}),求λk\lambda_k使得 f(x(k)+λkpk)=minλ0f(x(k)+λpk)f(x^{(k)}+\lambda_kp_k)=\min_{\lambda \ge0}f(x^{(k)}+\lambda p_k)

  4. x(k+1)=x(k)+λkpkx^{(k+1)}=x^{(k)}+\lambda_kp_k,计算f(x(k+1))f(x^{(k+1)})f(x(k+1))f(x(k))<ε\left \| f(x^{(k+1)}) -f(x^{(k)}) \right \|<\varepsilonx(k+1)x(k)<ε\left \| x^{(k+1)}-x^{(k)} \right \|<\varepsilon时,停止迭代,令x=x(k+1)x^*=x^{(k+1)}

  5. 否则,置k=k+1k=k+1,转第3步

只有当目标函数是凸函数时,梯度下降得到的才是全局最优解

《机器学习》的视角

梯度下降是一阶(first order)(只用一阶导,不用高阶导数)优化方法,是求解无约束优化问题最简单、最经典的方法之一。

考虑无约束优化问题minxf(x)\min_xf(x)f(x)f(x)是连续可微函数,如果能构造一个序列x0,x1,x2,...x^0,x^1,x^2,...满足

f(xt+1)<f(xt), t=0,1,2,...f(x^{t+1}) < f(x^t),\ t=0,1,2,...

那么不断执行这个过程,就可以收敛到局部极小点,根据泰勒展开有:

f(x)=f(x(k))+f(x(k))T(xx(k))f(x+Δx)=f(x(k))+f(x(k))T(x+Δxx(k))=f(x(k))+f(x(k))T(xx(k))+f(x(k))TΔx=f(x)+f(x(k))TΔx\begin{aligned} f(x)&=f(x^{(k)})+\nabla f(x^{(k)})^T(x-x^{(k)}) \\ f(x+\Delta x) &= f(x^{(k)})+\nabla f(x^{(k)})^T(x+\Delta x-x^{(k)}) \\ & = f(x^{(k)})+\nabla f(x^{(k)})^T(x-x^{(k)}) + \nabla f(x^{(k)})^T\Delta x \\ & = f(x) + \nabla f(x^{(k)})^T\Delta x \end{aligned}

ablaf(x(k))TΔxabla f(x^{(k)})^T\Delta x是一个标量,其转置等于自己,所以

f(x+Δx)=f(x)+ΔxTf(x(k))f(x+\Delta x)=f(x) + \Delta x^T\nabla f(x^{(k)})

想要让f(x+Δx)<f(x)f(x+\Delta x)< f(x),只需要令:

Δx=γf(x)\Delta x=-\gamma \nabla f(x)

其中的步长γ\gamma是一个小常数

如果f(x)f(x)满足LL-Lipschitz条件,也就是说对于任意的xx,存在常数LL,使得f(x)L\left \| \nabla f(x) \right \|\le L成立,那么设置步长为12L\frac{1}{2L}就可以确保收敛到局部极小点。

同样地,当目标函数是凸函数时,局部极小点就对应全局最小点,此时梯度下降可以确保收敛到全局最优解。

牛顿法

二阶导基本性质

对于点x=x0x=x_0

  • 一阶导f(x0)=0f'(x_0)=0时,如果二阶导f(x0)>0f''(x_0)>0,那么f(x0)f(x_0)是极小值,x0x_0是极小点

  • 一阶导f(x0)=0f'(x_0)=0,如果二阶导f(x0)<0f''(x_0)<0,那么f(x0)f(x_0)是极大值,x0x_0是极大点

  • 一阶导f(x0)=0f'(x_0)=0,如果二阶导f(x0)=0f''(x_0)=0,那么x0x_0是鞍点

证明:

对于任意x1x_1,根据二阶泰勒展开,有

f(x1)=f(x0)+f(x0)(x1x0)+12f(x0)(x1x0)2+...+Rn(x1)f(x_1)=f(x_0)+f'(x_0)(x_1-x_0)+\frac{1}{2}f''(x_0)(x_1-x_0)^2+...+R_n(x_1)

因为f(x0)>0f''(x_0)>0f(x0)=0f'(x_0)=0,所以,不论x1>x0x_1> x_0还是x1<x0x_1< x_0,总有f(x1)>f(x0)f(x_1)>f(x_0),也就是周围的函数值都比f(x0)f(x_0)大,而x0x_0又是极值点,所以是极小点。

牛顿法

对于矩阵形式,xx是一个nx1的列向量,H(x)H(x)f(x)f(x)的海赛矩阵,即二阶导,shape是n×nn\times n

f(x)=f(x(x))+gkT(xx(k))+12(xx(k))TH(x(k))(xx(k))f(x)=f(x^{(x)})+g^T_k(x-x^{(k)})+\frac{1}{2}(x-x^{(k)})^TH(x^{(k)})(x-x^{(k)})

函数f(x)f(x)有极值的必要条件是在极值点处一阶导为0,特别地,当H(x(k))H(x^{(k)})是正定矩阵时(二阶导大于0),是极小值。

牛顿法利用极小点的必要条件ablaf(x)=0abla f(x)=0,每次迭代从点x(k)x^{(k)}开始,求目标函数极小点,作为第k+1k+1次迭代值x(k+1)x^{(k+1)},具体地,假设ablaf(x(k+1))=0abla f(x^{(k+1)})=0,有

f(x)=f(x(x))+gkT(xx(k))+12(xx(k))TH(x(k))(xx(k))=f(x(x))+[gkT+12(xx(k))TH(x(k))](xx(k))=f(x(x))+[gk+12H(x(k))(xx(k))]T(xx(k))\begin{aligned} f(x) &= f(x^{(x)})+g^T_k(x-x^{(k)})+\frac{1}{2}(x-x^{(k)})^TH(x^{(k)})(x-x^{(k)}) \\ &=f(x^{(x)})+[g^T_k+\frac{1}{2}(x-x^{(k)})^TH(x^{(k)})](x-x^{(k)}) \\ &=f(x^{(x)})+[g_k+\frac{1}{2}H(x^{(k)})(x-x^{(k)})]^T(x-x^{(k)}) \\ \end{aligned}

把其中的gk+12H(x(k))(xx(k))g_k+\frac{1}{2}H(x^{(k)})(x-x^{(k)})看成一阶导,则上式就是一阶泰勒展开。记Hk=H(x(k))H^k=H(x^{(k)}),令x=x(k+1)x=x^{(k+1)},令一阶导为0:

gk+12Hk(x(k+1)x(k))=0gk=12Hk(x(k+1)x(k))2Hk1gk=x(k+1)x(k)x(k+1)=2Hk1gk+x(k)\begin{matrix} g_k+\frac{1}{2}H^k(x^{(k+1)}-x^{(k)})=0 \\ g_k=-\frac{1}{2}H^k(x^{(k+1)}-x^{(k)}) \\ -2H^{-1}_kg_k=x^{(k+1)}-x^{(k)} \\ x^{(k+1)} = -2H^{-1}_kg_k+x^{(k)} \end{matrix}

可以无视这个2,变成:

x(k+1)=x(k)Hk1gkx^{(k+1)}=x^{(k)}-H^{-1}_kg_k

或者

x(k+1)=x(k)+pkx^{(k+1)}=x^{(k)}+p_k

其中,

Hkpk=gkH_kp_k=-g_k

牛顿法

输入:目标函数f(x)f(x),梯度g(x)=f(x)g(x)=\nabla f(x),海赛矩阵H(x)H(x),精度要求ε\varepsilon

输出:f(x)f(x)的极小点xx^*

  1. 取初始点x(0)x^{(0)},置k=0k=0

  2. 计算gk=g(x(k))g_k=g(x^{(k)})

  3. gk<ε\left \| g_k \right \|<\varepsilon,则停止计算,得到近似解x=x(k)x^*=x^{(k)}

  4. 计算Hk=H(x(k))H_k=H(x^{(k)}),并求pkp_k,满足 Hkpk=gkH_kp_k=-g_k

  5. x(k+1)=x(k)+pkx^{(k+1)}=x^{(k)}+p_k

  6. k=k+1k=k+1,转到第2步

其中的步骤4,求pkp_k时,pk=Hk1gkp_k=-H^{-1}_kg_k需要求解Hk1H^{-1}_k很复杂。

拟牛顿法的思路

基本想法就是通过一个nn阶矩阵Gk=G(x(k))G_k=G(x^{(k)})来近似代替H1(x(k))H^{-1}(x^{(k)})

DFP(Davidon-Fletcher-Powell)

BFGS(Broydon-Fletcher-Goldfarb-Shanno)

拉格朗日对偶性

信息论相关

凸集

假设SS为在实或复向量空间的集。若对于所有x,ySx,y\in S和所有的t[0,1]t\in [0,1]都有tx+(1t)yStx+(1-t)y\in S,则SS称为凸集。

也就是说,SS中任意两点间的线段都属于SS

性质:

如果SS是凸集,对于任意的u1,u2,...,urSu_1,u_2,...,u_r\in S,以及所有的非负λ1,λ2,...,λr\lambda_1,\lambda_2,...,\lambda_r满足λ1+λ2+...+λr=1\lambda_1+\lambda_2+...+\lambda_r=1,都有k=1rλkukS\sum_{k=1}^{r}\lambda_ku_k\in S。这个组合称为u1,u2,...,uru_1,u_2,...,u_r的凸组合。

凸函数

函数是凸函数:曲线上任意两点x和y所作割线(与函数图像有两个不同交点的直线,如果只有一个交点,那就是切线)一定在这两点间的函数图象的上方:

tf(x)+(1t)f(y)f(tx+(1t)y),0t1tf(x)+(1-t)f(y)\ge f(tx+(1-t)y),0\le t \le 1

有如下几个常用性质:

  • 一元可微函数在某个区间上是凸的,当且仅当它的导数在该区间上单调不减

  • 一元连续可微函数在区间上是凸的,当且仅当函数位于所有它的切线的上方:对于区间内的所有xxyy,都有 f(y)f(x)+f(x)(yx)f(y) ≥ f(x) + f '(x) (y − x)(右边就是一阶泰勒展开)。特别地,如果f(c)=0f '(c) = 0,那么f(c)f(c)f(x)f(x)的最小值。

  • 一元二阶可微的函数在区间上是凸的,当且仅当它的二阶导数是非负的;这可以用来判断某个函数是不是凸函数。如果它的二阶导数是正数,那么函数就是严格凸的,但反过来不成立。

  • 多元二次可微的连续函数在凸集上是凸的,当且仅当它的黑塞矩阵在凸集的内部是半正定的

KL散度

熵的小结:https://blog.csdn.net/haolexiao/article/details/70142571

相对熵(relative entropy)又称为KL散度(Kullback–Leibler divergence,简称KLD),信息散度(information divergence)信息增益(information gain)

KL散度是两个概率分布P和Q差别的非对称性的度量。 KL散度是用来度量使用基于Q的编码来编码来自P的样本 平均所需的额外的位元数。 典型情况下,P表示数据的真实分布Q表示数据的理论分布,模型分布,或P的近似分布

注意:DKL(PQ)D_{KL}(P||Q)是指的用分布Q近似数据的真实分布P,先写P再写Q,公式里没有-ln的时候,就是p/q

对于离散随机变量:

DKL(PQ)=iP(i)lnP(i)Q(i)=iP(i)lnQ(i)P(i)D_{KL}(P||Q)=\sum_{i}P(i)ln\frac{P(i)}{Q(i)}=-\sum _{i}P(i)ln\frac{Q(i)}{P(i)}

KL散度仅当概率P和Q各自总和均为1,且对于任何i皆满足Q(i)>0Q(i)>0P(i)>0P(i)>0时,才有定义。如果出现0ln00ln0,当做0

对于连续随机变量:

DKL(PQ)=p(x)lnp(x)q(x)dxD_{KL}(P||Q)=\int ^{\infty }_{-\infty}p(x)ln\frac{p(x)}{q(x)}dx

性质:

KL散度大于等于0

证明:

先了解一下Jensen不等式:

如果φ\varphi是一个凸函数,那么有:

φ(E(x))E(φ(x))\varphi(E(x))\le E(\varphi(x))

对于离散随机变量i=1nλi=1\sum_{i=1}^n \lambda_i=1

φ(i=1ng(xi)λi)i=1nφ(g(xi))λi\varphi(\sum _{i=1}^n g(x_i)\lambda_i)\le \sum^{n}_{i=1}\varphi(g(x_i))\lambda_i

当我们取g(x)=xg(x)=xλi=1/n\lambda_i=1/nφ(x)=log(x)\varphi(x)=log(x)时,就有

log(i=1nxin)i=1nlog(xi)nlog(\sum^{n}_{i=1}\frac{x_i}{n})\ge \sum^{n}_{i=1}\frac{log(x_i)}{n}

对于连续随机变量,如果f(x)是非负函数,且满足(f(x)是概率密度函数):

f(x)dx=1\int^{\infty}_{-\infty}f(x)dx=1

如果φ\varphi在g(x)的值域中是凸函数,那么

φ(g(x)f(x)dx)φ(g(x))f(x)dx\varphi(\int^{\infty}_{-\infty}g(x)f(x)dx)\le \int ^{\infty}_{-\infty}\varphi(g(x))f(x)dx

特别地,当g(x)=xg(x)=x时,有

φ(xf(x)dx)φ(x)f(x)dx\varphi(\int^{\infty}_{-\infty}xf(x)dx)\le \int ^{\infty}_{-\infty}\varphi(x)f(x)dx

回到这个问题中,g(x)=q(x)p(x)g(x)=\frac{q(x)}{p(x)}φ(x)=logx\varphi(x)=-logx是一个严格凸函数,那么φ(g(x))=logq(x)p(x)\varphi(g(x))=-log\frac{q(x)}{p(x)},所以

DKL(PQ)=p(x)lnp(x)q(x)dx=p(x)(lnq(x)p(x))dx=(lnq(x)p(x))p(x)dxln(q(x)p(x)p(x)dx)ln(q(x)dx)=ln1=0\begin{aligned} D_{KL}(P||Q)&=\int ^{\infty }_{-\infty}p(x)ln\frac{p(x)}{q(x)}dx \\ &=\int ^{\infty }_{-\infty}p(x)(-ln\frac{q(x)}{p(x)})dx \\ &=\int ^{\infty }_{-\infty}(-ln\frac{q(x)}{p(x)})p(x)dx \\ &\ge -ln(\int^{\infty }_{-\infty}\frac{q(x)}{p(x)}p(x)dx) \\ &\ge -ln(\int^{\infty }_{-\infty}q(x)dx) \\ &=-ln1=0 \end{aligned}

采样

假设有一个很复杂的概率密度函数p(x)p(x),求解随机变量基于此概率下的某个函数期望,即:

Exp(x)[f(x)]E_{x\sim p(x)}[f(x)]

求解有两种方法:

解析法:

将上式展开成积分,并通过积分求解:

xp(x)f(x)dx\int _x p(x)f(x)dx

对于简单的分布,可以直接这么做。但dnn就不可行了。

蒙特卡洛法:

根据大数定理,当采样数量足够大时,采样的样本就可以无限近似地表示原分布:

1Nxip(x),i=1Nf(xi)\frac{1}{N}\sum ^N_{x_i\sim p(x),i=1}f(x_i)

拒绝采样

拒绝采样又叫接受/拒绝采样(Accept-Reject Sampling),对于目标分布p(x)p(x),选取一个容易采样的参考分布q(x)q(x),使得对于任意xx都有p(x)Mq(x)p(x)\le M\cdot q(x),则可以按如下过程采样:

  • 从参考分布q(x)q(x)中随机抽取一个样本xix_i

  • 从均匀分布U(0,1)U(0,1)产生一个随机数uiu_i

  • 如果uip(xi)Mq(xi)u_i\le \frac{p(x_i)}{Mq(x_i)},则接受样本xix_i,否则拒绝。

重新进行如上3个步骤,直到新产生的样本xix_i被接受

其中的第三步是因为p(x)Mq(x)p(x)\le M\cdot q(x),所以p(xi)Mq(xi)1\frac{p(x_i)}{Mq(x_i)}\le 1,说明只有函数值在p(xi)p(x_i)下方的xix_i才接受,所以xip(x)x_i\sim p(x)。相当于为目标分布p(x)p(x)选一个包络函数Mq(x)M\cdot q(x),包络函数紧,每次采样时样本被接受的概率越大,采样效率越高。实际应用中还有自适应拒绝采样等。

重要性采样

强化学习经常使用重要性采样。重要性采样主要应用在一些难以直接采样的数据分布上。

我们令待采样分布为p(x)p(x),有另一个简单可采样且定义域和p(x)p(x)相同的概率密度函数为p~(x)\tilde {p}(x),可以得到:

Exp(x)[f(x)]=xp(x)f(x)dx=xp~(x)p(x)p~(x)f(x)dx=Exp~(x)[p(x)p~(x)f(x)]1Nxip~(x),i=1Np(x)p~(x)f(x)\begin{aligned} E_{x\sim p(x)}[f(x)]&=\int _x p(x)f(x)dx \\ &=\int _x \tilde{p}(x) \frac{p(x)}{\tilde {p}(x)}f(x)dx \\ &=E_{x\sim \tilde{p}(x)} [\frac{p(x)}{\tilde {p}(x)}f(x)] \\ & \simeq \frac{1}{N} \sum ^N_{x_i\sim \tilde{p}(x),i=1}\frac{p(x)}{\tilde{p}(x)}f(x) \end{aligned}

因此,只需要从简单分布p~(x)\tilde {p}(x)中采样,然后分别计算p(x)p(x)p~(x)\tilde {p}(x)f(x)f(x)就可以了。

最好选择一个和原始分布尽量接近的近似分布进行采样。

例如,要对一个均值1,方差1的正态分布进行采样,有两个可选的分布:均值1方差0.5和均值1方差2。

从图像上可以看到方差为0.5的过分集中在均值附近,而且方差为2的与原分布重合度较高,所以应该选方差为2的。

马尔可夫蒙特卡洛采样(MCMC)

如果是高维空间的随机向量,拒绝采样和重要性采样经常难以找到合适的参考分布,采样效率低下(样本的接受概率低或者重要性权重低),此时可以考虑马尔可夫蒙特卡洛采样法。

蒙特卡洛法指基于采样的数值近似求解方法,马尔可夫链用于进行采样。

基本思想是:

  • 针对待采样的目标分布,构造一个马尔可夫链,使得该马尔可夫链的平稳分布就是目标分布;

  • 然后从任何一个初始状态出发,沿着马尔可夫链进行状态转移,最终得到的状态转移序列会收敛到目标分布

  • 由此可以得到目标分布的一系列样本

核心点是如何构造合适的马尔可夫链,即确定马尔可夫链的状态转移概率。

Metropolis-Hastings采样

对于目标分布p(x)p(x),首先选择一个容易采样的参考条件分布q(xx)q(x^*|x),并令

A(x,x)=min{1,p(x)q(xx)p(x)q(xx)}A(x,x^*)=\min \{1,\frac{p(x^*)q(x|x^*)}{p(x)q(x^*|x)}\}

然后根据如下过程采样:

  • 随机选一个初始样本x(0)x^{(0)}

  • For t = 1, 2, 3, ... :

    • 根据参考条件分布q(xx(t1))q(x^*|x^{(t-1)})抽取一个样本xx^*

    • 根据均匀分布U(0,1)U(0,1)产生随机数uu

    • u<A(x(t1),x)u < A(x^{(t-1)},x^*),则令x(t)=xx^{(t)}=x^*【接受新样本】,否则令x(t)=x(t1)x^{(t)}=x^{(t-1)}【拒绝新样本,维持旧样本】

可以证明,上述过程得到的样本序列{...,x(t1),x(t),...}\{...,x^{(t-1)},x^{(t)},...\}最终会收敛到目标分布p(x)p(x)

吉布斯采样

吉布斯采样是Metropolis-Hastings算法的一个特例,核心思想是只对样本的一个维度进行采样和更新。对于目标分布p(x)p(x),其中x=(x1,x2,...,xd)x=(x_1,x_2,...,x_d)是一个多维向量,按如下过程采样:

  • 随机选择初始状态x(0)=(x1(0),x2(0),...,xd(0))x^{(0)}=(x^{(0)}_1,x^{(0)}_2,...,x^{(0)}_d)

  • For t = 1, 2, 3, ... :

    • 对于前一步产生的样本x(t1)=(x1(t1),x2(t1),...,xd(t1))x^{(t-1)}=(x^{(t-1)}_1,x^{(t-1)}_2,...,x^{(t-1)}_d),依次采样和更新每个维度的值,即依次抽取分量x1(t)p(x1x2(t1),x3(t1),...,xd(t1))x^{(t)}_1\sim p(x_1|x^{(t-1)}_2,x^{(t-1)}_3,...,x^{(t-1)}_d), x2(t)p(x1x1(t1),x3(t1),...,xd(t1))x^{(t)}_2\sim p(x_1|x^{(t-1)}_1,x^{(t-1)}_3,...,x^{(t-1)}_d), ..., xd(t)p(x1x1(t1),x2(t1),...,xd1(t1))x^{(t)}_d\sim p(x_1|x^{(t-1)}_1,x^{(t-1)}_2,...,x^{(t-1)}_{d-1})

    • 形成新的样本x(t)=(x1(t),x2(t),...,xd(t))x^{(t)}=(x^{(t)}_1,x^{(t)}_2,...,x^{(t)}_d)

同样可以证明,上述过程得到的样本序列{...,x(t1),x(t),...}\{...,x^{(t-1)},x^{(t)},...\}会收敛到目标分布p(x)p(x)。另外,步骤2中对样本每个维度的抽样和更新操作,不是必须按下标顺序进行的,可以是随机顺序。

注意点:

  • 拒绝采样中,如果某一步中采样被拒绝,则该步不会产生新样本,需要重新采样。但MCMC采样每一步都会产生一个样本,只是有时候这个样本与之前的样本一样而已。

  • MCMC采样是在不断迭代过程中逐渐收敛到平稳分布的。实际应用中一般会对得到的样本序列进行"burn-in"处理,即截除掉序列中最开始的一部分样本,只保留后面的样本。

MCMC得到的样本序列中相邻的样本是不独立的,因为后一个样本是由前一个样本根据特定的转移概率得到的。如果仅仅是采样,并不要求样本之间相互独立。

如果确实需要产生独立同分布的样本,可以:

  • 同行运行多条马尔可夫链,这样不同链上的样本是独立的

  • 或者在同一条马尔可夫链上每隔若干个样本才选取一个,这样选出来的样本也是近似独立的。

最后更新于