纪念我第一个博客的碎碎念
先前我花了四五个月的业余时间学习了Ng的机器学习公开课,学习的过程中我就在想,如果我能把这个课程啃完,就开始写一些博客,把自己的所得记录下来,现在是实现的时候了。也如刘未鹏的《暗时间》里所说,哪怕更新频率很低,也应该坚持(从现在开始)写博客,记录有价值的东西(思考的产物),好处多多。我没有大神的气场,只是觉得,就算没人看,作为自己的备忘也不错,侥幸能坚持很久的话,日积月累下也会非常赏心悦目。
为了完成这个课程,我捡回了概率和线代以及统计学的知识,英语对我而言也从一门符号学真正地变成了一门语言,虽然我现在的英语水平还是不怎么样,但是我真正意识到了它的重要性,现在英语对我是有吸引力的,我想要接触和学习它,而不是先前那样我对它有一种先入为主的排斥感。
对我来说学习过程中仅有讲义和授课视频还是不够的,为了寻找更多的资料,我先是找到张雨石的公开课笔记,在惊叹与佩服中学习着,然后在一次无意中点开了他的其他文章,阅读着他平稳中带着精彩的经历,发现这么厉害一个人也在佩服着其他人,也在不断地学习、思考、进步着,那可真是一山更比一山高,这或许也正是他们出色的原因。
在开始这个课程的学习的同时,我也开始在阅读一些使心智变成熟与时间管理这样的书,而在张雨石一年又一年的书单中,我记下了我打算看的书的名字并标上了星级,其中我对《暗时间》这本书的名字十分感兴趣,立刻就下载来看了,结果确实不失所望,作者刘未鹏也是一个非常棒的人,他的个人经历、一些方法论以及后面的对一些实在的算法的分析解读,都让我看得非常开心和佩服,而看完这本《暗时间》之后,我的书单又加长了一些…并且不出所料,他也有比较推崇的人,是一位博客名前缀为g9的大大,而在搜到g9大大的博客,看到“负暄琐话”4个字和他给出的需要解码的email,以及他数量众多的博客之后,我再次窃喜:又是一位牛大啊…
这个课程的意义于我而言不仅仅是获得了课程中的知识那么简单。我仿佛因此打开了许多大门:我不再排斥英语,我开始阅读,我开始思考自身,我主动地学习,我透过博客与书籍接触到了很棒很有趣的人,并为之感到欢喜。我希望自己能够走近他们,最终成为他们之中的一员。
前言
首先是CS229这个课程的安排:http://cs229.stanford.edu/schedule.html
课程分为4个部分,分别是监督学习、学习理论、无监督学习、强化学习与控制。
以及课程资料:http://cs229.stanford.edu/materials.html
资料中包含了课程讲义(相当于课本),作业及解答,复习资料(这里的复习资料指的是此课程所需要的预备知识,包括贯穿全部的线性代数、概率论,SVM部分需要的凸优化知识,还有高斯进程,隐马尔科夫模型),与一些追加材料,大概是一些算法应用的例子和一些要注意的东西。
其中我想在一开始就说明的是,强化学习与控制部分是最后5个教学视频里的内容,但奇怪的是这个部分的讲义却只有一份,对应着第16个视频的内容,仅是马尔科夫过程的简介,往后4个视频里讲解的内容都没有讲义了,我在网上没找到这部分的讲义,后来我小心翼翼地写了一份邮件想要发到上面网页中提供的邮箱请他们把后面的讲义发给我,奈何在翻了墙和换了邮箱的情况下,邮件还是被退回,只好作罢。(求后面的讲义啊…)
接下来是网易公开课的视频:http://open.163.com/special/opencourse/machinelearning.html
这个课程有中文字幕,等于在听老师讲课。这个网页里打包好的课件中讲义、作业与解答是全的,复习资料少了三个,没有追加材料,如果要下载,建议把这个资料包下了,缺的又觉得有必要看的,再到CS229官网去下。
除了课本与授课老师之外,支撑着我完成已有讲义的学习的,还有这两份参考书,两份笔记都与讲义与授课内容高度统一,各有特点
张雨石的机器学习笔记专栏:http://blog.youkuaiyun.com/column/details/ml-ng-record.html
这份笔记注重授课内容的复述,可以看成是授课内容的文字版
还有这位JerryLead:http://www.cnblogs.com/jerrylead/default.html?page=3
这份笔记的突出之处在于这里有笔者自己的思考(比如对一些看似很自然的地方提出疑问或者把一些微妙的断点补上)与更高的完成度(Ng有时候会因为时间原因或者课程安排而没有把一些内容说满,这里完成了不少这样的内容)
他们都是我难望项背的人,这两份笔记也都给了我很大的帮助,放在这里没用对比优劣的意思,只是想把这些都记录下来。
1、线性回归(Linear Regression)
1.1、线性回归模型与解决方案
考虑下面的情况,这里给了一个房屋面积和价格的数据表:
我们的目标⟹\Longrightarrow⟹找到一条直线预测新来的房屋价格
怎么画这条直线⟹\Longrightarrow⟹参数θ\thetaθ决定
怎么对θ\thetaθ取值使预测结果尽量准确⟹\Longrightarrow⟹最小化J(θ)J(\theta)J(θ),使每一个房屋的预测值与实际值之差的平方和最小,即误差最小
这样我们的目标就已经转移到了minJ(θ)minJ(\theta)minJ(θ)。
1.2、方案可靠性研究——最小二乘的概率解释
这里有个问题,既然是要让误差最小,为什么不是计算预测值与实际值之差的绝对值∣hθ(x(i))−y(i)∣\left|h_\theta{(x^{(i)})}-y^{(i)}\right|hθ(x(i))−y(i)之和呢?
一个是因为比较难算,还有一个就是因为这里有个对J(θ)J(\theta)J(θ)的概率解释。
首先,让我们承认误差的存在,当预测值加上一个误差时,才能得到实际值:
y(i)=hθ(x(i))+ϵ(i)=θTx(i)+ϵ(i)(4)y^{(i)}=h_\theta(x^{(i)})+\epsilon^{(i)}=\theta^Tx^{(i)}+\epsilon^{(i)}\tag{4}y(i)=hθ(x(i))+ϵ(i)=θTx(i)+ϵ(i)(4)
第二,我们假设误差ϵ\epsilonϵ服从正态分布(Normal distribution),也称为高斯分布(Gaussian distribution):ϵ(i)\epsilon^{(i)}ϵ(i)~N(0,σ2)N(0,\sigma^2)N(0,σ2)
误差服从正态分布这个假设有两个原因,一个是影响误差的因素有很多,这些因素都是随机分布的,但是它们在整体上会趋向于正态分布,另一个是因为在把误差假设为服从正态分布后,相应的工作一般都能取得比较好的效果,虽然它们还是没有非常精确,但是已经足够了。
通过这个假设以及正态分布的公式f(x)=12πσexp(−(x−μ)22σ2)f(x)=\frac{1}{\sqrt{2\pi}\sigma}\exp(-\frac{(x-\mu)^2}{2\sigma^2})f(x)=2πσ1exp(−2σ2(x−μ)2),我们可以得到ϵ(i)\epsilon^{(i)}ϵ(i)的概率密度:
p(ϵ(i))=12πσexp(−(ϵ(i)−0)22σ2)=12πσexp(−(ϵ(i))22σ2)(5)p(\epsilon^{(i)})=\frac{1}{\sqrt{2\pi}\sigma}\exp(-\frac{(\epsilon^{(i)}-0)^2}{2\sigma^2})=\frac{1}{\sqrt{2\pi}\sigma}\exp(-\frac{(\epsilon^{(i)})^2}{2\sigma^2})\tag{5}p(ϵ(i))=2πσ1exp(−2σ2(ϵ(i)−0)2)=2πσ1exp(−2σ2(ϵ(i))2)(5)
再把(4)式代进来,我们可以得到:
p(y(i)∣x(i);θ)=12πσexp(−(y(i)−θTx(i))22σ2)(6)p(y^{(i)}\mid x^{(i)};\theta)=\frac{1}{\sqrt{2\pi}\sigma}\exp(-\frac{(y^{(i)}-\theta^Tx^{(i)})^2}{2\sigma^2})\tag{6}p(y(i)∣x(i);θ)=2πσ1exp(−2σ2(y(i)−θTx(i))2)(6)
p(y(i)∣x(i);θ)p(y^{(i)}\mid x^{(i)};\theta)p(y(i)∣x(i);θ)可以读作“给定参数θ\thetaθ时,在x(i)x^{(i)}x(i)发生了的情况下,y(i)y^{(i)}y(i)发生的概率是多少”。注意里面的分号,分号是用来区别条件与参数的,如果把分号改成了逗号,如p(y(i)∣x(i),θ)p(y^{(i)}\mid x^{(i)},\theta)p(y(i)∣x(i),θ),此时就应该读作“在θ\thetaθ与x(i)x^{(i)}x(i)同时发生的情况下,y(i)y^{(i)}y(i)发生的概率是多少”
第三,假设ϵ(i)\epsilon^{(i)}ϵ(i)是独立同分布(Independent and identical distribution)随机变量,这样我们就可以引出似然函数(Likelihood function):
\begin{align}
L(\theta) & = L(\theta;X,\vec{y})\
& = p(\vec{y}\mid X;\theta) \
& =\prod_{i=1}^m p(y^{(i)}\mid x^{(i)};\theta)\
& = \prod_{i=1}m\frac{1}{\sqrt{2\pi}\sigma}\exp(-\frac{(y{(i)}-\thetaTx{(i)})2}{2\sigma2})\tag{7}
\end{align}
y⃗\vec{y}y是实际值集合的向量,XXX是特征的集合,m是样本的数量:
y⃗=[y(1)y(2)⋯y(m)](8)\vec{y}=\left[ \begin{matrix}y^{(1)}&y^{(2)}&\cdots&y^{(m)}\end{matrix}\right]\tag{8}y=[y(1)y(2)⋯y(m)](8)
X=[x0(1)x1(1)⋯xn(1)x0(2)x1(2)⋯xn(2)⋮⋮⋱⋮x0(m)x1(m)⋯xn(m)]=[—(x(1))T——(x(2))T—⋮—(x(m))T—](9)X=\left[ \begin{matrix}x_0^{(1)}&x_1^{(1)}&\cdots&x_n^{(1)}\\x_0^{(2)}&x_1^{(2)}&\cdots&x_n^{(2)}\\\vdots&\vdots&\ddots&\vdots\\x_0^{(m)}&x_1^{(m)}&\cdots&x_n^{(m)}\end{matrix}\right]=\left[ \begin{matrix}—(x^{(1)})^T—\\—(x^{(2)})^T—\\\vdots\\—(x^{(m)})T—\end{matrix}\right]\tag{9}X=x0(1)x0(2)⋮x0(m)x1(1)x1(2)⋮x1(m)⋯⋯⋱⋯xn(1)xn(2)⋮xn(m)=—(x(1))T——(x(2))T—⋮—(x(m))T—(9)
我们回过头来看L(θ)=L(θ;X,y⃗)L(\theta) = L(\theta;X,\vec{y})L(θ)=L(θ;X,y)这个式子的含义,我读作“XXX和y⃗\vec{y}y被观测到时,参数为θ\thetaθ的概率是多少“,它与$ p(\vec{y}\mid X;\theta)所描述的是同一件事情,只是换了个角度来描述而已。在这个表达下,所描述的是同一件事情,只是换了个角度来描述而已。
在这个表达下,所描述的是同一件事情,只是换了个角度来描述而已。在这个表达下,X和和和\vec{y}都是常量,表示被观测到的值,比如一开始就提供的房屋面积与相应价格的数据表,同时参数都是常量,表示被观测到的值,比如一开始就提供的房屋面积与相应价格的数据表,同时参数都是常量,表示被观测到的值,比如一开始就提供的房屋面积与相应价格的数据表,同时参数\theta成为了变量,即在观测到成为了变量,即在观测到成为了变量,即在观测到X和和和\vec{y}的前提下,每一个的前提下,每一个的前提下,每一个\theta$都有一个相对应的概率值。
举个例子:有两个箱子,甲箱里有99个白球1个黑球,乙箱里有1个白球99个黑球。参数θ\thetaθ是箱子里黑白球的配置比例,XXX是摸球这件事,y⃗\vec{y}y是指摸到白球。
概率的问法是:从箱子里摸出白球的概率是多少?
答:$ p(\vec{y}\mid X;\theta_甲)=0.99,,,p(\vec{y}\mid X;\theta_乙)=0.01$
用似然性的问法是:摸出来一个白球,猜是从哪个箱子里摸出来的?
答:甲箱。为什么?因为甲箱里摸到白球的概率更大。
现在,在XXX和y⃗\vec{y}y被观测到的情况下(摸到白球),猜是从哪条由θ\thetaθ画出来的直线上观测到的(哪个箱子)?选择使XXX和y⃗\vec{y}y出现的概率最大的那个θ\thetaθ,这也是最大似然估计(Maximum likelihood)的原则。
所以接下来的工作是对L(θ)L(\theta)L(θ)求最大值,同时为了计算方便,我们选择最大化单调递增的logL(θ)logL(\theta)logL(θ):
\begin{align}
l(\theta)&=\log L(\theta) \& = \log\prod_{i=1}m\frac{1}{\sqrt{2\pi}\sigma}\exp(-\frac{(y{(i)}-\thetaTx{(i)})2}{2\sigma2})\
& =\sum_{i=1}m\log\frac{1}{\sqrt{2\pi}\sigma}\exp(-\frac{(y{(i)}-\thetaTx{(i)})2}{2\sigma2}) \
& =\sum_{i=1}m\left[\log\frac{1}{\sqrt{2\pi}\sigma}+\log\exp(-\frac{(y{(i)}-\thetaTx{(i)})2}{2\sigma2})\right] \
& = m\log\frac{1}{\sqrt{2\pi}\sigma}-\frac{1}{\sigma2}\cdot\frac{1}{2}\sum_{i=1}m(y{(i)}-\thetaTx{(i)})2\tag{10}
\end{align}
因为(10)式左边是常数,右边有负号的关系,最大化l(θ)l(\theta)l(θ)就是最小化12(y(i)−θTx(i))2\frac{1}{2}(y^{(i)}-\theta^Tx^{(i)})^221(y(i)−θTx(i))2,也就是第一节中提到的minJ(θ)minJ(\theta)minJ(θ)。当J(θ)J(\theta)J(θ)最小时,L(θ)L(\theta)L(θ)的概率最大,即我们观测到的数据出现的概率最大,取这个时候的参数θ\thetaθ,即可画出一条最合理的直线用于预测。
下面提供两种求解最小二乘参数θ\thetaθ的方法,一种是梯度下降,另一种是最小二乘的正规方程。
1.3、最小二乘解法一——梯度下降(Gradient descent)
对于梯度下降,Ng在课上给了一个比喻:想象你正站在一个山坡上,你环顾四周,找到一个坡度最陡的方向,往那个方向走一步,然后再往坡度最陡的方向走出相同长度的一步,当你用同样的方式走了很多步的时候,你最终会到达一个最低点。
这里有两个地方需要注意,第一,这个方法让你走到的不一定是整座山的最低点,可能是某个山洼,即局部最小,对初始值敏感;第二,因为步长的关系,你可能会在最低点附近徘徊,注意调整步长。
所以梯度下降是一个需要进行迭代的,求局部最小值的一个算法。下面直接给出迭代规则:
θj:=θj−α∂∂θjJ(θ)(11)\theta_j:=\theta_j-\alpha\frac{\partial}{\partial\theta_j}J(\theta)\tag{11}θj:=θj−α∂θj∂J(θ)(11)
把J(θ)J(\theta)J(θ)带入并求θj\theta_jθj的偏导数即可得到更新规则,:=:=:=在这里表示赋值,当只有一个样本的时候,更新规则如下:
\begin{align}
\frac{\partial}{\partial\theta_j}J(\theta)&=\frac{\partial}{\partial\theta_j}\frac{1}{2}(h_\theta{(x)}-y)^2 \
& = 2\cdot\frac{1}{2}(h_\theta{(x)}-y)\cdot\frac{\partial}{\partial\theta_j}(h_\theta{(x)}-y)\
& = (h_\theta{(x)}-y)\cdot\frac{\partial}{\partial\theta_j}(\sum_{i=0}^n \theta_ix_i-y)\
& = (h_\theta{(x)}-y)\cdot\frac{\partial}{\partial\theta_j}(\theta_0x_0+\theta_1x_1+\cdots+\theta_jx_j+\cdots-y)\
& = (h_\theta{(x)}-y)x_j\tag{12}
\end{align}
带入(11),得到:
θj:=θj−α(hθ(x(i))−y(i))xj=θj+α(y(i)−hθ(x(i)))xj(13)\theta_j:=\theta_j-\alpha(h_\theta{(x^{(i)})}-y^{(i)})x_j=\theta_j+\alpha(y^{(i)}-h_\theta{(x^{(i)})})x_j\tag{13}θj:=θj−α(hθ(x(i))−y(i))xj=θj+α(y(i)−hθ(x(i)))xj(13)
每做一次求导都表示找到了当前位置梯度最陡的方向,α\alphaα代表学习速率,即步长,下山的时候每一步跨出的长度。下面提供两种具体的实施方式。
一种叫做批量梯度下降(Batch gradient descent),重复如下操作直至收敛:
θj:=θj+α∑i=1m(y(i)−hθ(x(i)))xj(14)\theta_j:=\theta_j+\alpha\sum_{i=1}^m(y^{(i)}-h_\theta{(x^{(i)})})x_j\tag{14}θj:=θj+αi=1∑m(y(i)−hθ(x(i)))xj(14)
每一次迭代都需要把所有的样本(共m个)取到,把每一个样本的真实值减去预测值,再把这个差值乘上该样本的第i个特征。重复迭代p次之后收敛了(当两次迭代的值几乎不发生变化时,可判断收敛),获得了参数之一(如θ0\theta_0θ0),然后又重复如上步骤获得θ1\theta_1θ1、θ2\theta_2θ2、……θn\theta_nθn,最后取得最终的参数向量θ\thetaθ,去画出我们的预测直线。
还有一种方法叫做随机梯度下降(Stochastic gradient descent),不断地随机取值,使其一直循环:
θj:=θj+α(y(i)−hθ(x(i)))xj(15)\theta_j:=\theta_j+\alpha(y^{(i)}-h_\theta{(x^{(i)})})x_j\tag{15}θj:=θj+α(y(i)−hθ(x(i)))xj(15)
没错这跟公式(13)是同个公式。随机梯度下降会很快接近最小值但是很难真正收敛,会在最小值附近震荡,它在实际应用中效果不错,在数据量很大的时候,比起批量梯度下降,更喜欢用到它。
1.4、最小二乘解法二——正规方程组(The normal equations)
上面提到梯度下降是一个需要进行迭代的,求局部最小值的一个算法。
相反的,正规方程组是一个不需要进行迭代的,利用观测到的全部数据直接计算出全局最小值的方法。
由公式(8)和公式(9)(这是被观测到的所有数据)可以得到:
\begin{align}
X^T\theta-\vec{y}&=\left[ \begin{matrix}(x{(1)})T\theta\\vdots\(x{(m)})T\theta\end{matrix}\right]-\left[ \begin{matrix}y{(1)}\\vdots\y{(m)}\end{matrix}\right]\
& = \left[ \begin{matrix}(x{(1)})T\theta-y{(1)}\\vdots\(x{(m)})T\theta-y{(m)}\end{matrix}\right]\tag{16}
\end{align}
这里有一个公理,对于向量zzz,有zTz=∑izi2z^Tz=\sum_iz_i^2zTz=∑izi2,又因为上式是一个向量,有:
\begin{align}
\frac{1}{2}(X\theta-\vec{y})T(X\theta-\vec{y})=\frac{1}{2}\sum_{i=1}m(h_\theta{(x{(i)})}-y{(i)})^2&=J{(\theta)}\tag{17}
\end{align}
由上面可知J(θ)J(\theta)J(θ)是一个常量,同时∇θ\nabla_\theta∇θ表示对θ\thetaθ求偏导数,trA=∑i=1nAiitrA=\sum_{i=1}^nA_{ii}trA=∑i=1nAii,下面开始对正规方程组求解:
\begin{align}
\nabla_\theta J{(\theta)}&=\nabla_\theta trJ{(\theta)}\
&=\nabla_\theta tr\frac{1}{2}(X\theta-\vec{y})^T(X\theta-\vec{y})\
&=\nabla_\theta tr\frac{1}{2}(\thetaTXT-\vec{y}^T)(X\theta-\vec{y})\
&=\frac{1}{2}\nabla_\theta tr(\thetaTXTX\theta-\thetaTXT\vec y-\vec yX\theta+\vec y^T\vec y)\
&=\frac{1}{2}\nabla_\theta tr\left[(\thetaTXTX\theta-(\thetaTXT\vec y)^T-\vec yX\theta\right]\
&=\frac{1}{2}\nabla_\theta tr(\thetaTXTX\theta)-\nabla_\theta tr\vec yX\theta\
&=\frac{1}{2}(XTX\theta+XTX\theta)-X^T\vec y\
&=XTX\theta-XT\vec y
\tag{18}\end{align}
第一行tra=atra=atra=a;第二行带入;第三行求括号内转置;第四行打开;第五行去掉与θ\thetaθ无关的y⃗Ty⃗\vec y^T\vec yyTy,并因为θTXTy⃗\theta^TX^T\vec yθTXTy是常数,有aT=aa^T=aaT=a;第六行把12∇θtr\frac{1}{2}\nabla_\theta tr21∇θtr乘进去;第七行左边多项式利用公∇ATtrABATC=BTATCT+BATC\nabla_{A^T}trABA^TC=B^TA^TC^T+BA^TC∇ATtrABATC=BTATCT+BATC,并令AT=θA^T=\thetaAT=θ,B=BT=XTXB=B^T=X^TXB=BT=XTX,C=IC=IC=I;第七行右边多项式利用trAB=trBAtrAB=trBAtrAB=trBA与∇AtrAB=BT\nabla_AtrAB=B^T∇AtrAB=BT;第八步结束。
得到结果后另倒数为0,有:
XTXθ=XTy⃗⟹θ=(XTX)−1XTy⃗(19)X^TX\theta=X^T\vec y\Longrightarrow\theta=(X^TX)^{-1}X^T\vec y\tag{19}XTXθ=XTy⟹θ=(XTX)−1XTy(19)
这样,我们就利用所有的数据计算出了参数θ\thetaθ。这就是最小二乘正规方程组的解法。
顺便给出一个拟合完成的图:
局部加权线性回归选择和预测点相近的点来做线性回归,通过加权的方式来忽略远处的点对预测的影响,一点一点地画直线的方式,近似地达到右图的效果。它的目标函数是加权的最小二乘:
J(θ)=12∑i=1mw(i)(hθ(x(i))−y(i))2(20)J{(\theta)}=\frac{1}{2}\sum_{i=1}^mw^{(i)}(h_\theta{(x^{(i)})}-y^{(i)})^2\tag{20}J(θ)=21i=1∑mw(i)(hθ(x(i))−y(i))2(20)
权值w(i)w^{(i)}w(i)要达到的效果是,离预测点越近的点权重越大,离预测点越远的点权重越小。下面是一个常用的函数:
w(i)=exp(−(x(i)−x)22τ2)(21)w^{(i)}=\exp (-\frac{(x^{(i)}-x)^2}{2\tau^2})\tag{21}w(i)=exp(−2τ2(x(i)−x)2)(21)
在这里,它是这么一个指数函数:w(i)=e−tw^{(i)}=e^{-t}w(i)=e−t,其中t=(x(i)−x)2t=(x^{(i)}-x)^2t=(x(i)−x)2,该指数函数图像如下图:
x(i)x^{(i)}x(i)是样本点,xxx是预测点,这两个点离得越近,t=(x(i)−x)2t =(x^{(i)}-x)^2t=(x(i)−x)2越小,www的值就越大(因为t=(x(i)−x)2≥0t =(x^{(i)}-x)^2\geq0t=(x(i)−x)2≥0,故www满权重为1);两个点离得越远,t=(x(i)−x)2t =(x^{(i)}-x)^2t=(x(i)−x)2越大,www的值就会越逼近0,那么这个样本点在拟合直线时几乎就不起作用,乃至被忽略。
原式中τ\tauτ是衰减速率,τ\tauτ越大衰减越慢,τ\tauτ越小衰减越快。式(21)跟正态分布长的很像,但是两者之间没有关系,只要能达到越远的点越被忽略这个效果,www也可以是另外的式子。
整个算法有点微分的意思,x轴上的每个点都根据附近的样本点画出了一条尽可能小的曲线,当这些曲线连接起来的时候,就达到了上面右图中的效果。
但是相应的,每次预测时都需要调用所有的样本来拟合那段曲线,否则它是不知道哪些点会对它的拟合产生实质性影响的,所以在数据量巨大的时候,它的代价会非常高。
2、logistic回归(Logistic regression)
2.1、logistic回归模型
logistic回归也称为逻辑回归,与线性回归这样输出是连续的、具体的值(如具体房价123万元)不同,逻辑回归的输出是0~1之间的概率,但可以把它理解成回答“是”或者“否”(即离散的二分类)的问题。回答“是”可以用标签“1”表示,回答“否”可以用标签“0”表示。
比如,逻辑回归的输出是“某人生病的概率是多少”,我们可以进一步理解成“某人是否生病了”。设置一个阈值如0.5,如果输出的是“某人生病的概率是0.2”,那么我们可以判断“此人没有生病”(贴上标签“0”)。
下面直接给出其假设:
hθ(x)=g(θTx)=11+e−θTx(22)h_\theta(x)=g(\theta^Tx)=\frac{1}{1+e^{-\theta^Tx}}\tag{22}hθ(x)=g(θTx)=1+e−θTx1(22)
其中:
g(z)=11+e−z(23)g(z)=\frac{1}{1+e^{-z}}\tag{23}g(z)=1+e−z1(23)
这个函数被称为逻辑函数(Logistic function)或者是S形函数(Sigmoid function),它是由伯努利分布通过广义线性模型求解得到的,不是凭空捏造的。下面是它的函数图像:
下面直接给出一个样本下求导后的结果:
∇θl(θ)=(y−hθ(x))xj(28) \nabla_\theta l(\theta)=(y-h_\theta(x))x_j\tag{28} ∇θl(θ)=(y−hθ(x))xj(28)
相应的迭代规则就是:
θj:=θj+α(y(i)−hθ(x(i)))xj(29) \theta_j:=\theta_j+\alpha(y^{(i)}-h_\theta(x^{(i)}))x_j\tag{29} θj:=θj+α(y(i)−hθ(x(i)))xj(29)
没错它跟式(13)长得一模一样,但是要注意的是,因为hθ(x)h_\theta(x)hθ(x)不一样,所以实际上它们是两个不同的算法,只是迭代规则一样而已,这几乎是一种通用的规则。同样的,我想在真正训练的时候,也是可以使用批量或随机梯度法的。
2.4、感知器算法(The perceptron learning algorithm)
逻辑回归通过输出0~1的概率值,拐了个弯来做分类,感知器算法则直接输出0或1的离散值,下面是其假设:
Expected node of symbol group type, but got node of type cr
通过式(29)同样的迭代,即可得到感知器算法。感知器算法是神经网络的基础,它会被作为分析的起点。
2.5、牛顿方法(Newton’s method)
在这里提出牛顿方法是因为它跟梯度法一样,可以用来求解最大化似然函数l(θ)l(\theta)l(θ)的问题,而且往往牛顿法的收敛速度更快。
牛顿方法本身是一个求解多项式实根的迭代法,如f(x)=4x3+2x2+x+1f(x)=4x^3+2x^2+x+1f(x)=4x3+2x2+x+1是一个难解的高次多项式,如果要求其f(x)=0f(x)=0f(x)=0的实根,牛顿法是一个不错的选择。下面配图以说明牛顿法的思想:
下面推导其迭代规则:令Δ=θ(0)−θ(1)\Delta=\theta^{(0)}-\theta^{(1)}Δ=θ(0)−θ(1),我们可以得到f(θ)f(\theta)f(θ)在θ(0)\theta^{(0)}θ(0)处切线的斜率为f′(θ(0))=f(θ(0))Δf{'}(\theta^{(0)})=\frac{f(\theta^{(0)})}{\Delta}f′(θ(0))=Δf(θ(0)),由此可以得到Δ=f(θ(0))f′(θ(0))\Delta=\frac{f(\theta^{(0)})}{f{'}(\theta^{(0)})}Δ=f′(θ(0))f(θ(0)),从而有θ(1)=θ(0)−Δ=θ(0)−f(θ(0))f′(θ(0))\theta^{(1)}=\theta^{(0)}-\Delta=\theta^{(0)}-\frac{f(\theta^{(0)})}{f{'}(\theta^{(0)})}θ(1)=θ(0)−Δ=θ(0)−f′(θ(0))f(θ(0)),化为更一般的情况我们就有了迭代规则:
θ(n+1)=θ(n)−f(θ(n))f’(θ(n))(31)\theta^{(n+1)}=\theta^{(n)}-\frac{f(\theta^{(n)})}{f{’}(\theta^{(n)})}\tag{31}θ(n+1)=θ(n)−f’(θ(n))f(θ(n))(31)
那么我们如何把牛顿法应用在最大化l(θ)l(\theta)l(θ)的问题上呢?
上面我们说到,牛顿法是用来求解f(θ)=0f(\theta)=0f(θ)=0的解的,相对应的,我们原来是如何使l(θ)l(\theta)l(θ)最大化的?
之前的方法都是对l(θ)l(\theta)l(θ)求导并使其导数为0,即“l’(θ)=0l’(\theta)=0l’(θ)=0”。
所以我们只要令f(θ)=l’(θ)f(\theta)=l’(\theta)f(θ)=l’(θ),并求解出θ\thetaθ,即可使l’(θ)l’(\theta)l’(θ)最大化。此时的迭代规则为:
θ(n+1)=θ(n)−l′(θ)(n)l′′(θ)(n)(32)\theta^{(n+1)}=\theta^{(n)}-\frac{l'(\theta)^{(n)}}{l''(\theta)^{(n)}}\tag{32}θ(n+1)=θ(n)−l′′(θ)(n)l′(θ)(n)(32)
上面是θ\thetaθ是实数的情况,当θ\thetaθ是向量时,迭代规则如下:
θ(n+1)=θ(n)−H−1∇θl(θ)(33)\theta^{(n+1)}=\theta^{(n)}-H^{-1}\nabla_\theta l(\theta) \tag{33}θ(n+1)=θ(n)−H−1∇θl(θ)(33)
其中HHH是函数l(θ)l(\theta)l(θ)的二次导数矩阵,被称为HessianHessianHessian矩阵,矩阵中每一个元素的值为:
Hij=∂2l(θ)∂θi∂θj(34)H_{ij}=\frac{\partial^2l(\theta)}{\partial\theta_i\partial\theta_j} \tag{34}Hij=∂θi∂θj∂2l(θ)(34)
HHH是一个n∗nn*nn∗n矩阵,nnn是向量参数θ\thetaθ的长度。
牛顿法相比起梯度往往法收敛速度更快,特别是迭代值距离收敛值比较近的时候,每次迭代都能使误差变成原来的平方,但是在高维时矩阵HHH的逆计算会非常耗时。
3、广义线性模型(Generalized Linear Models)
3.1、构造广义线性模型
上面我们已经看到了一个线性回归的例子,和一个逻辑回归分类的例子,这两个例子都是由广义线性模型推导出来的。接下来我们就从广义线性模型的假设开始,把这两个模型推导出来。下面是GLM的三个假设:
- y∣x;θy\mid x;\thetay∣x;θ~ExponentialFamily(η)ExponentialFamily(\eta)ExponentialFamily(η):固定参数θ\thetaθ,在给定xxx的情况下,yyy服从指数分布族(The exponential family)中以η\etaη为参数的某个分布。
- 给定一个xxx,我们需要的目标函数为hθ(x)=E[T(y)∣x;θ]h_\theta(x)=E\left[T(y)\mid x;\theta\right]hθ(x)=E[T(y)∣x;θ],后者为该分布的期望。
- 令η=θTx\eta=\theta^Txη=θTx。
其中,所有可以表示成如下形式的概率分布,都属于指数分布族:
p(y;η)=b(y)exp(ηTT(y)−a(η))(35)p(y;\eta)=b(y)\exp(\eta^TT(y)-a(\eta)) \tag{35}p(y;η)=b(y)exp(ηTT(y)−a(η))(35)
在式(35)中,η\etaη成为分布的自然参数(Nature parameter);T(y)T(y)T(y)是充分统计量(Sufficient statistic),通常情况下T(y)=yT(y)=yT(y)=y。当参数a、b、T都固定的时候,久定义了一个以η\etaη为参数的函数族。
3.2、线性回归模型在广义线性模型下的推导
在上文中我们是通过对高斯分布的概率解释来获得确信的J(θ)J(\theta)J(θ)的,那我们就尝试着将高斯分布放在广义线性模型下推导,看看能做出来什么模型。
首先我们看向第1点,将高斯分布表示成式(35)的形式:
\begin{align}
p(y;\mu)&=\frac{1}{\sqrt{2\pi}\sigma}\exp(-\frac{(y-\mu)2}{2\sigma2})\
&=\frac{1}{\sqrt{2\pi}\sigma}\exp(-\frac{1}{2}y^2)\cdot \exp(\mu y-\frac{1}{2}\mu ^2)
\tag{36}\end{align}
其中:
\begin{align}
\eta&=\eta^T=\mu\
T(y)&=y\
a(\eta)&=\frac{1}{2}\mu ^2=\frac{1}{2}\eta ^2\
b(y)&=\frac{1}{\sqrt{2\pi}\sigma}\exp(-\frac{1}{2}y^2)
\tag{37}\end{align}
这就正好表示成了式(35)的形式,符合了第1点假设的要求,接下来我们看向第2点,使我们的假设函数等于高斯分布的期望:
\begin{align}
h_\theta(x)&=E\left[T(y)\mid x;\theta\right]\
&=E\left[y\mid x;\theta\right]\
&=\mu
\tag{38}\end{align}
最后我们看向第3点,令η=θTx\eta=\theta^Txη=θTx,又因为我们在式(37)中有η=μ\eta=\muη=μ,则:
\begin{align}
h_\theta(x)&=\mu=\eta=\theta^Tx
\tag{39}\end{align}
这就是最开始的式(2)的由来,惊不惊喜,刺不刺激。
3.3、逻辑回归模型在广义线性模型下的推导
逻辑回归的假设里是用到了伯努利分布的,我们同样地将伯努利分布放在广义线性模型下进行推导,看看是不是能得到逻辑回归的hθ(x)h_\theta(x)hθ(x),如式(22)。
还是从第一点看起,把伯努利分布表示成指数分布族的形式:
\begin{align}
p(y;\phi)&=\phiy(1-\phi){1-y}\
&=(e{\log\phi})y(e{\log(1-\phi)}){1-y}\
&=\exp(y\log\phi+(1-y)\log(1-\phi))\
&=\exp((\log\frac{\phi}{1-\phi})y+\log(1-\phi))
\tag{40}\end{align}
相应的参数为:
\begin{align}
\eta&=\log\frac{\phi}{1-\phi}\Longrightarrow\phi=\frac{1}{1+e^{-\eta}}\
T(y)&=y\
a(\eta)&=-\log(1-\phi)=\log((1-\phi){-1})=\log(1+e\eta)\
b(y)&=1
\tag{41}\end{align}
接下来执行第二点与第三点,令假设函数hθ(x)h_\theta(x)hθ(x)等于伯努利函数的期望,并且令η=θTx\eta=\theta^Txη=θTx:
\begin{align}
h_\theta(x)&=E\left[T(y)\mid x;\theta\right]\
&=E\left[y\mid x;\theta\right]\
&=\phi\
&=\frac{1}{1+e^{-\eta}}\
&=\frac{1}{1+e{-\thetaTx}}
\tag{42}\end{align}
这就是逻辑回归的假设函数式(22)的由来了。
除了高斯分布与伯努利分布,大多数的概率分布都能表示成指数分布族的形式,如多项式分布(Multinomial),对有K个离散结果的事件建模;泊松分布(Poisson),对计数过程进行建模,如网站访问量的计数问题;指数分布(Exponential),对有间隔的证书进行建模,如预测公交车的到站时间的问题;等等,然后通过进一步的推导,就能得到各自的线性模型,这大大扩展了线性模型可解决问题的范围。
3.4、多分类回归模型(Softmax regression)在广义线性模型下的推导
3.4.1、多项式分布
接下来再来利用多项式分布进行建模,得到解决多分类问题的模型,它是逻辑回归二分类模型的扩展。
首先我们给出多项式分布的一些公式:
\begin{align}
p(y=i)&=\phi_i\
\sum_{i=1}^k\phi_i&=1
\tag{43}\end{align}
其中$y \in KaTeX parse error: Expected '}', got 'EOF' at end of input: { 1,2,3,\cdots,k KaTeX parse error: Expected 'EOF', got '}' at position 1: }̲,又因为\phi_i累加为1,我们可以只保留累加为1,我们可以只保留累加为1,我们可以只保留k-1$个参数:
ϕk=1−∑i=1k−1ϕi(44)\phi_k=1-\sum_{i=1}^{k-1}\phi_i\tag{44}ϕk=1−i=1∑k−1ϕi(44)
3.4.2、从两点分布出发
为了有个直观的理解,我们先来看k=2k=2k=2时的情况。
此时p(y=1)=ϕ1p(y=1)=\phi_1p(y=1)=ϕ1,p(y=2)=ϕ2p(y=2)=\phi_2p(y=2)=ϕ2,并且ϕ1+ϕ2=1\phi_1+\phi_2=1ϕ1+ϕ2=1,于是有ϕ2=1−ϕ1\phi_2=1-\phi_1ϕ2=1−ϕ1
\begin{align}
p(y=1)&=\phi_1\
p(y=2)&=\phi_2\
\phi_1+\phi_2&=1 \Longrightarrow \phi_2=1-\phi_1\
p(y\mid \phi_1)&=\phi_1{2-y}\phi_2{y-1}=\phi_1{2-y}(1-\phi_1){y-1}\
&\Longrightarrow p(y\mid \phi)=\phi{2-y}(1-\phi){y-1}
\tag{45}\end{align}
当y=1y=1y=1时,只取ϕ1\phi_1ϕ1;当y=2y=2y=2时,只取ϕ2\phi_2ϕ2,即1−ϕ11-\phi_11−ϕ1,又因为只有一个参数ϕ1\phi_1ϕ1了,干脆就直接令ϕ1=ϕ\phi_1=\phiϕ1=ϕ。这与两点分布所要表达的意义是一模一样的,只不过两点分布只做二分类,就采取了取0和1的方式,使得整个式子更简洁(参考式(40)的第一行)。
3.4.3、多项式概率分布
那么在更一般的多分类问题中,其概率分布为:
p(y;ϕ)=ϕ1I(y=1)ϕ2I(y=2)⋯ϕkI(y=k)(46)p(y;\phi)=\phi_1^{I(y=1)}\phi_2^{I(y=2)}\cdots\phi_k^{I(y=k)}\tag{46}p(y;ϕ)=ϕ1I(y=1)ϕ2I(y=2)⋯ϕkI(y=k)(46)
这个式子所要表达的意思是,当yyy的值取到iii的时候,只取相应的ϕi\phi_iϕi的值。因为当kkk的取值不同之时,ϕi\phi_iϕi的指数表示形式会不一样(如上面的例子,同是两点分布,仅是取值不同就会使得其指数表达不一样)。
为了更方便的表达,引入了上面的指示器函数(Indicator function),如果判断条件为真,输出1,反之输出0:
\begin{align}
I\begin{Bmatrix} Ture \end{Bmatrix}&=1\
I\begin{Bmatrix} False \end{Bmatrix}&=0
\tag{47}\end{align}
如I{3=3}=1I\begin{Bmatrix} 3=3 \end{Bmatrix}=1I{3=3}=1,I{3=2}=0I\begin{Bmatrix} 3=2 \end{Bmatrix}=0I{3=2}=0。
所以,当y=1y=1y=1时,式(47)的第一个乘子ϕ1\phi_1ϕ1的指数为1,其他所有的乘子的指数全部为0,进一步推广到所有取值的情况,就起到了当y的值取到i的时候,只取相应的ϕi的值的作用。
但是这个指示器函数是一个抽象的东西,如何才能确实地实现这个作用呢?我们给出如下所示的T(y)T(y)T(y):
\begin{align}
T(1)=\left[ \begin{matrix}1\0\0\\vdots\0\end{matrix}\right],
T(2)=\left[ \begin{matrix}0\1\0\\vdots\0\end{matrix}\right],\cdots,
T(k-1)=\left[ \begin{matrix}0\0\0\\vdots\1\end{matrix}\right],
T(k)=\left[ \begin{matrix}0\0\0\\vdots\0\end{matrix}\right]
\tag{48}\end{align}
这样T(y)T(y)T(y)中的某个元素可以表示成:
\begin{align}
T(y)_i=I\begin{Bmatrix} y=i \end{Bmatrix}
\tag{49}\end{align}
举例来说,当y=2y=2y=2时,T(2)T(2)T(2)中的第一个元素为:T(2)1=I{2=1}=0T(2)_1=I\begin{Bmatrix} 2=1 \end{Bmatrix}=0T(2)1=I{2=1}=0;T(2)T(2)T(2)中的第二个元素为:T(2)2=I{2=2}=1T(2)_2=I\begin{Bmatrix} 2=2 \end{Bmatrix}=1T(2)2=I{2=2}=1,这就解决了指数上标什么时候取0什么时候取1的问题。
于是我们从式(46)中进一步得到具体的多项式分布的概率分布:
\begin{align}
p(y;\phi)&=\phi_1^{I\begin{Bmatrix} y=1\end{Bmatrix}}\phi_2^{I\begin{Bmatrix} y=2\end{Bmatrix}}\cdots\phi_k^{I\begin{Bmatrix} y=k \end{Bmatrix}}\
&=\phi_1^{I\begin{Bmatrix} y=1\end{Bmatrix}}\phi_2^{I\begin{Bmatrix} y=2\end{Bmatrix}}\cdots\phi_k{1-\sum_{i=1}{k-1}I\begin{Bmatrix} y=i\end{Bmatrix}}\
&=\phi_1{T(y)_1}\phi_2{T(y)_2}\cdots\phi_k{1-\sum_{i=1}{k-1}T(y)_i}
\tag{50}\end{align}
3.4.4、在广义线性模型下的推导
还是从广义线性模型的第一点开始做起,将式(50)的概率分布表示成指数分布族的形式:
\begin{align}
p(y;\phi)&=\phi_1{T(y)_1}\phi_2{T(y)2}\cdots\phi_k{1-\sum_{i=1}{k-1}T(y)i}\
&=\exp(T(y)1\log\phi_1+T(y)2\log\phi_2+\cdots+((1-\sum{i=1}^{k-1}T(y)i)\log\phi_k)\
&=\exp(T(y)1\log\frac{\phi_1}{\phi_k}+T(y)2\log\frac{\phi_2}{\phi_k}+\cdots+T(y){k-1}\log\frac{\phi{k-1}}{\phi_k}+\log\phi_k)\
&=b(y)\exp(\eta^TT(y)-a(\eta))
\tag{51}\end{align}
我们有:
\begin{align}
\eta&=\left[ \begin{matrix}\log\frac{\phi_1}{\phi_k}\\log\frac{\phi_2}{\phi_k}\\vdots\\log\frac{\phi{k-1}}{\phi_k}\end{matrix}\right]\
T(y)&=式(48)\
a(\eta)&=-\log\phi_k\
b(y)&=1
\tag{52}\end{align}
观察式(52)我们可以发现,向量η\etaη的第iii个元素为:
ηi=logϕiϕk⟹eηi=ϕiϕk(53)\eta_i=\log\frac{\phi_i}{\phi_k}\Longrightarrow e^{\eta_i}=\frac{\phi_i}{\phi_k}\tag{53}ηi=logϕkϕi⟹eηi=ϕkϕi(53)
令其累加,得到:
∑i=1keηi=∑i=1kϕiϕk=1ϕk⟹ϕk=1∑i=1keηi(54)\sum_{i=1}^ke^{\eta_i}=\frac{\sum_{i=1}^k\phi_i}{\phi_k}=\frac{1}{\phi_k}\Longrightarrow \phi_k=\frac{1}{\sum_{i=1}^ke^{\eta_i}}\tag{54}i=1∑keηi=ϕk∑i=1kϕi=ϕk1⟹ϕk=∑i=1keηi1(54)
把ϕk\phi_kϕk带入式(53),即可得到多项式分布的期望:
ϕi=eηi∑i=1keηi(55)\phi_i=\frac{e^{\eta_i}}{\sum_{i=1}^ke^{\eta_i}}\tag{55}ϕi=∑i=1keηieηi(55)
那么接下来的第二步与第三步,令假设函数hθ(x)h_\theta(x)hθ(x)等于多项式分布的期望,并且令η=θTx\eta=\theta^Txη=θTx,就显得非常地顺理成章了:
\begin{align}
h\theta(x)&=E\left[T(y)\mid x;\theta\right]\
&=E\left[ \begin{matrix}p(y=1\mid x;\theta)\p(y=2\mid x;\theta)\\vdots\p(y={k-1}\mid x;\theta)\end{matrix}\right]\
&=\left[ \begin{matrix}\phi_1\\phi_2\\vdots\\phi{k-1}\end{matrix}\right]\
&=\left[ \begin{matrix}\frac{e{\eta_1}}{\sum_{i=1}ke{\eta_i}}\\frac{e{\eta_2}}{\sum{i=1}ke{\eta_i}}\\vdots\\frac{e{\eta_{k-1}}}{\sum_{i=1}ke^{\eta_i}}\end{matrix}\right]
=\left[ \begin{matrix}\frac{e{\theta_1Tx}}{\sum_{i=1}ke{\theta_iTx}}\\frac{e{\theta_2Tx}}{\sum_{i=1}ke{\theta_iTx}}\\vdots\\frac{e{\theta_{k-1}Tx}}{\sum_{i=1}ke{\theta_i^Tx}}\end{matrix}\right]
\tag{56}\end{align}
这就是softmax回归的模型了。
如果要求取参数θ\thetaθ,还是和线性回归与逻辑回归一样,取其log似然函数做最大似然估计,具体方法可以选择梯度法或者牛顿法。