Logistic模型统计分析的R软件实现.doc

Logistic模型统计分析的R软件实现2011年12月第4期连云港师范高等专科学校JournalofLianyungangTeachersCollegeDecember,2011No.4Logistic模型统计分析的R软件实现光琳(江苏联合职业技术学院连云港财经分院,江苏连云港222003)摘要随着回归分析理论的不断发展,Logistic回归模型已成为目前重要的研究课题,并在医学,工程技术等领域得到了广泛应用.它的数据计算依赖于统计软件,文章结合实例研究了如何在R软件中实现对Logistic模型的统计分析.关键词Logistic模型;R软件;统计诊断中图分类号o212.1文献标识码A文章编号10097740(2011)040106031引言Logistic模型最早由比利时生物学家Verhulst于1838年导出,该模型形式简单但在生物,经济,环境科学等领域都有着广泛应用,特别是在人口数量预测,耐用消费品销售量预测等方面有着重要的应用.它主要研究二分变量与诸多自变量的相互关系.随着理论研究的不断深入,Logistic模型也生发出多种形式,其中如下形式:,n,logit(P)=logIJ=flo+X1+X2+ppl上,X(*)在工程以及医学,流行病学,生物学上都有广泛的应用价值,13益受到人们的重视.由于Logistic回归模型涉及较复杂的数学理论,数据统计分析的计算往往较为复杂,一般无法用手工计算,L2只能在计算机上实现.目前的统计软件有很多,如SAS,SPSS,STAT,SPLUS等,但大多是商业性的,购买需支付昂贵的费用.而R软件是一款免费的统计软件,它提供了有弹性的,互动的环境来分析,展示数据,_3且提供若干统计程序包以及一些集成统计计算工具和函数,使得用户可以灵活机动地进行数据分析,简化了数据分析过程.它可以完成大多数模型的统计计算,并帮助用户分析计算结果.本文将结合实例,展示如何在R软件中实现对Logistic模型的统计分析.2Logistic模型的拟合回归模型(*)的拟合主要是求模型中的参数估计值,Logistic模型的参数估计一般采用极大似然估计,对于得到的似然方程xT(YP)10:=0,可采用GaussNewton迭代法口进行计算.鉴于Logistic模型是基于二项分布族的广义线性模型,因此在R软件中可通过glm语句建立回归关系,再用summary语句得到其详细结果.3Logistic模型的诊断对得到的模型拟合结果,R软件还可以对其进行诊断,可用influence.measures语句进行影响分析.通过该语句可以直接得到像cook距离,帽子矩阵对角元值hii等基本统计量.再通过函数运算可得到其他统计量,进而进行诊断分析.4数值实例资料为5O位急性白血病人的生存时间Y,入院后取得的淋巴浸润等级x以及出院后有无巩固治疗x2的关系.b该数据曾被多次引用,周雁曾利用Splus和Matlab软件进行拟合,而谢书培J也只用了30例数据,得到的结果不够详细.根据数据,用R软件拟合Logistic模型,具体过程如下:X<一c(0,2,2,2,2,0,0,2,2,0,2,0,0,0,0,0,0,0,0,2,0,0,0,0,0,2,0,0,0,0,0,2,0,2,0,0,2,0,2,0,0,0,2,0,0,2,O,0,0,0)x2<一c(0,0,0,0,0,0,1,0,0,0,1,0,0,0,1,0,1,1,l,1,1,0,1,1,1,0,0,0,0,0,0,0,1,0,1,1,1,0,1,0,0,1,1,1,1,1,1,1,1,1)Y<一C(0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,收稿日期20110908作者简介光琳(1980一),男,河南灵宝人,江苏联合职业技术学院连云港财经分院讲师,扬州大学数学科学学院硕士生.1,1,1,1,1,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,1,1,1,1,1,1)life<一data.frane(X1,x2,Y)log<一glm(Y一?Xl+x2,family=binomial,data=life)summary(1og)根据拟合的结:果,得到l0g()=一1.64l90.7070X1+2.7844X2,即.exp(一1.一64190.7070X1+2.7844X2)一1+exp(一1.64190.7070X1+2.7844x2)R软件可以用来进行预测分析:-3pre<一predict(1og,data.frame(Xl=2,x2=0)P<一exp(pre)/(1+exp(pre);P得结果0.04496518;pre<一predict(1og,data.frame(X1=2,x2=1)P<一exp(pre)l,(1+exp(pre);P得结果0.4325522.比较上述两个结果可见,当入院后取得的淋巴浸润等级X=2时,出院后巩固治疗(x2=1)的存活概率是没有巩固治疗(X2=O)存活概率的9.62倍.接下来,再对模型进行诊断,在程序窗口中输入influence.measures(】0g),回车后得到一些结果,其中cooks,hii这两列分别为cook距离,帽子矩阵对角元值,该两项是统计诊断中的基cook距离),学生化残差ri=J.ti.,/1一hii通过下面的程序来计算各统计量:P<一exp(一1.64190.7070*Xl+2.7844*X2)/(1+exp(一1.64190.7070*xl+2.7844*X2);W<一P(1一P);rPi<一(Yp)/sqrt(w);rPi得到pearson残差rP;hii<一hatvaues(1og);SCi<一rPi/(1一hii);SCi得到Score统计量SC;CD<一cooks.distance(1og);LDi<一3*CDi;LDi得到似然距离LDi;学生化残差<一sqrt(SC;);学生化残差得到学生化残差.为了看得清楚,我们可用如下语句将所得到的这些诊断统计量写成矩阵的形式:matrix(c(1:50,rPi,h.,CDi,LDi,SCi,学生化残差),nrow=5o).根据得到的这些统计诊断量(限于篇幅,得到的结果不具体列出,读者可按上述程序自行得到结果),便可以进行诊断分析.这里之所以用学生化残差,而未用ri,仅仅是为了表明在R软件的语句中,不但可以用英文作变量名,还可以用中文作变量名.R软件还可以用来作图,用plot语句作散点图,tj可再用lines语句进行连线.如想作出cook距离和标准化pearson残差的散点图,用如下程序:序号<一1:50plot(序号,CDi);lines(序号,CDi)pot(序号,学生化残差);lines(序号,学生化残差)可分别得到如下两个诊断统计量的散点图:图2标准化pearson残差l07从图1中可较清楚地看到第l6,22和4l号点为强影响点,第20,43和46号点也有较大影响;从图2中可看到16,22和4l号点为异常点.还可以用text(i,Yi,labels=i,adj=1)将图中的第i个点标记出来,也可以用title(main=“)为图加上标题.其他散点图亦可类似上述方法得到.5小结通过上述实例分析可见,利用R软件可以很方便地计算出各种所需量,从而实现对Logistic模型的统计分析,这里只简单列举参数估计和统计诊断量的计算,其他如参数检验,局部影响分析等,R软件也可方便地实现,从而大大减少了计算量.如果数据较多,而且是Excel表格形式或者其他类型数据库的,R软件也可方便地导人,而不用逐个输入.总之,利用R软件进行Logistic模型的统计分析方便有效,在数理统计学的教学中也可以加入这方面的内容.参考文献1韦博成,鲁国斌.统计诊断引论M.南京:东南大学出版社,1991:442459.2王怀亮.基于R软件的logistic回归模型建模J.黑龙江对外经贸,2010(11):117118.3薛毅,陈立萍.统计建模与R软件M.北京:清华大学出版社,2OOr7:367372.4周雁.广义线性模型的诊断与实例分析J.四)l1大学学报,200-/,44(6):l1631168.5谢书培,韩俊林.Logistic回归模型的统计诊断与实例分析J.聊城大学,2011,24(1):2731.6关彦辉.R软件在概率统计教学中的应用J.现代计算机,2oo9(12):879o.StatisticalAnalysisofLogisticModelsintheRSofewareGUANGLin(LianyungangBranchofFinanceandEconomics,JiangsuUnionTechnicalInstitute,Lianyungang222003,China)Abstract:Withthedevelopmentofregressionanalysistheory,logisticregressionmodelbecomesanimportanttopicinmedicine.engineeringandotherfi