山西醫(yī)科大學(xué)衛(wèi)生統(tǒng)計(jì)教研室（030001） 喬 楠 梁洪川 張海霞 趙俊琴 王 聰 王 彤

半相依Poisson回歸模型及其在衛(wèi)生服務(wù)需求調(diào)查中的應(yīng)用*

山西醫(yī)科大學(xué)衛(wèi)生統(tǒng)計(jì)教研室（030001） 喬 楠 梁洪川 張海霞 趙俊琴 王 聰 王 彤△

半相依回歸（seemingly unrelated regression，SUR）也稱為似乎不相關(guān)回歸，它基于多元回歸模型，在參數(shù)估計(jì)過程中考慮了方程間的相關(guān)性信息，改進(jìn)了方程的估計(jì)效率。1962年Zellner首次將其運(yùn)用在線性模型框架中［1］，1975年Gallant將半相依線性回歸擴(kuò)展到非線性［2］，許多學(xué)者的不斷深入研究使其在經(jīng)濟(jì)、醫(yī)學(xué)、工業(yè)和社會(huì)科學(xué)等多個(gè)領(lǐng)域得到廣泛運(yùn)用和發(fā)展。

在國(guó)內(nèi)醫(yī)學(xué)領(lǐng)域中，關(guān)于SUR已經(jīng)在應(yīng)變量與自變量間為線性或非線性關(guān)系時(shí)進(jìn)行了理論方法的探討和應(yīng)用［3-4］。但在實(shí)際應(yīng)用中有時(shí)會(huì)遇到結(jié)果變量為多個(gè)離散變量計(jì)數(shù)資料的情況，例如在經(jīng)濟(jì)計(jì)量學(xué)研究中，管理者可能會(huì)對(duì)經(jīng)濟(jì)主體同時(shí)作出的多元離散經(jīng)濟(jì)決策感興趣；在醫(yī)療衛(wèi)生領(lǐng)域內(nèi)，醫(yī)務(wù)工作者希望了解某種疾病的多個(gè)離散指標(biāo)的相互關(guān)系及其影響因素。統(tǒng)計(jì)領(lǐng)域內(nèi)對(duì)因變量為多元連續(xù)分布資料的研究由來已久，但對(duì)因變量為多元離散分布資料的應(yīng)用相對(duì)較少。本文針對(duì)多元離散變量計(jì)數(shù)資料，以澳大利亞衛(wèi)生服務(wù)需求問題為例，擬合半相依Poisson回歸模型，闡明半相依廣義線性模型的原理與方法，同時(shí)完成R軟件實(shí)現(xiàn)，為多元離散數(shù)據(jù)分析提供一種參考。

模型原理與方法

1.基本模型

（1）半相依廣義線性模型

半相依廣義線性模型（seemingly unrelated generalized linearmodel）從模型結(jié)構(gòu)上與多元廣義線性模型相比，只在于半相依模型中某些方程的某些解釋變量對(duì)應(yīng)的回歸系數(shù)為0。由此，我們可以考慮對(duì)多元廣義線性模型的參數(shù)進(jìn)行限制，使得對(duì)應(yīng)的回歸系數(shù)為0，從而將多元廣義線性模型推廣到半相依廣義線性模型。

其表示形式，記為

其中yi為（n×1）向量，對(duì)應(yīng)于yi的n個(gè)觀測(cè)值；Z為由n個(gè)觀測(cè)和k個(gè)不同的解釋變量構(gòu)成的（n×k）矩陣；Hi是以0或1為元素的（k×ki）階選擇矩陣（或稱為限制矩陣）；βi為（ki×1）未知參數(shù)向量；ei＝（e1i，e2i，…，eni）為（n×1）隨機(jī)誤差向量。

（2）Poisson回歸模型

Poisson回歸屬?gòu)V義線性模型指數(shù)分布族，其模型連接函數(shù)多采用對(duì)數(shù)函數(shù)。Poisson回歸模型的一般形式可寫為

其中因變量yi服從參數(shù)為λi的Poisson分布：xi表示某一事件觀測(cè)的發(fā)生數(shù)，βi是解釋向量xi對(duì)應(yīng)的回歸系數(shù)。

（3）半相依Poisson回歸模型

