Logistic回歸模型

1、logistic回歸模型1 logistic回歸模型的基本知識1.1 logistic模型簡介主要應(yīng)用在研究某些現(xiàn)象發(fā)生的概率，比如股票漲還是跌，公司成功或失敗的概率，以及討論概率與那些因素有關(guān)。顯然作為概率值，一定有，因此很難用線性模型描述概率與自變量的關(guān)系，另外如果接近兩個極端值，此時(shí)一般方法難以較好地反映p的微小變化。為此在構(gòu)建與自變量關(guān)系的模型時(shí)，變換一下思路，不直接研究，而是研究的一個嚴(yán)格單調(diào)函數(shù)，并要求在接近兩端值時(shí)對其微小變化很敏感。于是logit變換被提出來： （1）其中當(dāng)從時(shí)，從，這個變化范圍在模型數(shù)據(jù)處理上帶來很大的方便，解決了上述面臨的難題。另外從函數(shù)的變形可得如下等價(jià)的

2、公式： （2）模型(2)的基本要求是，因變量是個二元變量，僅取0或1兩個值，而因變量取1的概率就是模型要研究的對象。而，其中表示影響的第個因素，它可以是定性變量也可以是定量變量，。為此模型(2)可以表述成： （3）顯然，故上述模型表明是的線性函數(shù)。此時(shí)我們稱滿足上面條件的回歸方程為logistic線性回歸。logistic線性回歸的主要問題是不能用普通的回歸方式來分析模型，一方面離散變量的誤差形式服從伯努利分布而非正態(tài)分布，即沒有正態(tài)性假設(shè)前提；二是二值變量方差不是常數(shù)，有異方差性。不同于多元線性回歸的最小二乘估計(jì)法則(殘差平方和最小)，logistic變換的非線性特征采用極大似然估計(jì)的迭代方

3、法尋求最佳的回歸系數(shù)。因此評價(jià)模型的擬合度的標(biāo)準(zhǔn)變?yōu)樗迫恢刀请x差平方和。定義1 稱事件發(fā)生與不發(fā)生的概率比為 優(yōu)勢比(比數(shù)比 odds ratio 簡稱or)，形式上表示為or= （4）定義2 logistic回歸模型是通過極大似然估計(jì)法得到的，故模型好壞的評價(jià)準(zhǔn)則有似然值來表征，稱-2為估計(jì)值的擬合似然度，該值越小越好，如果模型完全擬合，則似然值為1，而擬合似然度達(dá)到最小，值為0。其中表示的對數(shù)似然函數(shù)值。定義3 記為估計(jì)值的方差-協(xié)方差矩陣，為的標(biāo)準(zhǔn)差矩陣，則稱 （5）為的wald統(tǒng)計(jì)量，在大樣本時(shí)，近似服從分布，通過它實(shí)現(xiàn)對系數(shù)的顯著性檢驗(yàn)。定義4 假定方程中只有常數(shù)項(xiàng)，即各變量的系數(shù)

4、均為0，此時(shí)稱 （6）為方程的顯著性似然統(tǒng)計(jì)量，在大樣本時(shí)，近似服從分布。1.2 logistic模型的分類及主要問題根據(jù)研究設(shè)計(jì)的不同，logistic回歸通常分為成組資料的非條件logistic回歸和配對資料的條件logistic回歸兩種大類。還兼具兩分類和多分類之分，分組與未分組之分，有序與無序變量之分。具體如下：兩分類非條件logistic回歸：分組數(shù)據(jù)的logistic回歸，未分組數(shù)據(jù)的logistic回歸；多分類非條件logistic回歸：無序變量logistic回歸，無序變量logistic回歸；條件logistic回歸：1:1型、1:m型和m:n型logistic回歸。關(guān)于lo

