广义线性模型(GLM)
一、广义线性模型概念
在讨论广义线性模型之前,先回顾一下基本线性模型,也就是线性回归。
在线性回归模型中的假设中,有两点需要提出:
(1)假设因变量服从高斯分布:$Y={{\theta }^{T}}x+\xi $,其中误差项$\xi \sim N(0,{{\sigma }^{2}})$,那么因变量$Y\sim N({{\theta }^{T}}x,{{\sigma }^{2}})$。
(2)模型预测的输出为$E[Y]$,根据$Y={{\theta }^{T}}x+\xi $,$E[Y]=E[{{\theta }^{T}}x+\xi ]={{\theta }^{T}}x$,记$\eta ={{\theta }^{T}}x$,则$\eta =E[Y]$
广义线性模型可以认为在以上两点假设做了扩展:
(1)因变量分布不一定是高斯分布,服从一个指数分布族(见下文)即可。
(2)模型预测输出仍然可以认为是$E[Y]$(实际上是$E[T(Y)]$,许多情况下$T(Y)=Y$),但是$Y$的分布不一定是高斯分布,$E[Y]$和$\eta ={{\theta }^{T}}x$也不一定是简单的相等关系,它们的关系用$\eta =g(E[Y])$描述,称为连接函数,其中$\eta $称为自然参数。
由于以上两点的扩展,广义线性模型的应用比基本线性模型广泛许多。对于广义线性这个术语,可以理解为广义体现在因变量的分布形式比较广,只要是一指数分布族即可,而线性则体现在自然参数$\eta ={{\theta }^{T}}x$是$\theta $的线性函数。
二、广义线性模型的构建
上文提到指数分布族,它是广义线性模型的基础,所以先简单了解一下指数分布族。
对于变量$y$,如果其分布可写成$p(y;\eta )=b(\eta )\exp ({{\eta }^{T}}T(y)-a(\eta ))$的形式,则称$y$服从一个指数分布族,自然参数$\eta $是分布的参数。为什么这样定义是牛逼的数学家弄的,咱就看看它在广义线性模型中怎么用的吧~
实际中的许多分布都是一个指数分布族,如高斯分布,二项分布,泊松分布,多项分布等等,所以之前写的线性回归、逻辑回归实际上都是个广义线性模型。以逻辑回归为例来看看。
逻辑回归假设$y$服从参数为$\phi $伯努利分布,$p(y)={{\phi }^{y}}{{(1-\phi )}^{1-y}}$,$E[y]=\phi $。
下面将其写出指数分布族的形式:
$\begin{align} p(y)&={{\phi }^{y}}{{(1-\phi )}^{1-y}} \\& =\exp (\log ({{\phi }^{y}}{{(1-\phi )}^{1-y}})) \\& =\exp (y\log (\phi )+(1-y)\log (1-\phi )) \\ & =\exp (y\log (\frac{\phi }{1-\phi })+\log (1-\phi )) \\\end{align}$
与指数分布族的一般形式对比可发现:
$\begin{align}& b(\eta )=1 \\& \eta =\log (\frac{\phi }{1-\phi })\Rightarrow \phi =\frac{1}{1+{{e}^{-\eta }}} \\& T(y)=y \\& a(\eta )=-\log (1-\phi ) \\\end{align}$
可见这是符合假设(1)的。根据假设(2),我们预测的是
$\begin{align}{{h}_{\theta }}(x)&=E[T(y)] \\& =E[y](=\phi ) \\& =\frac{1}{1+{{e}^{-\eta }}} \\& =\frac{1}{1+{{e}^{-{{\theta }^{T}}x}}} \\\end{align}$
这正是之前逻辑回归的模型。
可以看出以上推导过程中$E[y]=\frac{1}{1+{{e}^{-\eta }}}$这一步比较重要,起到了连接预测值$E[y]$和自然参数$\eta $的作用,这就是连接函数的作用。
回顾假设(1)(2)和以上逻辑回归推导过程,可以看出,构建一个广义线性模型需要两个步骤:
(1)确定预测变量$y$的分布是一个指数分布族
(2)确定连接函数。连接函数可以是任意的,但根据上文可以看出,一但步骤(1)中的分布形式给定,就可以推导出一个连接函数,这个根据分布推导出的分布称为标准连接函数,也是通常默认采用的。所以一般步骤(1)中分布的形式给定,步骤(2)也就默认确定了。
经过以上两个步骤,模型就建立好了,接下来就是写出似然函数,最大化似然函数估计模型参数。对于广义线性模型的参数估计,有个专门的算法IRWLS(iteratively weighted least squares),感兴趣的可以查阅相关文献。另外,关于模型的假设检验,也不写了。写不出来。。。看来几天数学课本,头炸了。。还是安静的做个程序员吧~
三、广义线性模型应用
广义线性模型的应用最广泛的的是逻辑回归和泊松回归。逻辑回归将因变量建模为伯努利分布,输出是二值的,通常用来做二分类。泊松回归将因变量的分布建模为泊松分布,一般用来预测类似顾客数目、一个时间段内给定事件发生数目的问题。
另外,对于多分类问题,将因变量建模为多项分布也是一个广义线性模型。
之前在逻辑回归中,没有提到广义线性模型,现在可以直接用R中提供的广义线性模型来拟合。
1 x<-read.table("q1x.txt");
2 x<-as.matrix(x);
3 y<-scan("q1y.txt");
4 y<-matrix(y,ncol=1);
5
6 gfit<-glm(y~x,family=binomial());
7 print(coef(gfit));
可以看出这和之前拟合的参数基本一致。
泊松回归和多项回归没弄得数据,不代码模拟了~
附:R中GLM相关
R中用glm(formula,family...)函数来做广义线性模型,并且提供了一下指数族分布:
binomial(link = "logit")
gaussian(link = "identity")
Gamma(link = "inverse")
inverse.gaussian(link = "1/mu^2")
poisson(link = "log")
quasi(link = "identity", variance = "constant")
quasibinomial(link = "logit")
quasipoisson(link = "log")
用上选项指定glm()中family参数,就可以得到不同的模型。
另外,robust包中的glmRob()函数可用来拟合稳健的广义线性模型,包括稳健Logistic回归,稳健泊松回归等。当拟合Logistic回归模型数据出现离群点和强影响点时,稳健Logistic回归便可派上用场;
多项分布回归:若响应变量包含两个以上的无序类别,便可使用mlogit包中的mlogit()函数拟合多项Logistic回归;
序数Logistic回归:若响应变量是一组有序的类别(比如,信用风险为差/良/好),便可使用rms包中的lrm()函数拟合序数Logistic回归
参考资料
[1]Andrew Ng 机器学习视频讲义:http://cs229.stanford.edu/
[2]统计之都:http://cos.name/2011/01/how-does-glm-generalize-lm-assumption/
[3]《R语言实战》 /(美)科巴科弗(Kabacoff,R.I.)著;高涛,肖楠,陈钢译.人民邮电出版社,2013.1