張春霞 蔣紅衛(wèi) 尹 平

傳統(tǒng)的統(tǒng)計(jì)分析方法通常是建立在獨(dú)立觀測值假定基礎(chǔ)上的。然而，在遇到空間數(shù)據(jù)時，獨(dú)立觀測值在現(xiàn)實(shí)研究中并不是普遍存在的〔1〕。對于具有地理空間屬性的數(shù)據(jù)，一般認(rèn)為在空間上，離得近的變量之間比離得遠(yuǎn)的變量之間具有更加密切的關(guān)系〔2〕。隨著遙感技術(shù)、地理信息系統(tǒng)和全球定位系統(tǒng)滲透到公共衛(wèi)生研究各個領(lǐng)域之中，空間測量技術(shù)在衛(wèi)生研究中得到了廣泛的應(yīng)用，空間統(tǒng)計(jì)學(xué)正在逐漸成為衛(wèi)生研究領(lǐng)域的熱點(diǎn)〔3〕。

由于空間相關(guān)性的存在，傳統(tǒng)的計(jì)數(shù)回歸模型難以較好地處理具有空間聚集性的數(shù)據(jù)，經(jīng)常出現(xiàn)模型擬合精度不高，過度離散或離散過小等情形。較典型的是，作為血吸蟲唯一中間宿主——釘螺的空間分布格式研究中，除了草叢密度，溫度等因素的影響之外，釘螺的分布往往具有明顯的空間聚集特性。蘇德隆早在1963年提出水網(wǎng)型釘螺空間服從負(fù)二項(xiàng)分布規(guī)律〔4〕。但是，也有研究者認(rèn)為釘螺空間分布并不服從單一的負(fù)二項(xiàng)分布〔5－7〕。由于存在大量因素影響釘螺的空間分布，僅僅通過強(qiáng)行擬合負(fù)二項(xiàng)分布，無法較好地解決釘螺空間分布與預(yù)測等問題，因而，為了解決此類存在空間相關(guān)性的數(shù)據(jù)分析問題，必須引入空間計(jì)數(shù)回歸模型。本文將在引入空間數(shù)據(jù)形式的基礎(chǔ)上，系統(tǒng)地研究空間負(fù)二項(xiàng)回歸模型及其參數(shù)估計(jì)，并利用釘螺小范圍的空間分布實(shí)例研究以說明空間負(fù)二項(xiàng)回歸模型。

空間數(shù)據(jù)基本結(jié)構(gòu)

空間數(shù)據(jù)主要包括三個部分:

(1)空間位置:s1，s2，…，sn;

(2)反應(yīng)變量:Y(s1)，Y(s2)，…，Y(sn);

(3)解釋變量:X(s1)，X(s2)，…，X(sn)。

其中，空間位置通常為二維空間域，si=(，)∈D?R2，與分別為空間位置的橫坐標(biāo)與縱坐標(biāo)，可

其條件均值為E(Y(s)|b(s))=μ(s)，方差V(Y(s)|b(s))= μ(s)+μ(s)2/k。

(3)連接函數(shù)

(4)回歸模型形式以為離散型變量，也可以為連續(xù)型變量;反應(yīng)變量Y與p維解釋變量X={X1，X2，…，Xp}對應(yīng)于每一個空間位置，表示區(qū)域化的反應(yīng)變量與解釋變量。

空間負(fù)二項(xiàng)回歸模型

空間負(fù)二項(xiàng)回歸給定以下四個基本假定:(1)反應(yīng)變量服從負(fù)二項(xiàng)分布;(2)反應(yīng)變量的期望值通過自然對數(shù)變換后，與各解釋變量之間呈線性關(guān)系;(3)回歸系數(shù)在所研究的空間區(qū)域內(nèi)具有一致性(即為常數(shù));(4)空間隨機(jī)部分呈高斯平穩(wěn)過程。根據(jù)這四個基本假定，構(gòu)建如下空間負(fù)二項(xiàng)回歸基本模型，用于反映空間數(shù)據(jù)的特征與變化規(guī)律。

1．回歸模型形式

假定{b(si)，si∈R2}，i=1，2，…，n 為潛在的不可觀測的空間隨機(jī)過程。其中，b(si)表示空間位置Si上無法被解釋變量所能解釋的隨機(jī)部分?？臻g負(fù)二項(xiàng)回歸模型可定義如下:

(1){b(s)，s∈R2}是一個高斯平穩(wěn)過程，具有期望E(b(s))=0，方差 －協(xié)方差函數(shù) cov(b(s+h)，b(s))=C(h)，其中，s，h∈R2，方差函數(shù) C(h)=C(θ)，θ∈Rv，既是距離的函數(shù)，又是方向的函數(shù)，但對于高斯平穩(wěn)過程，其必須滿足二階平穩(wěn)性與各向同性。

(2)每個位置上的反應(yīng)變量{Y(s)，s∈R2}，條件依賴于高斯平穩(wěn)過程{b(s)，s∈R2}，即反應(yīng)變量的概率密度函數(shù)為

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

采用最大似然法來估計(jì)空間負(fù)二項(xiàng)回歸模型的參數(shù)。假定高斯平穩(wěn)過程依賴于正態(tài)各向同性協(xié)方差函數(shù)

根據(jù)空間負(fù)二項(xiàng)回歸模型的設(shè)定，其似然函數(shù)為

將相應(yīng)的分布函數(shù)代入公式(5)中，可以獲取似然函數(shù)最大時各參數(shù)的估計(jì)。由于上式相當(dāng)復(fù)雜，是高維的非線性積分函數(shù)，直接求解極為困難，因此，采用蒙特卡羅最大似然法(MCML)〔8－9〕來進(jìn)行參數(shù)估計(jì)。在蒙特卡羅最大似然法中，常采用Gibbs抽樣的方法，具體計(jì)算可在WinBUGS軟件下進(jìn)行。按照空間負(fù)二項(xiàng)回歸模型編寫程序(見附錄)，給定待估計(jì)參數(shù)的先驗(yàn)(理論上，通常應(yīng)采用無信息先驗(yàn)，但鑒于winBUGS中dflat()函數(shù)易導(dǎo)致復(fù)雜模型難以收斂的問題，各參數(shù)的具體設(shè)定見附錄中的 winBUGS程序)，就可以對待估計(jì)參數(shù)及其可信區(qū)間等統(tǒng)計(jì)指標(biāo)加以估計(jì)。首先，通過SAS9獲得常規(guī)負(fù)二項(xiàng)回歸模型回歸系數(shù)的估計(jì)值μβi及其標(biāo)準(zhǔn)誤σβi;其次，令回歸模型參數(shù)初始值服從正態(tài)分布N(μβi，100*σ2βi)，隨機(jī)取值，獲得50組回歸系數(shù)初始值;其次，采用Gibbs抽樣法，迭代20000次，其中前5000次作為預(yù)熱，對后15000次作為待估參數(shù)后驗(yàn)分布的模擬，使用50組參數(shù)初始值，從而，獲得50條MCMC鏈，利用方差比法，判斷MCMC鏈?zhǔn)欠袷諗?最后，合并各鏈后15000次迭代數(shù)據(jù)，得到空間負(fù)二項(xiàng)回歸模型的參數(shù)估計(jì)值及其標(biāo)準(zhǔn)誤。

應(yīng)用實(shí)例

血吸蟲的終宿主是廣譜的，幾乎所有哺乳類動物都可以是;而中間寄主釘螺則是唯一的。了解釘螺特性，快速獲得釘螺的地理空間分布信息，并采取有效的滅螺措施是消滅血吸蟲病的一項(xiàng)重要舉措。本實(shí)例的釘螺樣本數(shù)據(jù)來自某市2009年4月一個小區(qū)域釘螺分布調(diào)查的部分。具體調(diào)查方案是，在洲灘上縱深設(shè)線設(shè)點(diǎn)調(diào)查，每個樣本點(diǎn)的尺寸大小是0．33米×0．33米予以制定，每個樣點(diǎn)按照水線走向設(shè)定位置坐標(biāo)，對樣本點(diǎn)中的活釘螺予以計(jì)數(shù)，并記錄與之相關(guān)的草叢密度，分為高密度(記為0)與低密度(記為1)兩類。原始數(shù)據(jù)見表1。

表1 釘螺調(diào)查數(shù)據(jù)

釘螺分布的均數(shù)為2.15，標(biāo)準(zhǔn)差為3.42，中位數(shù)為1，四分位數(shù)間距為2，可見釘螺的空間分布呈現(xiàn)出明顯的正偏態(tài)分布，且可能存在離散過大趨勢。Moran's I空間相依指數(shù)為0.188，P＜0.0001，表明釘螺的空間分布具有一定的空間相關(guān)性。為了更直觀地研究釘螺的空間相關(guān)性，對其作半方差圖，見圖1。