5、gistic回歸，主要研究的內(nèi)容包括：1 模型參數(shù)的估計(jì)及檢驗(yàn)2 變量模型化及自變量的選擇3 模型評價(jià)和預(yù)測問題4 模型應(yīng)用2 logistic模型的參數(shù)估計(jì)及算法實(shí)現(xiàn)2.1 兩分類分組數(shù)據(jù)非條件logistic回歸應(yīng)變量(反應(yīng)變量)分為兩類，取值有兩種，設(shè)事件發(fā)生記為y=1，不發(fā)生記為 y=0，設(shè)自變量是分組數(shù)據(jù)，取有限的幾個值；研究事件發(fā)生的概率與自變量的關(guān)系，其logistic回歸方程為： 或 例2.1.1 分組數(shù)據(jù)1 在一次住房展銷會上，與房地產(chǎn)商簽訂初步購房意向書的有n=325人，在隨后的3個月時(shí)間內(nèi)，只有一部分顧客購買了房屋。購買房屋的顧客記為1，否則記為0。以顧客的年家庭收入(萬

6、元)作為自變量,對數(shù)據(jù)統(tǒng)計(jì)后如表2.1.1所示，建立logistic回歸模型。表2.1.1 購房分組數(shù)據(jù)序號年家庭收入x(萬元)簽訂意向人數(shù)實(shí)際購買人數(shù)11.525822.5321333.5582644.5522255.5432066.5392277.5281688.5211299.51510例2.1.2 藥物療效數(shù)據(jù)2 為考察某藥物療效，隨機(jī)抽取220例病人并分配到治療組合對照組，治療組采用治療藥物，對照組采用安慰劑。治療一段時(shí)間后觀察病人的療效，得到表2.1.2數(shù)據(jù)。設(shè)y為療效指標(biāo)(y=1 有效，y=0無效)，為治療組指標(biāo)(1為治療組，0為對照組)，為年齡組指標(biāo)(1為>45歲，0為其

7、他)。表2.1.2 藥物療效數(shù)據(jù)序號治療分組年齡分組有療效無效合計(jì)111321850210402060301213152400184058上述兩個例子數(shù)據(jù)都是經(jīng)過統(tǒng)計(jì)加工后的分組數(shù)據(jù)，對此類數(shù)據(jù)進(jìn)行l(wèi)ogistic回歸，首先要明確應(yīng)變量對應(yīng)事件的發(fā)生概率如何確定和進(jìn)行l(wèi)ogit變換，其次才能建立logistic回歸。為便于數(shù)據(jù)處理，我們將此類數(shù)據(jù)的格式作個約定，排列格式為（組序號，自變量，該組事件發(fā)生數(shù)，該組總例數(shù)）。表2.1.3 分組數(shù)據(jù)的標(biāo)準(zhǔn)格式序號年家庭收入x(萬元)實(shí)際購買人數(shù)簽訂意向總?cè)藬?shù)11.582522.5133233.5265844.5225255.5204366.522397

8、7.5162888.5122199.51015表2.1.1 改造表序號治療分組年齡分組有效例數(shù)觀察例數(shù)1113250210406030121524001858表2.1.2 改造表經(jīng)過改造后，可得我們關(guān)心的事件的發(fā)生的頻率為 。其中為分組數(shù)，然后作logit變換，即。變換后的數(shù)據(jù)，形式上已經(jīng)可以采用一般的線性回歸的處理方式來估計(jì)回歸參數(shù)了。此時(shí)方程變?yōu)椋寒?dāng)然這樣處理并沒有解決異方差性，當(dāng)較大時(shí)，的近似方差為： （7）所以選擇權(quán)重 ，最后采用加權(quán)最小二乘法估計(jì)參數(shù)。注意，分組數(shù)據(jù)的logistic回歸只適用于大樣本分組數(shù)據(jù)，對小樣本的為分組數(shù)據(jù)不適用，并且以組數(shù)為回歸擬合的樣本量，明顯降低了擬合精

9、度，在實(shí)際應(yīng)用中必須謹(jǐn)慎。求解算法及步驟：1依據(jù)分組數(shù)據(jù)的標(biāo)準(zhǔn)格式，計(jì)算頻率、logit變換和權(quán)重2構(gòu)建加權(quán)最小二乘估計(jì)： （8）令 ，則方程又變成一般的線性回歸模型： （9）3構(gòu)造增廣矩陣?yán)孟シǖ镁仃?，得到估?jì)其中為殘差平方和, 回歸方差各系數(shù)檢驗(yàn)采用 總平方和，回歸平方和總平方和求解相當(dāng)于擬合方程的殘差平方和，故得上式st所以方程的檢驗(yàn)為例2.1.1的求解過程如下（由lllstat統(tǒng)計(jì)軟件計(jì)算）：表2.1.4 數(shù)據(jù)logit變換及權(quán)重 家庭年收入x實(shí)際購買mi簽訂意向ni比例pi邏輯變換logit權(quán)重ni*pi(1-pi)1.500000 8250.320000-0.753772 5.

