《多元統(tǒng)計(jì)分析-基于R（第3版）》課件 第5、6章-多元線性模型、廣義線性模型

5.1多元正態(tài)分布5.2多元線性模型5.3變量選擇5.4回歸診斷5.5回歸預(yù)測第5章多元線性模型2主編：費(fèi)宇2024/7/315.1.1多元正態(tài)分布的定義5.1多元正態(tài)分布2024/7/313主編：費(fèi)宇5.1.1多元正態(tài)分布的定義5.1多元正態(tài)分布2024/7/314主編：費(fèi)宇5.1.2多元正態(tài)分布的性質(zhì)5.1多元正態(tài)分布2024/7/315主編：費(fèi)宇5.1.2多元正態(tài)分布的性質(zhì)5.1多元正態(tài)分布2024/7/316主編：費(fèi)宇5.2多元線性模型5.2.1模型定義其中x1,…,xp是非隨機(jī)的自變量,y是隨機(jī)的因變量,β0是常數(shù)項(xiàng),β1,…,βp是回歸系數(shù),ε是隨機(jī)誤差項(xiàng).7主編：費(fèi)宇2024/7/315.2.1.模型定義模型(5.1)的樣本形式其中i=1,…,n,表示有n組觀測值.8主編：費(fèi)宇2024/7/315.2.1.模型定義模型(5.2)的矩陣形式其中9主編：費(fèi)宇2024/7/31例5.1(數(shù)據(jù)文件為eg5.1)10主編：費(fèi)宇表5-1抽樣調(diào)查得到的36個學(xué)生的相關(guān)成績2024/7/31yx1x2x3x4x5yx1x2x3x4x5858386909076456065608678909288879280768081758075787076738590888582868580807281829088828081868790………………………………928385908580627865608588788482739083878083857883模型（5.3）的參數(shù)的最小二乘估計(jì)的最小二乘估計(jì)

5.2.2模型的參數(shù)估計(jì)和檢驗(yàn)2024/7/3111主編：費(fèi)宇1.回歸方程的顯著性檢驗(yàn)其中是回歸平方和，而

是殘差平方和，拒絕域?yàn)?/p>

5.2.2模型的參數(shù)估計(jì)和檢驗(yàn)2024/7/3112主編：費(fèi)宇5.2.2模型的參數(shù)估計(jì)和檢驗(yàn)2024/7/3113主編：費(fèi)宇例5.1續(xù)1(數(shù)據(jù)文件為eg5.1)14主編：費(fèi)宇建立y關(guān)于x1、x2、x3、x4和x5的線性回歸方程,并對方程和回歸系數(shù)進(jìn)行顯著性檢驗(yàn).2024/7/31yx1x2x3x4x5yx1x2x3x4x5858386909076456065608678909288879280768081758075787076738590888582868580807281829088828081868790………………………………928385908580627865608588788482739083878083857883#例5.1回歸分析：全變量回歸setwd("C:/data")#設(shè)定工作路徑d5.1<-read.csv("exam5.1.csv",header=T)#將exam5.1.csv數(shù)據(jù)讀入到d5.1中l(wèi)m.exam<-lm(y~x1+x2+x3+x4+x5,data=d5.1)#建立y關(guān)于x1,x2,x3,x4和x5的線性回歸方程,數(shù)據(jù)為d5.1summary(lm.exam)#給出回歸系數(shù)的估計(jì)和顯著性檢驗(yàn)等例5.1續(xù)1(數(shù)據(jù)文件為eg5.1)15主編：費(fèi)宇R程序?yàn)?Call:lm(formula=y~x1+x2+x3+x4+x5,data=d5.1)Residuals:Min1QMedian3QMax-10.0696-1.7983-0.15352.93616.8726

Coefficients:EstimateStd.ErrortvaluePr(>|t|)(Intercept)-32.7353415.35701-2.1320.0413*x10.162710.150311.0820.2877x20.227840.138351.6470.1100x30.881160.111087.9337.46e-09***x4-0.051360.15476-0.3320.7423x50.168870.143761.1750.2494回歸分析結(jié)果為2024/7/3116主編：費(fèi)宇---Signif.codes:0‘***’0.001‘**’0.01‘*’0.05‘.’0.1‘’1

