R语言实战广义线性模型

Posted 2020-09-02 gy_jerry

本文对应《R语言实战》第13章：广义线性模型

广义线性模型扩展了线性模型的框架，包含了非正态因变量的分析。

两种流行模型：Logistic回归（因变量为类别型）和泊松回归（因变量为计数型）

 

glm()函数的参数

分布族	默认的连接函数
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”)

 

连用的函数

函数	描述
summary()	展示拟合模型的细节
coefficients(), coef()	列出拟合模型的参数（截距项和斜率）
confint()	给出模型参数的置信区间（默认为95%）
residuals()	列出拟合模型的残差值
anova()	生成两个拟合模型的方差分析表
plot()	生成评价拟合模型的诊断图
predict()	用拟合模型对新数据集进行预测

 

模型拟合和回归诊断：

#诊断图（以下model都是glm拟合的模型）
plot(predict(model, type = “response”), residuals(model, type = “deviance”))
 
#帽子值（hat value），学生化残差，Cook距离统计量近似值
plot(hatvalues(model))
plot(rstudent(model))
plot(cooks.distance(model))

#综合性诊断图
library(car)
influencePlot(model)

　　

当响应变量有许多时，诊断图非常有用；而当响应变量只有有限个值时（比如Logistic回归），诊断图的功效就会降低很多。 

 

 

Logistic回归

一般过程：

首先将所有变量作为预测变量拟合模型，通过回归系数的显著性，筛选对方程贡献显著的变量，再次拟合模型。使用anova()函数对两个嵌套模型的拟合优度进行检验，广义线性回归时，可采用卡方检验，在卡方值不显著的时候，表示少变量的模型与多变量的模型拟合效果没有差异。（具体例子详见书上285-288页）

解释模型参数：

Logistic回归中，响应变量是Y=1的对数优势比(log)。回归系数含义是当其他预测变量不变时，一单位预测变量的变化可引起的响应变量对数优势比的变化。由于对数优势比解释性差，所以常将结果进行指数化后进行解释，也可以使用confint()函数获取系数的置信区间。

exp(confint(fit.reduced))

　　

评价预测变量对结果概率的影响： 

由于使用概率的方式思考比使用优势比更直观，所以可以先创建一个包含感兴趣预测变量值的虚拟数据集，然后对这个数据集使用predict()函数，以预测这些值的结果概率。（具体例子详见书289-290页）

过度离势：

抽样于二项分布的期望方差是，n为观测数，pie为属于Y=1组的概率。

所谓过度离势，即观测到的响应变量的方差大于期望的二项分布的方差。过度离势会导致奇异的标准误检验和不精确的显著性检验。

当出现过度离势时，仍可使用glm()函数拟合Logistic回归，但此时需要将二项分布改为类二项分布（quasibinomial distribution）。

检测过度离势的一种方法是比较二项分布模型的残差偏差与残差自由度，如果比值：

 

比1大很多，便可认为存在过度离势。

具体检验方式：拟合两次模型，第一次使用family = “binomial”, 第二次使用family = “quasibinomial”, 记第一次返回的对象为fit，第二次返回的对象为fit.od，

pchsq(summary(fit.od)$dispersion * fit$df.residual, fit$df.residual, lower = F)

　　

 提供的p值即可对零假设H0：Φ=1 与备择假设 H1：Φ≠1进行检验，若p很小就可以拒绝H0.

Logistic回归的扩展：

稳健Logistic回归：robust包中的glmRob()函数，解决离群点和强影响点的问题

多项分布回归：mlogit包中的mlogit()函数，应用于响应变量包含两个以上的无序类别

序数Logistic回归：rms包中的lrm()函数，应用于响应变量是一组有序类别

 

 

泊松回归：

适用范围：通过一系列连续型和类别型预测变量来预测计数型结果变量。

分析过程与模型参数的解释与Logistic回归类似。

过度离势：

泊松分布的方差和均值相等。当响应变量观测的方差观测的方差比依据泊松分布预测的方差大时，泊松回归可能发生过度离势。处理计数型数据时经常发生过度离势，因此需要格外注意。

可能的原因：

1. 遗漏了某个重要的预测变量；
2. 可能因为时间相关；
3. 在纵向数据分析中，重复测量的数据由于内在群聚特性可导致过度离势。

如果存在过度离势，在模型中将无法进行解释，可能会得到很小的标准误和置信区间，并且显著性检验也过于宽松（也就是说，将会发现并不真是存在的效应）。

与Logistic回归类似，如果残差偏差与残差自由度的比例远远大于1，表明存在过度离势。

检验方法：qcc包

library(qcc)
qcc.overdispersion.test(breslow.dat$sumY, type = “poisson”)

　　

p值显著（小于0.05）表明存在过度离势。 

与Logistic回归的过度离势处理类似，通过用family = “quasipoisson”替换family = “poisson”改进。需要注意的是，使用类泊松（quasi-Poisson）方法所得的参数估计与泊松方法相同，但标准误变大了许多。

 

泊松回归的扩展：

1. 时间段变化的泊松回归：

对于泊松回归的讨论，是将响应变量局限在一个固定长度时间段中进行测量，整个观测集中时间长度都是不变的。允许时间段变化时，假设结果变量是一个比率，也就是将

 

修改为：

 

或等价形式：

 

具体实现方式是，在glm()函数里增加offset参数，offset = log(time).

1. 零膨胀的泊松回归：

在数据集中0计数的数目时常比用泊松模型预测的数目多，这些0值称为结构零值。此时使用零膨胀的泊松回归分析数据，将同时拟合两个模型，可以看做是Logistic回归和泊松回归的组合。pscl包中的zeroinfl()函数可做零膨胀泊松回归。

1. 稳健泊松回归：

robust包中的glmRob()函数可以拟合稳健广义线性模型，包含稳健泊松回归。解决存在离群点和强影响点带来的问题。