10、440000 2.500000 13320.406250-0.379490 7.718750 3.500000 26580.448276 -0.207639 14.344828 4.500000 22520.423077 -0.310155 12.692308 5.500000 20430.465116 -0.139762 10.697674 6.500000 22390.564103 0.257829 9.589744 7.500000 16280.571429 0.287682 6.857143 8.500000 12210.571429 0.287682 5.142857 9.50000

11、0 10150.666667 0.693147 3.333333 表2.1.5 回歸模型基本信息 總樣本 9 求解方法 加權(quán)最小二乘僅常數(shù)項(xiàng)beta0 -0.095029 方程f統(tǒng)計(jì)量 51.982160 f分布自由度1,7 方程檢驗(yàn)p值0.000176 總平方和 8.798294 回歸平方和 7.754112 殘差平方和 1.044181 表2.1.6 分組logistic回歸系數(shù)檢驗(yàn) 序號均值回歸系數(shù)系數(shù)標(biāo)準(zhǔn)誤t統(tǒng)計(jì)量自由度df檢驗(yàn)p值常數(shù)項(xiàng) 2.837815-0.8488820.113578-7.47399470.000056家庭年收入x14.9011400.1493230.020711

12、7.20986570.000056表2.1.7 0.086479 -0.014517 -0.014517 0.002876 本例logistic模型的回歸方程：對于多分類無序自變量的logistic回歸，即某個自變量為m個水平的名義變量（如治療方法a,b,c），只需要引入m-1（2個）個啞變量，然后采用上述方法進(jìn)行分析。例2.1.3 研究三種治療方法對不同性別病人的治療效果2，數(shù)據(jù)如表2.1.4表2.1.4 性別和治療法對某病治愈情況的影響性別治療方法有效無效總例數(shù)男a7828106b10111112c6846114女a(chǎn)40545b54559c34640由于治療方法有三種，沒有等級關(guān)系，所以屬

13、于無序的名義變量，故引入兩個啞變量分別代表a和b療法，其中表示方法a, 表示方法b, 表示方法c，將上述數(shù)據(jù)轉(zhuǎn)化成標(biāo)準(zhǔn)格式，得表2.1.5。表2.1.5 性別和治療法對某病治愈情況的影響性別有效總例數(shù)1107810610110111210068114010404500154590003440對于分類數(shù)據(jù)，也可以采用極大似然法進(jìn)行參數(shù)估計(jì)，具體見2.2節(jié)最后部分內(nèi)容。2.2 兩分類未分組(連續(xù))非條件logistic回歸應(yīng)變量取值為0和1，設(shè)事件發(fā)生記為y=1，否則為0，設(shè)自變量，n組觀測數(shù)據(jù)記為,。記，則與的logistic回歸模型是： （10）易知，是均值為的0-1型分布，其分布律為 ，則的

14、似然函數(shù)和對數(shù)似然函數(shù)分別為： 代入,得 （11）記，選取的估計(jì)使得達(dá)到極大，這就是logistic回歸模型的極大似然估計(jì)，該過程的求解需要采用牛頓迭代法。構(gòu)造得分函數(shù)，共k+1個非線性方程組，令其=0求解，其中 (12 )構(gòu)造信息矩陣，即二階導(dǎo)矩陣的負(fù)矩陣，其中 (13 )很明顯，故是一個對稱矩陣。求解算法及步驟：1 根據(jù)公式(12 ) 計(jì)算得分函數(shù)，公式（13）計(jì)算信息矩陣給定初值, k =1 和精度，可取0.0000012 采用牛頓迭代式 , ，通過以下方式求解。構(gòu)造增廣矩陣=，通過對if矩陣作k+1次ij消去變換求解若 或者 或者 ，則轉(zhuǎn)3否則k = k +1，繼續(xù)執(zhí)行第2步3 此時(shí)就