在因變量為離散變量的半相依廣義線性回歸中，當(dāng)各方程自變量均相同時(shí)其廣義最小二乘估計(jì)仍優(yōu)于分別對(duì)各方程進(jìn)行估計(jì)的傳統(tǒng)方法，因此我們可首先使用多元Poisson回歸擬合模型，然后剔除參數(shù)為零的自變量，再擬合半相依廣義線性回歸模型。

2.參數(shù)估計(jì)

半相依回歸模型的參數(shù)估計(jì)方法是由多元回歸模型中推廣而來。在多元回歸模型中，給定樣本Y＝（y1，…，yn），我們可以寫出yi的記分函數(shù)

兩者有關(guān)系F（β）＝E（Fobs（β）），當(dāng)聯(lián)接函數(shù)取自然聯(lián)接函數(shù)時(shí)F（β）＝Fobs（β）。令記分函數(shù)s（β）＝0可得參數(shù)的漸近正態(tài)估計(jì)由于此處的似然方程通常都是非線性的，所以方程的解必須借助迭代法。給定某初值Fisher記分迭代法的公式為

對(duì)于半相依廣義線性模型，將模型誤差的方差-協(xié)方差矩陣Wi作為權(quán)重，采用迭代最小二乘估計(jì)方法，參數(shù)估計(jì)的準(zhǔn)則是極小化以下目標(biāo)函數(shù)

采用迭代最小二乘法可使我們很方便地使用線性回歸中的已有程序。在F（β）滿足正定的情況下，通常迭代只需進(jìn)行幾次即可收斂。當(dāng)不收斂的情況出現(xiàn)，一般是由于初值選擇不理想所致，或者是由于在假定的參數(shù)空間中不存在極大似然值。對(duì)于前一種情況可變換初值，多試幾次，或者采用修正迭代法依次取λ＝1，0.9，…，0.5，直到獲得收斂；對(duì)于后一種，則要懷疑模型設(shè)定的合理性了。

實(shí)例分析

數(shù)據(jù)來自于1977-1978年澳大利亞衛(wèi)生服務(wù)調(diào)查資料［7］，為探討衛(wèi)生服務(wù)需求情況，收集了5190位居民的兩周就診數(shù)和兩日處方數(shù)的數(shù)據(jù)，欲分析兩日就診數(shù)和兩日處方數(shù)與性別、年齡、年收入和是否參保的關(guān)系。

通常，兩日就診數(shù)和兩日處方數(shù)可假定服從Poisson分布，所以可以使用Poisson回歸模型擬合數(shù)據(jù)。從邏輯上兩日就診數(shù)和兩日處方數(shù)存在相關(guān)關(guān)系，因此，我們可以使用兩因變量的相關(guān)信息以獲得有效的估計(jì)。使用R軟件的VGAM程序包進(jìn)行分析，模型擬合結(jié)果見表1。

由多元Poisson回歸模型擬合結(jié)果可看到，是否參保對(duì)兩周就診數(shù)差別無影響，而年收入的多少對(duì)兩日處方數(shù)差別無影響。

將是否參保和年收入變量分別從方程1和方程2剔除，則擬合半相依Poisson回歸模型結(jié)果見表2。

從表2擬合結(jié)果可見，與多元Poisson回歸模型擬合結(jié)果相比較，是否參保對(duì)兩日處方數(shù)的影響無統(tǒng)計(jì)學(xué)意義。對(duì)于兩周就診數(shù)來說，女性就診數(shù)高于男性，其OR值為1.241；年齡越高，兩周就診數(shù)越高，每5歲的OR值為3.449；年收入越高，其兩周處方數(shù)反而越低，每1000澳元OR值為0.758。對(duì)于兩日處方數(shù)，女性處方數(shù)高于男性，其OR值為1.811；年齡越高，兩周就診數(shù)越高，每5歲的OR值為21.029。另外不同年收入對(duì)兩日處方數(shù)的影響無統(tǒng)計(jì)學(xué)意義，是否參保對(duì)兩周就診數(shù)和兩日處方數(shù)的影響均無統(tǒng)計(jì)學(xué)意義。這樣我們可以得到最終的預(yù)測(cè)模型。

表1 多元Poisson回歸模型參數(shù)估計(jì)結(jié)果