Residualstandarderror:4.021on30degreesoffreedomMultipleR-squared:0.8945,AdjustedR-squared:0.877F-statistic:50.89on5and30DF,p-value:9.359e-14例5.1續(xù)(數(shù)據(jù)文件為eg5.1)2024/7/3117主編：費(fèi)宇回歸方程的F值為50.89,相應(yīng)的p值為9.359

10-14,說明回歸方程是顯著的;但t檢驗(yàn)對應(yīng)的p值則顯示:常數(shù)項(xiàng)和x3是顯著的,而x1、x2、x4和x5不顯著.最優(yōu)模型一般滿足2個條件（1）模型反映了變量間的真實(shí)關(guān)系

（2）模型包含的變量盡量少例5.1（續(xù)2）建立y關(guān)于x1、x2、x3、x4和x5的線性回歸方程,并對方程和回歸系數(shù)進(jìn)行顯著性檢驗(yàn)（逐步回歸建立“最優(yōu)方程”）.

5.3變量選擇2024/7/3118主編：費(fèi)宇>#例5.1回歸分析：逐步回歸>#假設(shè)exam5.1.xls中的數(shù)據(jù)已經(jīng)讀入到d5.1中>lm.exam<-lm(y~x1+x2+x3+x4+x5,data=d5.1)#建立全變量回歸方程>lm.step<-step(lm.exam,direction="both")#進(jìn)行逐步回歸5.3變量選擇2024/7/3119主編：費(fèi)宇例5.1(續(xù)2)逐步回歸程序Start:AIC=105.63y~x1+x2+x3+x4+x5

DfSumofSqRSSAIC-x411.78486.83103.76-x1118.95503.99105.01-x5122.31507.36105.25<none>485.05105.63-x2143.85528.90106.74-x311017.441502.49144.335.3變量選擇2024/7/3120主編：費(fèi)宇例5.1(續(xù))回歸結(jié)果Step:AIC=103.76y~x1+x2+x3+x5

DfSumofSqRSSAIC-x1117.91504.73103.06-x5120.57507.40103.25<none>486.83103.76-x2142.99529.81104.80+x411.78485.05105.63-x311112.961599.79144.595.3變量選擇2024/7/3121主編：費(fèi)宇Step:AIC=103.06y~x2+x3+x5

DfSumofSqRSSAIC-x5117.40522.14102.28<none>504.73103.06+x1117.91486.83103.76+x410.74503.99105.01-x2170.76575.50105.78-x311848.492353.23156.485.3變量選擇2024/7/3122主編：費(fèi)宇Step:AIC=102.28y~x2+x3

DfSumofSqRSSAIC<none>522.14102.28+x5117.40504.73103.06+x1114.74507.40103.25+x410.25521.89104.26-x2166.64588.78104.60-x311953.302475.43156.305.3變量選擇2024/7/3123主編：費(fèi)宇Call:lm(formula=y~x2+x3,data=d5.1)Residuals:Min1QMedian3QMax-10.4395-2.5508-0.44592.73677.2345

Coefficients:EstimateStd.ErrortvaluePr(>|t|)(Intercept)-18.842907.58902-2.4830.0183*x20.249230.121442.0520.0481*x30.968040.0871311.1111.09e-12***---Signif.codes:0‘***’0.001‘**’0.01‘*’0.05‘.’0.1‘’1

Residualstandarderror:3.978on33degreesoffreedomMultipleR-squared:0.8865,AdjustedR-squared:0.8796F-statistic:128.8on2and33DF,p-value:2.566e-16回歸模型匯總信息:summary(lm.step)2024/7/3124主編：費(fèi)宇5.4.1殘差分析和異常點(diǎn)探測

殘差向量e是模型中隨機(jī)誤差項(xiàng)