15、是回歸系數(shù)的數(shù)值估計(jì)，k就是迭代次數(shù)，消去變換后的矩陣的前子陣就是方差-協(xié)方差矩陣的估計(jì)陣=v ，下面給出檢驗(yàn)有關(guān)計(jì)算：計(jì)算wald統(tǒng)計(jì)量 ,近似服從分布，檢驗(yàn)p值 標(biāo)準(zhǔn)誤, , 例2.2.1 公共交通調(diào)查數(shù)據(jù)1 在一次關(guān)于公共交通的社會調(diào)查中，調(diào)查項(xiàng)目為“是乘坐公共汽車上下班，還是騎自行車上下班”。因變量y=1表示乘坐公共汽車，y=0表示騎自行車。自變量是年齡，作為連續(xù)變量；是月收入(元)；是性別，=1表示男性，=0表示女性。調(diào)查對象為工薪族群體，數(shù)據(jù)如表2.2.1所示。表2.2.1 公共交通社會調(diào)查 序號年齡月收入性別交通 y11885000221120000323850014239500

16、1528120001631850007361500018421000019469500110481200001155180001125621000113581800011418850101520100010162512001017271300101828150010193095011203210001021331800102233100010233812001024411500102545180011264810001027521500112856180011以下計(jì)算結(jié)果采用lllstat 1.0 軟件得到：表2.2.2 主要計(jì)算結(jié)果序號均值回歸系數(shù)系數(shù)標(biāo)準(zhǔn)誤wald統(tǒng)計(jì)量自由度df檢驗(yàn)p值o

17、r=exp(b)常數(shù)項(xiàng)0.535714-3.6550162.0912233.05476610.0805010.025861年齡1273.2142860.0821680.0521192.48551610.1148991.085639月收入0.4642860.0015170.0018650.66146610.4160431.001518性別36.107143-2.5018441.1578184.66917510.0307090.081934表2.2.3 logistic模型基本信息總樣本28求解方法 極大似然法 & newton迭代迭代次數(shù)(僅beta0) 7(4) -2loglikeli

18、hood(beta) 25.970652 僅常數(shù)項(xiàng)beta0 -0.143101 -2loglikelihood(beta0) 38.673263 方程wald值(相減) 12.702611 方程自由度 4 方程檢驗(yàn)p值 0.012824 對于例2.1.3分組數(shù)據(jù)的極大似然估計(jì)法，主要過程如下： 代入,得 則有 ；其中,分別表示分組i中事件發(fā)生次數(shù)和總觀察數(shù)，如表2.1.4和2.1.5所示。然后可采用newton-raphson迭代法進(jìn)行求解。由lllstat計(jì)算得到如下結(jié)果。表2.2.4 性別和療法對某病治愈的影響(未分組logistic似然估計(jì)法)序號均值回歸系數(shù)系數(shù)標(biāo)準(zhǔn)誤wald統(tǒng)計(jì)量自

19、由度df檢驗(yàn)p值常數(shù)項(xiàng)1.0000001.4183990.29869022.55051310.000002性別0.500000-0.9616180.29979710.28847210.001339治療a0.3333330.5847450.2641084.90196610.026826治療b0.3333331.5607630.31596124.40099310.000001表2.2.5回歸系數(shù)方差矩陣v(beta)(信息矩陣i(beta)的逆矩陣)0.089215-0.072957-0.029931-0.030097-0.0729570.089878-0.0000780.000128-0.029

20、931-0.0000780.0697530.029993-0.0300970.0001280.0299930.0998312.3 條件logistic回歸2,3條件logistic回歸是配對設(shè)計(jì)(病例-對照)中常用的一種統(tǒng)計(jì)分析方法，通過配對方法收集資料：每一配對組可包括一個病例和一個或多個對照，有1：1型、1：m型配對。假設(shè)收集了如下數(shù)據(jù)：表2.3.1 n個1:m配對組,k個協(xié)變量的比例資料 配對組號病例組第1對照組第m個對照組12配對資料用配對的方法來控制影響因素的干擾，并且每個配對組都可以建立一個logistic回歸方程：為此需要估計(jì)的參數(shù)有n個常數(shù)項(xiàng)和k個回歸系數(shù)，配對數(shù)越多估計(jì)的參數(shù)