圖1 原始釘螺分布的半方差圖

從圖1可以看出，原始釘螺空間分布具有一定的空間聚集性。塊金值約為3.5，基臺約為16，變程約為8，表明釘螺空間分布的變異較大，具有各向同性，并且存在著較長程的空間相關(guān)性。

為了更好地分析小范圍釘螺分布的規(guī)律，擬合了常規(guī)的負(fù)二項(xiàng)回歸〔11〕以及本文介紹的空間負(fù)二項(xiàng)回歸，其參數(shù)估計(jì)與模型擬合見表2。

表2 常規(guī)負(fù)二項(xiàng)回歸與空間負(fù)二項(xiàng)回歸模型擬合與參數(shù)估計(jì)

由表2可以看出，無論擬合常規(guī)負(fù)二項(xiàng)回歸，還是擬合空間負(fù)二項(xiàng)回歸，各回歸參數(shù)估計(jì)均具有統(tǒng)計(jì)學(xué)意義。這也就是說，無論是否考慮釘螺的空間相關(guān)性，草叢密度對釘螺分布均是具有一定的影響，提示其與釘螺的發(fā)生與發(fā)展具有較強(qiáng)的關(guān)聯(lián)性。而從專業(yè)上講，釘螺分布具有明顯的空間聚集性，必須考慮其空間分布格局與聚集性特征。由表2可見，由于空間負(fù)二項(xiàng)回歸考慮了釘螺的空間相關(guān)性，因而，從模型擬合的角度來看，相對于常規(guī)負(fù)二項(xiàng)回歸而言，空間負(fù)二項(xiàng)回歸的空模型與含有草叢密度變量模型的-2log likelihood統(tǒng)計(jì)量均有所下降。這表明無論是空模型，還是考慮草叢密度的模型，空間負(fù)二項(xiàng)回歸模型擬合精度均有較明顯的提高。由此可見，空間相關(guān)性對負(fù)二項(xiàng)回歸模型存在著較強(qiáng)的影響。另外，從參數(shù)估計(jì)的角度來看，與常規(guī)負(fù)二項(xiàng)回歸相比較，空間負(fù)二項(xiàng)回歸模型截距與回歸系數(shù)的估計(jì)值相對較小，而其標(biāo)準(zhǔn)誤相對較大。這提示如果不考慮釘螺的空間相關(guān)性特征，那么，有可能會導(dǎo)致I型錯誤增大，使得負(fù)二項(xiàng)回歸模型參數(shù)估計(jì)失效。

由于具有空間聚集性的解釋變量納入到回歸模型，有可能會影響到反應(yīng)變量空間相關(guān)性的分析結(jié)果，故而，必須在排除這部分解釋變量的影響之后，重新審視反應(yīng)變量的空間相關(guān)性。本實(shí)例中，從專業(yè)上講，草叢也具有一定的空間聚集性，必須通過擬合含有草叢密度的空間負(fù)二項(xiàng)回歸模型，排除草叢密度的影響，再對釘螺的空間相關(guān)性進(jìn)行分析。擬合空間負(fù)二項(xiàng)回歸模型后，其模型殘差的 Moran's I空間相依指數(shù)為0.078，P＜0.0001。與擬合空間負(fù)二項(xiàng)回歸模型前的Moran's I空間相依指數(shù)相比較，該指數(shù)數(shù)值有所下降，這說明排除草叢密度的影響之后，釘螺的空間相關(guān)性有所降低。該統(tǒng)計(jì)量假設(shè)檢驗(yàn)仍然具有統(tǒng)計(jì)學(xué)意義，這表明釘螺的空間分布仍然具有一定的空間相關(guān)性。見圖2。

圖2 擬合空間負(fù)二項(xiàng)回歸模型排除草叢密度之后的釘螺分布的半方差圖

與圖1相比較，圖2與圖1存在著較強(qiáng)的相似性，半方差曲線走勢大致相同，變程也約為8，只是塊金值與基臺值較小，分別為1.0和1.9。這表明擬合含有草叢密度的空間負(fù)二項(xiàng)回歸模型，可以有效地降低釘螺空間分布的變異性，對釘螺的空間分布格局與規(guī)律沒有產(chǎn)生較大影響，使得其空間相關(guān)性的基本趨勢保持不變。

結(jié) 論