的估計(jì)，殘差分析可以診斷模型的基本假定是否成立。

5.4

回歸診斷2024/7/3125主編：費(fèi)宇例5.2

計(jì)算例5.1得到的逐步回歸模型lm.step的普通殘差和標(biāo)準(zhǔn)化殘差,判斷可能存在的異常點(diǎn),畫出相應(yīng)的殘差散點(diǎn)圖,并直觀判斷模型的基本假定是否成立.5.4

回歸診斷2024/7/3126主編：費(fèi)宇#例5.2#假設(shè)由例5.1已經(jīng)得到逐步回歸模型lm.stepy.res<-residuals(lm.exam)#計(jì)算模型lm.exam的普通殘差y.rst<-rstandard(lm.step)#計(jì)算回歸模型lm.step的標(biāo)準(zhǔn)化殘差print(y.rst)#輸出回歸模型lm.step的標(biāo)準(zhǔn)化殘差y.rsty.fit<-predict(lm.step)#計(jì)算回歸模型lm.step的預(yù)測值plot(y.res~y.fit)#繪制以普通殘差為縱坐標(biāo)，預(yù)測值為橫坐標(biāo)的散點(diǎn)圖plot(y.rst~y.fit)#繪制以標(biāo)準(zhǔn)化殘差為縱坐標(biāo)，預(yù)測值為橫坐標(biāo)的散點(diǎn)圖5.4

回歸診斷2024/7/3127主編：費(fèi)宇分別采用residuals(),rstandard()和rstudent()來計(jì)算普通殘差,標(biāo)準(zhǔn)化殘差和學(xué)生化殘差.123456-1.226479490.701233481.85465439-0.18487397-0.731575470.14591132789101112-0.651653781.37662024-0.28171298-0.96473838-0.798622470.81284419131415161718-0.48393343-1.176685880.913377160.564389020.658766891.49006874192021222324-2.871217390.527102680.81076269-0.668013511.20184149-1.040201892526272829300.32282704-0.04616114-0.159120010.21602487-0.21306706-0.23026109313233343536-0.24302334-2.03567204-0.33183300-0.073548931.804380090.737029325.4

回歸診斷2024/7/3128主編：費(fèi)宇回歸模型lm.step的標(biāo)準(zhǔn)化殘差y.rst如下第19號點(diǎn)是異常點(diǎn)2024/7/3129主編：費(fèi)宇圖5-1例5.2中的普通殘差圖(左)和標(biāo)準(zhǔn)化殘差圖(右)例5.3

通過方差穩(wěn)定變換來更新例5.1得到的逐步回歸模型lm.step,并計(jì)算更新后模型的標(biāo)準(zhǔn)化殘差,畫出相應(yīng)的殘差散點(diǎn)圖,并直觀判斷模型的基本假定是否成立.5.4

