在R语言中实现Bayes方法对logistic的回归分析

田 薇, 李晓毅, 付志慧

(沈阳师范大学 数学与系统科学学院, 沈阳 110034)



统计学

在R语言中实现Bayes方法对logistic的回归分析

田 薇, 李晓毅, 付志慧

(沈阳师范大学 数学与系统科学学院, 沈阳 110034)

对于logistic回归分析的处理办法,一直采用的都是极大似然估计的EM算法,由于计算方法的固定及计算过程的复杂性,例如,该算法对于初值的选取要求很高,否则收敛速度很慢。Gibbs抽样法作为一种高效灵活的估计方法广泛应用于广义线性回归模型,其中Probit回归模型由于联系函数为正态分布,使得回归系数的后验分布为共轭正态,从而抽样简单快捷。而Logit模型的后验分布比较复杂,无法直接抽取。本文基于增加数据的Gibbs抽样方法,通过引入Plya-Gamma分布族的潜在变量,使得模型中的回归系数参数的满条件分布为共轭正态分布,从而回归系数的马氏链很容易构造,回归系数的估计为后验均值估计。通过一组实际数据,分别调用R语言Glm包和 BayesLogit包,并对比2种方法的估计结果,二者差别不大,表明Plya-Gamma潜变量Bayes估计法在处理logistic回归模型时的可用性、准确性。

logistic回归模型; Plya-Gamma分布; MCMC; R语言

0 引 言

logistic回归模型常用在寻找危险因素、预测、判别3个方面。由于线性回归模型只能分析连续型数据,具有很强的局限性,而实际中的logistic回归用途是极为广泛的,logistic回归几乎已经成了流行病学和医学中最常用的分析方法,logit模型在处理属性数据或分类数据方面极具优势。在估计logistic回归模型时,一般采用极大似然法。若引入合理的先验分布函数,Bayes方法对于许多模型的参数估计问题一直特别有效,对提高统计推断质量具有实际意义。另外,忽略参数的先验信息,有时是一种浪费,甚至还会导致不合理的结论。logit回归模型联系函数为logistic分布函数,回归系数的后验分布无具体形式,需要采用Metropolis-Hastings抽样法,该方法在应用过程中若建议分布选取的不合理,会导致接受概率很小,因此马氏链收敛很慢。

本文尝试使用一种由Albert提出的数据添加新方法,该方法在一定程度上区别于常见的数据添加方式。引入一个来自Plya-Gamma分布的随机变量(具体的构造方法将在文中加以说明)。为了提高估计的质量,除了当前样本数据,还可以利用客观信息和经验累积的信息,先验信息的加入,参数估计更加稳定,也更合理和符合实际。回归系数的共轭分布仍为正态分布,使得抽样及后验估计很容易得到,更加方便计算。在R语言的BayesLogit程序包中,MCMC抽样及数据处理均可得以实现。最后通过一组实际数据,利用Glm包实现传统似然估计方法,并与MCMC估计方法相比较。

定义1 随机变量X是带有参数b(b>0),和c∈R的Plya-Gamma分布,其中变量Χ分布记为X～PG(b,c), 形式为

这里gk～Ga(b,1),是独立的伽玛随机变量。

其中:Vw=(XTΩX+B-1)-1;mω=Vω(XTκ+B-1b);κ=(y1-n1/2,…,yN-nN/2);Ω为ωi的对角线矩阵。

接下来开始对上面所得到的抽样方法进行演绎证明,先从一些定理和积分公式开始。

定理p(ω)为随机变量ω的密度函数,且ω～PG(b,0),b>0。对于所有的a∈R,有下列恒等式:

其中κ=a-b/2 。

对式子(2)的非正态化联合密度处理,可得ω的条件分布

通过定理和积分,可以得到第i个观测值的似然函数为

(4)

其中p(ωi|ni,0)为带参数随机变量(ni,0),服从Plya-Gamma分布的密度函数。

在n组数据下,β的后验条件分布为

从整理得到的分布形式知,β的后验条件分布服从正态分布,即P(β|ω,y)∝N(mω,Vω)。

其中:mω=Vω(XTκ+B-1b);Vw=(XTΩX+B-1)-1。其理论依据为正态分布(方差已知)的共轭先验还是正态分布。此处:z=(κ1/ω1,…,κN/ωN);Ω=diag(ω1,…,ωN)。由式(5)可知,β的条件分布为高斯似然,且先验p(β)也服从高斯分布,因此该线性模型得以简单计算。

2 模拟研究

在二项分布族中,logistic回归模型是最重要的模型。对于响应变量Y有p个自变量(或称为解释变量),记为X1,X2,…,Xp。在p个自变量作用下出现成功的条件概率为P=P{Y=1|X1,X2,…,Xp},那么其logistic回归模型可表示为

其中:称β0为截距;称β1,β2,…,βp为logistic回归模型系数。

对上式作logit变换,logistic回归模型可以变成下列线性形式:

可以使用线性回归模型对参数进行估计,这也是logistic回归模型属于广义线性模型的原因。