21、就越多，但是一般的數(shù)據(jù)量難以支撐這樣的估計(jì)，故一般的logistic回歸不適合配對資料。不過在參數(shù)估計(jì)時(shí)，常數(shù)項(xiàng)會被消去，所以方程組減少了n個常數(shù)項(xiàng)的估計(jì)，復(fù)雜度大大降低。對于回歸參數(shù)的估計(jì)采用條件似然函數(shù)替代一般的似然函數(shù)進(jìn)行。對于第i個配對組而言，共有m+1個觀察對象，記為，其中僅有一例發(fā)病，且正好是病例組a發(fā)病，而對照組均沒有發(fā)病的條件概率（類似bayes概率）可以表示成： （14）其中= ，而， (15 )故n個配對組的條件似然函數(shù)表示為： (16 )則對數(shù)似然函數(shù)為： （17）令 ，它是一個與第i個樣本點(diǎn)有關(guān)的k維向量，表示向量中的第g個元素，則有如下得分函數(shù)和信息矩陣：=注意此時(shí)的

22、，沒有常數(shù)項(xiàng)。至此(17)式中的參數(shù)可采用newton-raphson迭代法求解了，初值依然取為0向量。不過該方程的求解已經(jīng)相對復(fù)雜多了。方程似然度檢驗(yàn)和回歸系數(shù)的wald檢驗(yàn)同非條件logistic回歸。例2.3.1 研究肥胖、口服避孕藥雌激素與子宮內(nèi)膜癌的關(guān)系，隨機(jī)抽取20名患者，對于每名患者，在隨機(jī)抽取年齡相近的正常人作為對照。檢測患者與正常人的肥胖程度和雌激素服用情況3。表2.3.1 肥胖和雌激素與子宮內(nèi)膜癌關(guān)系病例-對照研究數(shù)據(jù)配對組病例肥胖病例雌激素對照1肥胖對照1雌激素對照2肥胖對照2雌激素111000021111013110111401000150010016110010711

23、0111811101191011111001010011010110120101011311001114110010151100011601010117010010181110011910010120110100例2.3.1 求解的主要結(jié)果，由lllstat軟件計(jì)算得到：表2.3.2 條件logistic回歸系數(shù)檢驗(yàn) 序號均值(病例)回歸系數(shù)系數(shù)標(biāo)準(zhǔn)誤wald統(tǒng)計(jì)量自由度df檢驗(yàn)p值肥胖0.6500001.8239140.54719211.11039010.000859雌激素0.8500001.5896210.45054412.44836710.000419表2.3.3條件logistic回歸

24、模型基本信息樣本量20求解方法極大似然+牛頓迭代迭代次數(shù)45-2loglikelihood(beta)33.306763-2loglikelihood(0)43.944492方程wald值(相減)10.637728方程自由度2方程檢驗(yàn)p值0.0048982.4 多分類有序反應(yīng)變量logistic回歸在實(shí)際應(yīng)用中，經(jīng)常遇到反應(yīng)變量為多分類有序變量的情況，例如評價(jià)指標(biāo)分為差、中、良、優(yōu)等，各等級之間是有序的。這種資料的logistic回歸分析通常稱為比例比數(shù)模型(累積概率模型) 4，它需要擬合m-1 (m為水平或等級個數(shù))個logistic回歸模型。有序累積概率logistic模型： 或 （18）

25、 （19）有序累積概率模型參數(shù)的極大似然估計(jì)就是尋找參數(shù)使得聯(lián)合概率實(shí)現(xiàn)最大化，由于觀測之間相互獨(dú)立，聯(lián)合概率被分解成邊緣概率之積。而觀測到的概率就是累積概率之差：第i個觀測值對應(yīng)似然值的貢獻(xiàn)取決于觀測到哪一個j值，因此對于次序響應(yīng)的每個j值，取所有的觀測之的乘積，有似然函數(shù)：，其中若，則，否則并且對于任一個觀測而言，只有一個等級事件發(fā)生，即，故有(19)式。其對數(shù)似然函數(shù)如下（對于分組數(shù)據(jù)，似然函數(shù)變?yōu)椋?，分組中各分類例數(shù)）。 （20）其中： （21）然后就可以通過極大似然法，就上newton-raphson方法加以求解參數(shù)了，注意的是。下面給出具體推導(dǎo)求解的詳細(xì)過程。對（20）式進(jìn)行化簡，