表2 半相依Poisson回歸模型參數(shù)估計(jì)結(jié)果

表3給出按傳統(tǒng)的處理方法對(duì)各方程分別用Poisson回歸模型擬合的結(jié)果。

從表1～3可見，對(duì)各方程分別用Poisson回歸模型擬合得到的參數(shù)的標(biāo)準(zhǔn)誤均大于多元Poisson回歸和半相依Poisson回歸模型參數(shù)的標(biāo)準(zhǔn)誤，說明這兩種模型的參數(shù)估計(jì)效率高于傳統(tǒng)方法。在本例，半相依Poisson回歸模型與多元Poisson回歸模型的參數(shù)估計(jì)結(jié)果并不一致，一是前者估計(jì)參數(shù)的標(biāo)準(zhǔn)誤小于后者估計(jì)參數(shù)的標(biāo)準(zhǔn)誤；另外，在多元Poisson回歸模型中，是否參保對(duì)兩日處方數(shù)的影響有統(tǒng)計(jì)學(xué)意義，但在半相依Poisson回歸模型，是否參保對(duì)兩日處方數(shù)的影響無統(tǒng)計(jì)學(xué)意義。從估計(jì)參數(shù)的標(biāo)準(zhǔn)誤大小來看，半相依Poisson回歸模型的參數(shù)估計(jì)效率高于多元Poisson回歸模型。本例中是否參保變量由兩種估計(jì)方法得到的95%可信區(qū)間均非常接近于0，其專業(yè)價(jià)值尚需進(jìn)一步研究確認(rèn)。

表3 各方程分別用Poisson回歸模型擬合的參數(shù)估計(jì)結(jié)果

討 論

本文從廣義多元線性模型出發(fā)，討論了當(dāng)因變量為多元分類變量時(shí)模型建模的一般理論。在指數(shù)分布族內(nèi)，導(dǎo)出了廣義多元線性模型的記分函數(shù)，從而得出參數(shù)估計(jì)的迭代最小二乘法；當(dāng)每一方程因變量由不同影響因素決定時(shí)，我們可以通過對(duì)自變量對(duì)應(yīng)參數(shù)施加限制，使其參數(shù)為零得到半相依廣義線性回歸模型。以1977-1978年澳大利亞衛(wèi)生服務(wù)調(diào)查數(shù)據(jù)為例，擬合半相依Poisson回歸模型，說明其擬合的基本過程及優(yōu)點(diǎn)。

由于算法復(fù)雜，目前能夠?qū)崿F(xiàn)半相依廣義線性回歸分析的軟件非常少，一般的軟件只是提供了其中的一種或幾種模型，如Gauss軟件的CML模塊只實(shí)現(xiàn)了半相依泊松回歸模型分析。在由Yee［6］等編寫的R程序包VGAM中，提供了眾多的離散多元因變量分析模型，包括廣義多元線性模型、半相依廣義線性回歸模型、向量廣義可加模型等。

1.Zellner A.An efficientmethod of estimating seem ingly unrelated regressions and tests for aggregation bias.JAm Statist Assoc，1962，57：348-368.

2.Gallant AR.seemingly unrelated nonlinear regressions.Journal of Econometics，1975，3：35-50.

3.梁洪川，韓宏，郎素萍，等.似乎不相關(guān)回歸模型及其在老年認(rèn)知問題中的應(yīng)用.中國(guó)衛(wèi)生統(tǒng)計(jì)，2005，22（6）：362-364.

4.趙俊康，梁洪川.非線性半相依回歸模型在生長(zhǎng)曲線研究中的應(yīng)用.中國(guó)衛(wèi)生統(tǒng)計(jì)，2012，29（3）：348-350.

5.Cameron A，Trivedi P.Regression Analysis of Count Data.Oxford University Press，1998.

6.Yee TW.VGAM：Vector Generalized Linear and Additive Models.R package version0.6-7.http：//www.stat.auckland.ac.nz/～yee，2005.

（責(zé)任編輯：劉 壯）

*：國(guó)家自然科學(xué)基金項(xiàng)目（81072385）；全國(guó)統(tǒng)計(jì)科研計(jì)劃重點(diǎn)項(xiàng)目（2009LZ033）

△通信作者：王彤