綜上所述，空間負(fù)二項(xiàng)回歸模型考慮了空間相關(guān)性的影響，是對常規(guī)的負(fù)二項(xiàng)回歸模型的一種改進(jìn)，可用于分析具有空間相關(guān)性的數(shù)據(jù)，能有效地降低空間相關(guān)性對回歸系數(shù)估計(jì)的影響，更好地分析空間分布格局與規(guī)律。在實(shí)際應(yīng)用中，如若數(shù)據(jù)存在著空間相關(guān)性時，必須將空間相關(guān)性納入到統(tǒng)計(jì)分析考量之中，提高回歸模型的擬合精度與效能。

1．空間相關(guān)性的處理

首先，對原始數(shù)據(jù)進(jìn)行有關(guān)的空間統(tǒng)計(jì)描述，如，Moran's I空間相依指數(shù)，半方差圖等，可以較為有效地發(fā)現(xiàn)數(shù)據(jù)可能存在的空間相關(guān)性的程度與模式進(jìn)行初步分析。其次，在空間統(tǒng)計(jì)描述的基礎(chǔ)上，結(jié)合資料的具體特點(diǎn)，擬合相應(yīng)的空間回歸模型，以提高回歸模型的擬合精度，并保證回歸系數(shù)估計(jì)的有效性。最后，在排除解釋變量(可能具有空間聚集性的解釋變量)的影響之后，通過殘差來分析反應(yīng)變量的空間聚集性與分布格局。

2．空間負(fù)二項(xiàng)回歸模型的進(jìn)一步改進(jìn)

本文實(shí)例研究表明，空間負(fù)二項(xiàng)回歸模型是一種可以用于分析具有空間相關(guān)性計(jì)數(shù)型反應(yīng)變量的廣義線性模型，可以處理因空間異質(zhì)而導(dǎo)致的空間相關(guān)性問題。值得注意的是，該回歸模型缺乏較好地處理因反應(yīng)變量自相關(guān)性而導(dǎo)致的空間相關(guān)性問題。因此，空間負(fù)二項(xiàng)回歸模型研究工作尚需從以下四個方面進(jìn)行展開:(1)進(jìn)一步探討空間負(fù)二項(xiàng)回歸模型的理論性質(zhì);(2)將反應(yīng)變量的滯后變量作為解釋變量納入到回歸模型之中;(3)逐步放松空間負(fù)二項(xiàng)回歸的某些基本假定，擴(kuò)展其適用范圍;(4)對空間負(fù)二項(xiàng)回歸模型進(jìn)行統(tǒng)計(jì)模擬研究，探討其統(tǒng)計(jì)效能，全面地了解空間負(fù)二項(xiàng)回歸模型的特性。根據(jù)這些方向，將在以后的研究工作中對空間負(fù)二項(xiàng)回歸進(jìn)行更為系統(tǒng)的模型改進(jìn)與理論研究。

1．Getis A，Mur J，Zoller GH．Spatial Econometrics and Spatial Statistics．New York，Palgrave Macmillan，2004．

2．Anselin L，Getis A．Spatial statistical analysis and geographic information systems．The annuals of Regional Science，1992，26:19-33．

3．Lawson BL．Statistical Methods in Spatial Epidemiology，2nd Edition．Chichester，UK:Wiley，2006．

4．蘇德?。K德隆教授論文選集．天津:天津科學(xué)技術(shù)出版社，1995:97-101．

5．陳峰，楊樹勤．疾病的復(fù)合分布模型研究-疾病的統(tǒng)計(jì)分布(一)．中國衛(wèi)生統(tǒng)計(jì)，1995，12(6):12-15．

6．趙安，鮑曙明．知識驅(qū)動的鄱陽湖區(qū)域釘螺分布預(yù)測模型研究初探．科學(xué)通訊，2007，52(22):2663-2670．

7．張志杰，彭文祥，ong Senghuat．廣義負(fù)二項(xiàng)分布對釘螺分布的擬合．中國衛(wèi)生統(tǒng)計(jì)，2008，25(1):2-6．

8．Geyer CJ．Markov chain Monte Carlo maximum likelihood．In Computing Science and Statistics:Proceedings of the 23rd Symposium on the Interface(ed．E．M．Keramidas)．Interface Foundation of North America，F(xiàn)airfax Station，VA．，1991:156-163．

9．Geyer CJ，Thompson EA．Constrained Monte Carlo maximum likelihood for dependent data(with Discussion)．Journal of the Royal Statistical Society B，1992，54:657-699．

10．Hilbe MJ．Negative Binomial Regression．Cambridge University Press，Cambridge，U．K．，2010．