回歸診斷2024/7/3130主編：費(fèi)宇#例5.3#假設(shè)由例5.1已經(jīng)得到逐步回歸模型lm.steplm.step_new<-update(lm.step,log(.)~.)#對模型進(jìn)行對數(shù)變換y.rst<-rstandard(lm.step_new)#計(jì)算lm.step_new的標(biāo)準(zhǔn)化殘差y.fit<-predict(lm.step_new)#計(jì)算lm.step_new的預(yù)測值plot(y.rst~y.fit)#繪制以標(biāo)準(zhǔn)化殘差為縱坐標(biāo),預(yù)測值為橫坐標(biāo)的散點(diǎn)圖（見圖5-2）采用對數(shù)變換來解決方差非齊問題2024/7/3131主編：費(fèi)宇對數(shù)變換后：第19號點(diǎn)是異常點(diǎn)2024/7/3132主編：費(fèi)宇圖5-2例5.3中的標(biāo)準(zhǔn)化殘差圖lm.exam<-lm(log(y)~x1+x2+x3+x4,data=d5.1[-c(19),])#去掉第19號觀測值再建立全變量回歸方程lm.step<-step(lm.exam,direction="both")#用一切子集回歸法來進(jìn)行逐步回歸y.rst<-rstandard(lm.step)#計(jì)算回歸模型lm.step的標(biāo)準(zhǔn)化殘差y.fit<-predict(lm.step)#計(jì)算回歸模型lm.step的預(yù)測值plot(y.rst~y.fit)#繪制以標(biāo)準(zhǔn)化殘差為縱坐標(biāo)，預(yù)測值為橫坐標(biāo)的散點(diǎn)圖去掉19號觀測值再回歸2024/7/3133主編：費(fèi)宇殘差幾乎全部落在[-2，2]區(qū)域內(nèi)2024/7/3134主編：費(fèi)宇圖5-3例5.3中的標(biāo)準(zhǔn)化殘差圖:去掉19號觀測值5.4.2回歸診斷：一般的方法殘差分析無法分析模型的影響點(diǎn),即探測哪些點(diǎn)對模型的推斷有重要影響,本節(jié)給出的回歸診斷方法,可以診斷模型的基本假定是否成立,哪些值是異常點(diǎn),哪些點(diǎn)是強(qiáng)影響點(diǎn).在R中,函數(shù)plot()和influence.measures()可以用來繪制診斷圖和計(jì)算診斷統(tǒng)計(jì)量.5.4

回歸診斷2024/7/3135主編：費(fèi)宇例5.4對例5.3得到的逐步回歸模型lm.step_new進(jìn)行回歸診斷分析.5.4回歸診斷2024/7/3136主編：費(fèi)宇#例5.4#假定由例5.3已經(jīng)獲得模型lm.step_newpar(mfrow=c(2,2))#在一個2×2網(wǎng)格中創(chuàng)建4個繪圖區(qū)plot(lm.step_new)#繪制模型診斷圖influence.measures(lm.step_new)#計(jì)算各個觀測值的診斷統(tǒng)計(jì)量運(yùn)行上述程進(jìn)行序可得回歸診斷圖(圖5-4)和如下36個觀測值對應(yīng)的診斷統(tǒng)計(jì)量的值.第11，19和33號觀測值被診斷為強(qiáng)影響點(diǎn).2024/7/3137主編：費(fèi)宇I(lǐng)nfluencemeasuresoflm(formula=log(y)~x2+x3，data=d5.1)：dfb.1_dfb.x2dfb.x3dffitcov.rcook.dhatinf10.172353-0.052013-1.36e-01-0.291711.0522.82e-020.06622-0.0789410.0622346.53e-030.110451.1604.17e-030.069130.196029-0.049262-1.17e-010.373080.8364.31e-020.0383…110.014319-0.1934602.46e-01-0.270291.3462.49e-020.2065*12-0.0378780.0103993.07e-020.060361.1881.25e-030.0806…19-2.0235661.0016117.40e-01-2.218990.2329.50e-010.1645*20-0.0066640.071366-8.14e-020.147121.0937.33e-030.0419…33-0.0808070.140284-1.13e-01-0.152341.4647.96e-030.2556*第3,19和35號觀測值可能是異常點(diǎn)和強(qiáng)影響點(diǎn).2024/7/3138主編：費(fèi)宇圖5-4例5.4的回歸診斷圖回歸預(yù)測分為點(diǎn)預(yù)測和區(qū)間預(yù)測兩種,可以采用函數(shù)predict()來實(shí)現(xiàn).例5.5