26、可得 （22） （23） （24） （25） （26） （27） （28） （29） （30） （31） （32）由此構(gòu)建信息矩陣和，并可迭代求解了。注：若為分組數(shù)據(jù)，上述每項(xiàng)乘以。例2.4.1 研究性別和兩種治療方法對某種疾病療效的影響3,將療效分成效果顯、有效和無效三個等級，根據(jù)試驗(yàn)調(diào)查，得到如下資料。表2.4.1 性別和兩種治療方法對某種疾病療效的影響性別治療方法顯著有效無效合計(jì)女新藥165627傳統(tǒng)671932男新藥52714傳統(tǒng)101011表2.4.2 多分類有序反應(yīng)變量數(shù)據(jù)格式行號性別治療方法頻數(shù)療效等級111161211523116341061510726101937015180

27、122901731000111100021200103計(jì)算結(jié)果，由lllstat統(tǒng)計(jì)軟件給出：表2.4.3 回歸系數(shù)方差矩陣v(beta)(信息矩陣i(beta)的逆矩陣) 0.3747330.324880-0.257757-0.1928230.3248800.323782-0.244457-0.169612-0.257757-0.2444570.2894880.069404-0.192823-0.1696120.0694040.236257表2.4.4 有序分類因變量logistic回歸系數(shù)檢驗(yàn) 序號回歸系數(shù)系數(shù)標(biāo)準(zhǔn)誤wald統(tǒng)計(jì)量自由度df檢驗(yàn)p值常數(shù)項(xiàng)a1-2.6935760.61215

28、519.36137710.000011常數(shù)項(xiàng)a2-1.8120400.56901810.14105910.001450性別1.0523520.5380413.82552810.050477治療方法2.1872720.48606320.24980010.000007表2.4.5 有序分類因變量logistic回歸模型基本信息 樣本分組數(shù)12求解方法極大似然+牛頓迭代迭代次數(shù)17注意：該結(jié)果與 sas, dps不一致。poisson回歸模型1 簡介一般情況下，單位容積水中的細(xì)菌數(shù)，單位時(shí)間內(nèi)某些事件發(fā)生的次數(shù)，單位面積上降落的灰塵的顆粒數(shù)等，都可以用poisson分布來描述。一般poisson分布

29、描述成隨機(jī)變量，概率分布律為：易知，通?？赡苁艿奖姸嘁蛩氐挠绊?，不妨假設(shè)這些因素為(自變量，協(xié)變量），令，對于分組數(shù)據(jù)，poisson分布的期望發(fā)生數(shù)假設(shè)為7： （1）其中為回歸參數(shù)，為第i組的總觀測數(shù)?；貧w模型的似然函數(shù)為poisson分布條件下各個格子概率函數(shù)的乘積，因此poisson分布的極大似然函數(shù)和對數(shù)似然函數(shù)具體形式分別為：代入，得 （2）令 （3）, （4）則可采用newton-raphson迭代法求解參數(shù)的極大似然估計(jì)了。對于僅有常數(shù)項(xiàng)的poisson模型，其估計(jì)值為，用于計(jì)算對數(shù)似然比。2案例分析例1 3 doll和hill(1966)研究英國男性醫(yī)生患冠心病與抽煙、年齡關(guān)系

30、。由于死亡與追蹤人數(shù)和追蹤時(shí)間有關(guān)，故用追蹤人數(shù)和追蹤時(shí)間的乘積(人年)作為觀察單位數(shù)。假定其目標(biāo)變量(死亡人數(shù))近似服從poisson分布，其調(diào)查取樣共74588調(diào)查單位，死亡598例。主要研究因素有抽煙（1為抽煙,0為不抽煙）；調(diào)查對象年齡分成4組(35-44歲,45-54歲,55-64歲,65-74歲)，此為多分類變量，需要設(shè)置三個變量加以區(qū)分，可將其中一個年齡組作為參照組，不妨取35-44歲，計(jì)算時(shí)不考慮對照組信息。表1 英國男性醫(yī)生患冠心病與抽煙、年齡關(guān)系分組抽煙34-44歲45-54歲55-64歲65-74歲死亡數(shù)總例數(shù)11100032523072101001044324831001020628612410001186126635010002187906001001210673700010285710800001282585由lllst