徐州醫(yī)學(xué)院流行病與衛(wèi)生統(tǒng)計(jì)教研室(221002) 王 婷 曾 平 何 鵬
Probit回歸最先由Bliss在1934年提出〔1〕,在進(jìn)行殺蟲(chóng)劑的毒理實(shí)驗(yàn)時(shí),Bliss發(fā)現(xiàn)了一個(gè)有趣的現(xiàn)象:無(wú)論他配制的殺蟲(chóng)劑濃度有多高,在用藥之后總會(huì)有一兩只昆蟲(chóng)還活著;無(wú)論他怎么稀釋殺蟲(chóng)劑,即便只是用裝過(guò)殺蟲(chóng)劑的容器,試驗(yàn)結(jié)果也總會(huì)有幾只昆蟲(chóng)死掉。Bliss原創(chuàng)性地采用概率分布這種新的數(shù)學(xué)思想來(lái)解決殺蟲(chóng)劑實(shí)驗(yàn)時(shí)所遭遇的困境。Probit回歸建立了“劑量”與“使用該劑量時(shí)一只動(dòng)物會(huì)死掉的概率”這兩者間的關(guān)系,因此稱(chēng)為Probit(單位概率)模型,現(xiàn)在已經(jīng)越來(lái)越多地應(yīng)用到二分類(lèi)數(shù)據(jù)的回歸分析中。本文將主要在貝葉斯統(tǒng)計(jì)框架內(nèi)討論P(yáng)robit回歸和參數(shù)后驗(yàn)分布的潛變量Gibbs抽樣。
1.Probit回歸和后驗(yàn)分布
設(shè)解釋變量為X,回歸系數(shù)向量為β,根據(jù)廣義線(xiàn)性模型原理〔2〕,可建立 Probit回歸:
這里假定事件發(fā)生服從參數(shù)為p的Bernoulli分布,n表示樣本量,Φ表示標(biāo)準(zhǔn)正態(tài)的累積概率分布函數(shù),如Φ(1.96)=0.975,Φ-1為累積分布函數(shù)的逆函數(shù),如Φ-1(0.975)=1.96,這樣通過(guò)Probit連接函數(shù)Φ-1將取值為0~1之間的p映射到了整個(gè)實(shí)數(shù)空間。似然。設(shè) g(β)為回歸參數(shù)β的先驗(yàn)分布(prior distribution),貝葉斯統(tǒng)計(jì)和頻率統(tǒng)計(jì)最大的區(qū)別之一就在于假定參數(shù)為隨機(jī)變量,當(dāng)有關(guān)于未知參數(shù)的歷史知識(shí)、主觀認(rèn)識(shí)或者專(zhuān)家意見(jiàn)時(shí),可以選擇有信息先驗(yàn)。當(dāng)對(duì)未知參數(shù)的信息一無(wú)所知、或先驗(yàn)分布有太多參數(shù)需要指定時(shí),認(rèn)為參數(shù)在其空間內(nèi)有等可能的取值概率而不特別偏向某些取值,則選取均勻分布作為先驗(yàn)分布,又稱(chēng)貝葉斯假定。根 據(jù) 貝 葉 斯 原 理, 后 驗(yàn) 分 布〔3〕(posterior distribution)p(β|Y,X)為:
3.潛變量Probit回歸的Gibbs抽樣
設(shè)存在一個(gè)潛在的連續(xù)變量 y*,y*稱(chēng)為特征(trait)或傾向得分(propensity score)〔5〕??紤]以下的模型:均數(shù)(-1.5897)和中位數(shù)(-1.5018)來(lái)看,男孩比女 孩更易發(fā)胖。
圖1 腰圍后驗(yàn)樣本的直方圖、核密度圖、軌跡圖和自相關(guān)圖
表1 參數(shù)的后驗(yàn)樣本描述
本文討論了醫(yī)學(xué)領(lǐng)域中二分類(lèi)數(shù)據(jù)分析的貝葉斯Probit回歸,像絕大部分貝葉斯模型一樣,貝葉斯Probit回歸參數(shù)的后驗(yàn)分布異常復(fù)雜,需要采用MCMC模擬抽樣,Gibbs抽樣是眾多MCMC算法中最常用的模擬方法。Gibbs抽樣中需要已知參數(shù)的滿(mǎn)條件分布,在此基礎(chǔ)上迭代抽樣生成參數(shù)的Markov鏈,但從Probit回歸的后驗(yàn)分布卻不能得到簡(jiǎn)單并且抽樣方便的滿(mǎn)條件分布,因此執(zhí)行Gibbs抽樣也就不大可行,本文通過(guò)增加潛在變量解決了這個(gè)問(wèn)題。
潛在變量并不能被觀察到,因此在貝葉斯統(tǒng)計(jì)中當(dāng)作未知量看待,則此時(shí)Probit回歸的后驗(yàn)分布為g(β,Y*|Y),是回歸參數(shù)和潛變量的聯(lián)合密度函數(shù)。在這個(gè)后驗(yàn)分布中,β,Y*各自的滿(mǎn)條件分布分別是g(β|Y*,Y)和 g(Y*|β,Y),在 Probit回歸中,前者為多元正態(tài)分布,后者為截尾正態(tài)分布,這兩個(gè)滿(mǎn)條件分布都比較簡(jiǎn)單而且容易直接進(jìn)行模擬抽樣,因此執(zhí)行Gibbs抽樣也就沒(méi)有困難。在給定回歸參數(shù)初始值后在這兩個(gè)滿(mǎn)條件分布之間反復(fù)迭代生成參數(shù)和潛變量的Markov鏈,在算法收斂后則可認(rèn)為生成的參數(shù)序列來(lái)自Probit回歸后驗(yàn)分布。構(gòu)造潛變量的Gibbs抽樣可以看作是一種數(shù)據(jù)擴(kuò)增技術(shù),通過(guò)在模型的后驗(yàn)分布中增加輔助變量使得Gibbs抽樣更加容易。另一個(gè)好處是,生成的潛在變量向量可以進(jìn)一步作為模型診斷的信息加以利用。
1.Salsburg D.The lady tasting tea:how statistics revolutionized science in the twentieth century.Holt Paperbacks,2002.
2.Dobson AJ,Barnett A.An introduction to generalized linear models,third edition.Chapman & Hall,2009.
3.Gelman A,Carlin JB,Stern HS,et al.Bayesian data analysis,2nd ed.London:Chapman & Hall,2004.
4.Albert J.Bayesian computation with R,2nd Ed.New York:Springer,2009.
5.Lynch SM.Introduction to applied bayesian.Statistics and Estimation for Social Scientists,New York:Springer,2009.
6.http://www.r-project.org/