給定解釋變量x2=80，x3=90，利用例5.1得到的回歸模型對y進(jìn)行點(diǎn)預(yù)測和區(qū)間預(yù)測(置信度為95%).5.5回歸預(yù)測2024/7/3139主編：費(fèi)宇#例5.5#假定由例5.1已經(jīng)獲得模型lm.steppreds<-data.frame(x2=80,x3=90)#給定解釋變量x2和x3的值predict(lm.step,newdata=preds,interval="prediction",level=0.95)#進(jìn)行點(diǎn)預(yù)測和區(qū)間預(yù)測5.5回歸預(yù)測2024/7/3140主編：費(fèi)宇點(diǎn)預(yù)測和區(qū)間預(yù)測的程序如下:>predict(lm.step,newdata=preds,interval="prediction",level=0.95)#區(qū)間預(yù)測fitlwrupr188.2190779.7974896.640675.5回歸預(yù)測2024/7/3141主編：費(fèi)宇運(yùn)行上述程序可得y的點(diǎn)預(yù)測和區(qū)間預(yù)測的結(jié)果如下:程序中選項(xiàng)interval="prediction"表示要給出預(yù)測區(qū)間,選項(xiàng)level=0.95表示置信水平是95%.計(jì)算結(jié)果y的點(diǎn)預(yù)測為88.22,預(yù)測區(qū)間為[79.80,96.64].2024/7/31主編：費(fèi)宇2024/7/31主編：費(fèi)宇43

2024/7/31主編：費(fèi)宇44第6章廣義線性模型6.1廣義線性模型的定義6.2Logistic模型6.3Probit模型6.4多項(xiàng)Logit模型6.5泊松對數(shù)線性模型6.6零膨脹計(jì)數(shù)模型6.7多項(xiàng)分布對數(shù)線性模型2024/7/31主編：費(fèi)宇456.1

廣義線性模型的定義第5章我們研究了多元線性模型,該模型的一個重要假定是因變量是連續(xù)型的變量(通常假定服從正態(tài)分布),但在許多情況下,這種假定并不合理,例如下面這兩種情況.(1)結(jié)果變量可能是類型變量.二值分類變量和多分類變量.(比如:是/否,差/一般/良好/優(yōu)秀等)顯然都不是連續(xù)型變量.2024/7/31主編：費(fèi)宇466.1廣義線性模型概述(2)結(jié)果變量可能是計(jì)數(shù)型變量(比如:一周交通事故的數(shù)目)這類變量都是非負(fù)的有限值,而且它們的均值和方差通常是相關(guān)的(一般線性模型假定因變量是正態(tài)變量,而且相互獨(dú)立).普通線性回歸模型(5.3)假定因變量y服從正態(tài)分布,其均值滿足關(guān)系式:μ=Xβ,這表明因變量的條件均值是自變量的線性組合.本章介紹六種常見的廣義線性模型：Logistic模型、Probit模型、多項(xiàng)Logit模型、泊松對數(shù)線性模型、零膨脹計(jì)數(shù)模型和多項(xiàng)分布對數(shù)線性模型.2024/7/31主編：費(fèi)宇476.1廣義線性模型的定義1.廣義線性模型的定義：(1)隨機(jī)成分:設(shè)y1,y2,…,yn是來自于指數(shù)分布族的隨機(jī)樣本,即yi的密度函數(shù)為其中ai(.),b(.),ci(.)是已知函數(shù),參數(shù)αi是典則參數(shù),?是散度參數(shù).2024/7/31主編：費(fèi)宇481.廣義線性模型的定義：(2)連接函數(shù):設(shè)yi的均值為μi而函數(shù)g(.)是單調(diào)可微的連接函數(shù),使得其中是協(xié)變量,是未知參數(shù)向量.指數(shù)分布族正態(tài)分布二項(xiàng)分布泊松分布2024/7/31主編：費(fèi)宇492024/7/31主編：費(fèi)宇502.正態(tài)線性回歸模型正態(tài)分布屬于指數(shù)分布族,其密度函數(shù)為與(6.1)對照可知2024/7/31主編：費(fèi)宇512.正態(tài)線性回歸模型只要取連接函數(shù)為,則正態(tài)線性回歸模型滿足廣義線性模型的定義.類似的,容易驗(yàn)證,二項(xiàng)分布和泊松分布都屬于指數(shù)分布族.下面介紹實(shí)際中應(yīng)用廣泛的Logistic模型、Probit模型、多項(xiàng)Logit模型、泊松對數(shù)線性模型、零膨脹計(jì)數(shù)模型和多項(xiàng)分布對數(shù)線性模型.2024/7/31主編：費(fèi)宇526.2Logistic模型1.模型定義 設(shè)yi服從參數(shù)為pi的二項(xiàng)分布,則μi=E(yi)=pi采用邏輯連接函數(shù),即這個廣義線性模型稱為Logistic模型.2024/7/31主編：費(fèi)宇53例6.1(數(shù)據(jù)文件為eg6.1)表6-1