当logistic回归模型的分布函数为

农村集体土地上不动产登记工作中，不动产登记权利人主体不一致的情况很常见，突出表现在以下3个方面：①规划审批手续的建房人与土地审批手续的使用人不同；②土地使用权人与登记簿中房屋所有人不同；③房屋所有权与土地所有权人的主体不一致。

其联合分布函数为

通常利用极大似然法,对该线性回归模型中的参数进行估计。

通过一组实际数据,用R语言中BayesLogit包对数据进行处理,得到spambase数据包的其他特征数,详见表1。并与原始方法得到的模型系数估计值进行比照,详见表2(logistic回归模型系数估计值表)。其结果表明利用BayesLogit方法处理logistic回归分析问题的可行性。

表1 spambase数据包其他特征数Tab.1 spambase packet number of other features

表2 logistic回归模型系数估计值表Tab.2 logistic regression coefficient estimates table

3 结 论

[ 1 ]GAMERMAN D. Sampling from the posterior distribution in generalized linear mixed models[J]. Statistics and Computing, 1997(7):57-68.

[ 2 ]HOLMAN R, GLAS C A W. Modeling non-ignorable missing data mechanisms with item response theory models[J]. BRIT J MATH STAT PSY, 2005,58(1):1-17.

[ 3 ]HAMBLETON R K. Fundamentals of item response theory[M]. NewYork:Sage Publication, 1991.

[ 4 ]RUBIN D B. Inference and missing data[J]. Biometrika, 1976,63(3):581-592.

[ 5 ]LITTLE R J A, RUBIN D B. Statistical analysis with missing data[M]. Manhattan:John Wiley&Sons, 2014.

[ 6 ]MASTERS G N.ARasch model for partial credit scoring[J]. Psychometrika, 1982,47(2):149-174.

[ 7 ]ALBERT J H. Bayesian estimation of normal ogive item response curves using Gibbs sampling[J]. J EDUCBEHAV STAT, 1992,17(3):251-269.

[ 8 ]JONES D H, NEDIAK M S. Item parameter calibration of LSAT items using MCMC approximation of Bayes posterior distribution[M]. Newtown:Law School Admission Coucil, 2005.

[ 9 ]GELMAN A, RUBIN D B. Inference from iterative simulation using multiple sequences[J]. STAT SCI, 1992:457-472.

[10]MARIS G,BECHGER T M. An introduction to the DAT Gibbs sampler for the two-parameter logistic(2PL) model and beyond[J]. International Journal of Methodology and Experimental Psychology, 2005,26(2):327-352.

[11]LUDLOW L H, O’LEARY M. Scoring omitted and not-reached items: practical data analysis implications[J]. EDUC PSYCHOL MEAS, 1999,59(4):615-630.

[12]HUISMAN M. Imputation of missing itemresponses:Some simple techniques[J]. QUAL QUANT, 2000,34(4):331-351.

[13]MURAKI E, BOCK R D. PARSCALE:IRT based test scoring and item analysis for graded open-ended exercises and performance tasks[M]. Scientific Software International, 1993.

[14]LORD F M. Maximum likelihood and Bayesian parameter estimation in item response theory[J]. J EDUC MEAS, 1986,23(2):157-162.

[15]MOUSTAKI I, KNOTT M. Weighting for item non-response in attitude scales by using latent variable models with covariates[J]. J R STAT SOC B, 2000,163(3):445-459.

Bayesian inference for logistic models in R Language

TIAN Wei, LI Xiaoyi, FU Zhihui

(College of Mathemetics and Systems Science, Shenyang Normal University, Shenyang 110034, China)

For the approach to logistic regression analysis, using a maximum likelihood estimation are the EM, due to the complexity and fixity of calculation, for example, the initial value of the algorithm is demanding, otherwise the convergence rate is slow. Gibbs sampling as an efficient and flexible estimation is widely used for generalized linear regression models, due to the contact function is normal in Probit model, so that the posterior distribution of the regression coefficients is Conjugated Normality and sampling is easier.The posterior of Logit model is complex, unable to directly extract, based on Gibbs to increase data by introducing latent variables Plya-Gamma distribution families, making the regression coefficient parameters of full conditional distribution Conjugated Normality, thereby Markov chains regression coefficient is easy to construct the estimated regression coefficients for the posterior mean estimate. Through a set of actual data, respectively, calling R language package of BayesLogit and Glm, and comparing the results of the two methods, the difference is small, indicating Plya-Gamma latent variable Bayesian estimation in dealing with the accuracy of logistic regression model.

logistic regression model; Plya-Gamma distribution; MCMC; R language

2016-04-16。

国家自然科学基金青年基金资助项目(11201313)。

田 薇(1990-),女,辽宁葫芦岛人,沈阳师范大学硕士研究生; 通信作者: 李晓毅(1956-),女,辽宁葫芦岛人,沈阳师范大学教授。

1673-5862(2016)03-0321-04

O212.8

A

10.3969/ j.issn.1673-5862.2016.03.014