某城市48個家庭的調(diào)查數(shù)據(jù)2024/7/31主編：費(fèi)宇542.模型的參數(shù)估計(jì)和檢驗(yàn)采用R軟件中的廣義線性模型過程glm()可以完成回歸系數(shù)的估計(jì),以及模型回歸系數(shù)的顯著性檢驗(yàn).程序如下：#例6.1廣義線性模型:Logistic模型setwd("C:/data")#設(shè)定工作路徑d6.1<-read.csv("exam6.1.csv",header=T)#將exam6.1.csv數(shù)據(jù)讀入到d6.1中g(shù)lm.logit<-glm(y~x1+x2,family=binomial(link=logit),data=d6.1)#建立y關(guān)于x1，x2的logistic回歸模型，數(shù)據(jù)為d6.1summary(glm.logit)#模型匯總注意邏輯連接函數(shù)是二項(xiàng)分布的典則連接函數(shù)，是默認(rèn)的連接函數(shù)，因此代碼中的(link=logit)可以省略.2024/7/31主編：費(fèi)宇55運(yùn)行以上程序可得如下結(jié)果：Call:glm(formula=y~x1+x2，family=binomial(link=logit)，data=d6.1)DevianceResiduals:Min1QMedian3QMax-2.30297-0.198320.022830.202511.59258Coefficients:EstimateStd.ErrorzvaluePr(>|z|)(Intercept)-7.531152.56352-2.9380.00331**x10.439560.138643.1700.00152**x2-0.081031.24747-0.0650.94821---Signif.codes:0‘***’0.001‘**’0.01‘*’0.05‘.’0.1‘’1(Dispersionparameterforbinomialfamilytakentobe1)Nulldeviance:61.105on47degreesoffreedomResidualdeviance:17.643on45degreesoffreedomAIC:23.643NumberofFisherScoringiterations:82024/7/31主編：費(fèi)宇56注意到x2對應(yīng)的p值(0.948)比較大，即x2不顯著，所以考慮采用逐步回歸.glm.step<-step(glm.logit)#逐步回歸summary(glm.step)#給出模型回歸系數(shù)的估計(jì)和顯著性檢驗(yàn)等運(yùn)行以上程序可得如下結(jié)果Start:AIC=23.64y~x1+x2

DfDevianceAIC-x2117.64721.647<none>17.64323.643-x1159.00863.008

Step:AIC=21.65y~x1

DfDevianceAIC<none>17.64721.647-x1161.10563.105>summary(glm.step)#給出模型回歸系數(shù)的估計(jì)和顯著性檢驗(yàn)等2024/7/31主編：費(fèi)宇57Call:glm(formula=y~x1,family=binomial(link=logit),data=d6.1)DevianceResiduals:Min1QMedian3QMax-2.28859-0.197030.022760.204001.60887

Coefficients:EstimateStd.ErrorzvaluePr(>|z|)(Intercept)-7.56822.5101-3.0150.00257**x10.43960.13873.1690.00153**---Signif.codes:0‘***’0.001‘**’0.01‘*’0.05‘.’0.1‘’1

(Dispersionparameterforbinomialfamilytakentobe1)

Nulldeviance:61.105on47degreesoffreedomResidualdeviance:17.647on46degreesoffreedomAIC:21.647

NumberofFisherScoringiterations:8模型預(yù)測如果要預(yù)測年收入為20萬元（x1=20）、家里有孩子（x2=1）的家庭有購買住房的可能性，可以采用以下命令：>yp<-predict(glm.step,data.frame(x1=20))>p.fit<-exp(yp)/(1+exp(yp));p.fit#估計(jì)x1=20時y=1的概率10.7728122容易看出，當(dāng)x1=20，x2=1時，估計(jì)y=1的概率約為0.77，即年收入為20萬元、家里有孩子的家庭有購買住房的可能性約為77%.2024/7/31主編：費(fèi)宇586.3Probit模型2024/7/31主編：費(fèi)宇596.3Probit模型:例6.1（續(xù)）2024/7/31主編：費(fèi)宇60運(yùn)行以上程序可得如下結(jié)果：Call:glm(formula=y~x1,family=binomial(link=probit),data=d6.1)

DevianceResiduals:Min1QMedian3QMax-2.2493-0.15220.00180.17681.6024

Coefficients:EstimateStd.ErrorzvaluePr(>|z|)(Intercept)-4.340281.27539-3.4030.000666***x10.249890.069443.5990.000320***---Signif.codes:0‘***’0.001‘**’0.01‘*’0.05‘.’0.1‘’1

(Dispersionparameterforbinomialfamilytakentobe1)

Nulldeviance:61.105on47degreesoffreedomResidualdeviance:17.349on46degreesoffreedomAIC:21.349

NumberofFisherScoringiterations:92024/7/31主編：費(fèi)宇612024/7/31主編：費(fèi)宇626.4多項(xiàng)Logit模型2024/7/31主編：費(fèi)宇63例6.2某城市48個家庭的調(diào)查數(shù)據(jù)2024/7/31主編：費(fèi)宇64問題：根據(jù)這個數(shù)據(jù)建立多項(xiàng)分布回歸模型并估計(jì)年收入為20萬元、家里有孩子的家庭有購買住房但還在還貸款的可能性.解:采用nnet程序包中的multinom()可以完成多項(xiàng)logit模型的擬合.#例6.2廣義線性模型:多項(xiàng)分布回歸模型library(nnet)setwd("C:/data")d6.2<-read.csv("exam6.2.csv",header=T)d6.2$x2<-as.factor(d6.2$x2)#將x2這一列因子化mlog<-multinom(y~x1+x2,data=d6.2)#建立模型summary(mlog)#查看所擬合的模型2024/7/31主編：費(fèi)宇65運(yùn)行以上程序可得如下結(jié)果：Call:multinom(formula=y~x1+x2，data=d6.2)Coefficients:(Intercept)x1x22-7.4438920.4329375-0.067896533-17.3785220.7438569-0.57429520Std.Errors:(Intercept)x1x222.5703380.13962821.24601334.4477300.18612381.704516ResidualDeviance:37.79579AIC:49.79579注意到x2對應(yīng)標(biāo)準(zhǔn)誤相對于x2的系數(shù)比較大，所以估計(jì)x2可能不顯著，采用step()函數(shù)對模型進(jìn)行逐步回歸.2024/7/31主編：費(fèi)宇66mlog.s<-step(mlog)#對mlog進(jìn)行逐步回歸summary(mlog.s)#查看所擬合的模型運(yùn)行以上程序可得如下結(jié)果：Call:multinom(formula=y~x1,data=d6.2)

Coefficients:(Intercept)x12-7.4794080.43324433-17.2933710.7313709

Std.Errors:(Intercept)x122.5180900.139753034.4241140.1834096

ResidualDeviance:37.98674AIC:45.986742024/7/31主編：費(fèi)宇672024/7/31主編：費(fèi)宇682024/7/31主編：費(fèi)宇692024/7/31主編：費(fèi)宇70估計(jì)48個家庭最可能屬于3類家庭中的哪一類？2024/7/31主編：費(fèi)宇716.5泊松對數(shù)線性模型1.模型的定義設(shè)y服從參數(shù)λ為的泊松分布,則μ=E(y)=λ,采用對數(shù)連接函數(shù),即這個廣義線性模型稱為泊松對數(shù)線性模型.2024/7/31主編：費(fèi)宇72例6.3(數(shù)據(jù)文件為eg6.3)表6-3

Breslow癲癇數(shù)據(jù)2024/7/31主編：費(fèi)宇73例6.3(數(shù)據(jù)文件為eg6.3)這個數(shù)據(jù)是robust包中的Breslow癲癇數(shù)據(jù)(Breslow,1993).我們討論在治療初期的八周內(nèi),癲癇藥物對癲癇發(fā)病數(shù)的影響,響應(yīng)變量為八周內(nèi)癲癇發(fā)病數(shù)(y),預(yù)測變量為前八周內(nèi)的基礎(chǔ)發(fā)病次數(shù)(x1),年齡(x2)和治療條件(x3),其中治療條件是二值變量,x3=0表示服用安慰劑,x3=1表示服用藥物.根據(jù)這個數(shù)據(jù)建立泊松對數(shù)線性模型并對模型的系數(shù)進(jìn)行顯著性檢驗(yàn).2024/7/31主編：費(fèi)宇742.模型的參數(shù)估計(jì)和檢驗(yàn)采用R軟件中的廣義線性模型過程glm()來建立泊松對數(shù)線性模型并對模型的系數(shù)進(jìn)行顯著性檢驗(yàn).程序如下：#例6.3廣義線性模型:泊松對數(shù)線性模型setwd("C:/data")d6.3<-read.csv("exam6.3.csv",header=T)#將exam6.3.scv數(shù)據(jù)讀入到d6.3中g(shù)lm.ln<-glm(y~x1+x2+x3,family=poisson(link=log),data=d6.3)#建立y關(guān)于x1,x2,x3的泊松對數(shù)線性模型summary(glm.ln)#模型匯總,給出模型回歸系數(shù)的估計(jì)和顯著性檢驗(yàn)等泊松分布的默認(rèn)連接函數(shù)是對數(shù)連接函數(shù)，因此代碼中的(link=log)可以省略.2024/7/31主編：費(fèi)宇75運(yùn)行以上程序可得如下結(jié)果：Call:glm(formula=y~x1+x2+x3,family=poisson(link=log),data=data6.3)DevianceResiduals:Min1QMedian3QMax-6.0569-2.0433-0.93970.792911.0061Coefficients:EstimateStd.ErrorzvaluePr(>|z|)(Intercept)1.94882590.135619114.370<2e-16***x10.02265170.000509344.476<2e-16***x20.02274010.00402405.6511.59e-08***x3-0.15270090.0478051-3.1940.0014**Signif.codes:0‘***’0.001‘**’0.01‘*’0.05‘.’0.1‘’1(Dispersionparameterforpoissonfamilytakentobe1)Nulldeviance:2122.73on58degreesoffreedomResidualdeviance:559.44on55degreesoffreedomAIC:850.71NumberofFisherScoringiterations:52024/7/31主編：費(fèi)宇762.模型的參數(shù)估計(jì)和檢驗(yàn)于是得回歸模型:從檢驗(yàn)結(jié)果可以看出:x1和x2的系數(shù)都顯著,說明基礎(chǔ)發(fā)病次數(shù)(x1),年齡(x2)和治療條件(x3)對八周內(nèi)癲癇發(fā)病數(shù)(y)重要影響.年齡(x2)的回歸系數(shù)為0.0227,表明保持其他預(yù)測變量不變,年齡增加1歲,癲癇發(fā)病數(shù)的對數(shù)均值將相應(yīng)的增加0.0227.2024/7/31主編：費(fèi)宇772.模型的參數(shù)估計(jì)和檢驗(yàn)在因變量的初始尺度(癲癇發(fā)病數(shù),而不是癲癇發(fā)病數(shù)的對數(shù))上解釋回歸系數(shù)比較容易,因此,指數(shù)化系數(shù):可以看出:保持其他預(yù)測變量不變,年齡增加1歲,癲癇發(fā)病數(shù)將乘以1.023;治療條件變化一個單位(即從安慰劑到藥物),癲癇發(fā)病數(shù)將乘以0.86,換言之,保持基礎(chǔ)癲癇發(fā)病數(shù)和年齡不變,服藥相