貝葉斯空間計(jì)量模型綜述：理論及其應(yīng)用

李麗輝，溫 濤，劉 達(dá)

(西南大學(xué) a.數(shù)學(xué)與統(tǒng)計(jì)學(xué)院；b.經(jīng)濟(jì)管理學(xué)院，重慶 400715)

一、引 言

空間計(jì)量模型起源于區(qū)域科學(xué)模型和計(jì)量經(jīng)濟(jì)模型的共同發(fā)展，研究如何在截面數(shù)據(jù)和面板數(shù)據(jù)模型中處理空間相互作用和空間結(jié)構(gòu)，是計(jì)量經(jīng)濟(jì)模型的一個(gè)分支?？傮w而言，空間計(jì)量模型的推斷理論主要包含估計(jì)理論、檢驗(yàn)理論及模型選擇理論三個(gè)方面。LeSage將一種非常特別且在某些領(lǐng)域極為有用的統(tǒng)計(jì)學(xué)方法——貝葉斯統(tǒng)計(jì)方法，運(yùn)用到空間計(jì)量建模中，建立了貝葉斯空間計(jì)量模型，突破了空間計(jì)量模型經(jīng)典方法的限制，適用面廣，在學(xué)術(shù)界產(chǎn)生了極大的影響，也是當(dāng)前空間計(jì)量模型理論和應(yīng)用領(lǐng)域的最新熱點(diǎn)[1]。貝葉斯空間計(jì)量模型的最新發(fā)展包括從截面到面板、從靜態(tài)到動(dòng)態(tài)、從解析到仿真，其應(yīng)用范圍也涉及經(jīng)濟(jì)學(xué)各個(gè)領(lǐng)域。

本文試圖通過對(duì)一些相關(guān)經(jīng)典文獻(xiàn)著作的簡單介紹及梳理，為理論工作者提供一定理論參考，為實(shí)證分析研究者提供一種分析類似問題的方法。首先，對(duì)空間計(jì)量模型的分類系統(tǒng)進(jìn)行梳理；其次，介紹貝葉斯空間計(jì)量模型的推斷理論及其應(yīng)用；再次，介紹貝葉斯空間計(jì)量模型選擇理論及其應(yīng)用；最后，指明貝葉斯空間計(jì)量模型的發(fā)展方向。

二、空間計(jì)量模型分類

空間計(jì)量經(jīng)濟(jì)學(xué)起源于區(qū)域科學(xué)和計(jì)量經(jīng)濟(jì)學(xué)的共同發(fā)展，研究如何在截面數(shù)據(jù)和面板數(shù)據(jù)中處理空間相互作用和空間結(jié)構(gòu)，是計(jì)量經(jīng)濟(jì)學(xué)的一個(gè)分支。最早由J.Paelinck在1972年荷蘭統(tǒng)計(jì)協(xié)會(huì)年會(huì)大會(huì)致詞時(shí)提出，后經(jīng)Anselin等人發(fā)展，最終形成了學(xué)科框架體系。自Krugman在“規(guī)模報(bào)酬與經(jīng)濟(jì)地理”中建立了“中心―外圍”模型并將空間因素引入主流經(jīng)濟(jì)學(xué)分析框架后，經(jīng)濟(jì)學(xué)家對(duì)空間因素的重視程度與日俱增，尤其是2008年Krugman因其“全面整合國際貿(mào)易與經(jīng)濟(jì)地理的研究成果”而獲得諾貝爾經(jīng)濟(jì)學(xué)獎(jiǎng)，使空間計(jì)量經(jīng)濟(jì)學(xué)成為舉世矚目的主流經(jīng)濟(jì)學(xué)?？臻g計(jì)量經(jīng)濟(jì)學(xué)被廣泛運(yùn)用至各個(gè)社會(huì)科學(xué)方面，包括社會(huì)學(xué)、犯罪學(xué)、政治學(xué)、經(jīng)濟(jì)學(xué)等。

在經(jīng)濟(jì)學(xué)領(lǐng)域，涉及的內(nèi)容包括空間溢出、城市發(fā)展和組群經(jīng)濟(jì)、貿(mào)易和經(jīng)濟(jì)增長等，其基本思想是對(duì)于某一區(qū)域的觀察值yi，由于其依賴于其鄰接區(qū)域某些變量值的加權(quán)和加上一些固定效應(yīng)或隨機(jī)噪聲，因此就用自回歸的方式模擬這樣的空間交互效應(yīng)，這使得模型能夠清楚地展示出空間依賴性并且在模型的誤差項(xiàng)中呈現(xiàn)出一種特殊的方差―協(xié)方差結(jié)構(gòu)。

空間計(jì)量模型按照數(shù)據(jù)分布類型主要分為以下幾種模型：

(一)高斯模型

這類模型以最常見的聯(lián)立自回歸模型(SAR)模型或者空間誤差模型(SEM)為代表，其基本形式可表述為：

y=Xβ+μ；μ=ρWμ+ε；ε～MVN(0,σ2In)

(1)

其中y=(y1,y2,…,yn)′為觀測值向量，X為一個(gè)n×p階設(shè)計(jì)矩陣，β=(β1,β2,…,βp)′為變量系數(shù)，In為n維單位矩陣，W為一個(gè)行標(biāo)準(zhǔn)化后的鄰接矩陣，ρ為空間自回歸參數(shù)。隨機(jī)誤差項(xiàng)被模擬為其鄰接區(qū)域的隨機(jī)誤差的加權(quán)和再加上一個(gè)隨機(jī)噪聲ε，其中隨機(jī)噪聲的分布可用多元標(biāo)準(zhǔn)正態(tài)分布(MVN)來刻畫。LeSage和Pace探討了誤差項(xiàng)的異方差結(jié)構(gòu)，在他們的推導(dǎo)下，SEM模型可表述為[2]25-41：

y=Xβ+ε′；

ε′～MVN(0,σ2(In-ρW)-1(In-ρW′)-1)

(2)

上述模型就是一個(gè)在誤差項(xiàng)上有著非對(duì)角線方差-協(xié)方差矩陣的一般的線性回歸模型。

空間計(jì)量模型中另外一種重要的模型形式就是空間滯后模型(SLM)，LeSage和Pace稱之為SAR模型。這種模型中的因變量被模擬為因變量的加權(quán)和再加上自變量的線性項(xiàng)與一個(gè)隨機(jī)誤差項(xiàng)，模型可表述為：

y=ρWy+Xβ+ε；ε～MVN(0,σ2In)

(3)

還有一種非常常見且應(yīng)用廣泛的空間計(jì)量模型是空間杜賓模型(SDM)，可表述為：

y=ρWy+Xβ+WXγ+ε；ε～MVN(0,σ2In)

(4)

其中γ是空間滯后自變量WX項(xiàng)的系數(shù)向量，即因變量不僅依賴于鄰接因變量的加權(quán)和，而且依賴于鄰接自變量的加權(quán)和。

這些模型共同的特點(diǎn)就是誤差項(xiàng)的設(shè)置都是高斯設(shè)置，即均值為0且方差―協(xié)方差矩陣結(jié)構(gòu)為SAR模型所設(shè)定的結(jié)構(gòu)。SLM和SDM模型在自變量的線性項(xiàng)上有更復(fù)雜的結(jié)構(gòu)，空間自相關(guān)參數(shù)ρ被限制在一定區(qū)域上，W作了行標(biāo)準(zhǔn)化。

(二)非高斯模型

當(dāng)響應(yīng)變量為二元離散選擇變量時(shí)，高斯模型不再適用。LeSage等研究了在2005年卡特里娜颶風(fēng)之后，新奧爾良的一位商人考察周邊同行是否重新營業(yè)，進(jìn)而決定自己是否重新開張的問題[3]。他們模擬響應(yīng)變量yi的形式為：

(5)

其中y*為潛變量，測量了潛在的凈收益(如果其大于0，則會(huì)重新開張)。潛變量可用SLM模型來模擬，即：

y*=ρWy*+Xβ+ε；ε～MVN(0,In)

(6)

(三)帶隨機(jī)效應(yīng)的廣義線性模型

Nelder和Wedderburn提出了廣義線性模型(GLM)，該模型能夠較好地模擬協(xié)變量與空間自相關(guān)關(guān)系[4]。在這種模型中，響應(yīng)變量yi設(shè)定為來自于指數(shù)族分布，均值為μi。自變量向量Xi的線性組合部分和μi之間的關(guān)系可通過下述連接函數(shù)g(·)建立：

g(μi)=ηi=Xiβ

(7)

如果引入相關(guān)隨機(jī)效應(yīng)ui，則式(7)變?yōu)椋?/p>

ηi=Xiβ+ui

(8)

u=(u1，u2，…，un)′為多元正態(tài)分布，均值為0，方差協(xié)方差矩陣為Σ，對(duì)矩陣Σ的不同設(shè)定就會(huì)產(chǎn)生不同的空間依賴性。比如在SAR模型設(shè)置中：

Σ=σ2(In-ρW)-1(In-ρW′)-1

(9)

對(duì)于這一類模型來說，其模型擬合和解釋是更容易的，這是因?yàn)檫@類模型的空間相關(guān)參數(shù)和自變量的效應(yīng)項(xiàng)是分開的。比如對(duì)于空間滯后模型，自變量部分可寫為：

(In-ρW)-1Xβ

(10)

可見X的效應(yīng)依賴于系數(shù)β和空間相關(guān)參數(shù)ρ。

第一類高斯模型和第二類非高斯模型都可以視為廣義線性模型。高斯模型就是限定ρ以后帶隨機(jī)效應(yīng)的高斯GLM，空間概率模型就是給定ρ的帶概率連接函數(shù)的二項(xiàng)GLM。

三、貝葉斯空間計(jì)量模型推斷理論及其應(yīng)用

貝葉斯統(tǒng)計(jì)方法在空間計(jì)量分析中的應(yīng)用可視為空間計(jì)量模型領(lǐng)域的一個(gè)前沿分支，由LeSage率先提出[1]。與基于經(jīng)典統(tǒng)計(jì)理論的空間計(jì)量模型相比，用貝葉斯方法進(jìn)行空間計(jì)量分析具有獨(dú)到的優(yōu)勢：一是該方法更為普遍，能運(yùn)用到更加廣泛的統(tǒng)計(jì)領(lǐng)域；二是該方法允許合理地利用先驗(yàn)信息，更利于直接明確地分析具體問題；三是該方法得到的結(jié)果不僅是一個(gè)預(yù)測值，而是一個(gè)完整的未來經(jīng)濟(jì)結(jié)果的概率分布，這種預(yù)測方法比其它預(yù)測方法產(chǎn)生的結(jié)果更加有用；四是相對(duì)于傳統(tǒng)統(tǒng)計(jì)方法，該方法能更加明確地合理處理不確定因素的問題。空間計(jì)量模型按照誤差項(xiàng)的分布主要可分為高斯模型、非高斯模型和廣義線性模型，而貝葉斯方法在空間計(jì)量模型中的應(yīng)用主要涉及模型的推斷與模型比較兩個(gè)方面，本部分將首先介紹貝葉斯方法應(yīng)用于這幾種空間模型的相關(guān)推斷理論及其應(yīng)用。

(一)高斯模型

貝葉斯方法的一個(gè)特征是在建模過程中引入先驗(yàn)信息，且與數(shù)據(jù)分布相結(jié)合以生成用于推斷的后驗(yàn)分布?；镜呢惾~斯方法會(huì)在后驗(yàn)分布中創(chuàng)造一個(gè)關(guān)于樣本均值和先驗(yàn)信息的加權(quán)矩陣，但權(quán)重會(huì)受到與先驗(yàn)信息相關(guān)的樣本數(shù)據(jù)可得信息數(shù)量的強(qiáng)烈影響，在大樣本情形下，這能對(duì)后驗(yàn)分布的解析提供一個(gè)簡化的方法。

對(duì)于高斯模型，LeSage提出了兩類處理的方法：解析法與模擬仿真算法[1]。下面以模型(3)為例進(jìn)行闡述：

1.解析法

對(duì)于模型(3)，其似然概率為：

(11)

其中A=(In-ρW)代表矩陣的行列式。為了產(chǎn)生模型參數(shù)的后驗(yàn)分布，需要指定參數(shù)的先驗(yàn)分布。LeSage采用正態(tài)逆伽瑪(NIG)分布作為參數(shù)β和σ2的先驗(yàn)分布[1]；參數(shù)ρ的可行范圍(1/λmin，1/λmax)，其中λmin和λmax分別代表空間權(quán)重矩陣的最小、最大特征根。若給定σ、β服從多元正態(tài)分布N(c,σ2T)，σ的邊緣分布服從給定(a,b)逆伽瑪分布，則模型(3)可寫為：

y=ρWy+Xβ+ε；ε～MVN(0,σ2In)

(12)

π(β,σ2)～NIG(c,T,a,b)=π(β|σ2)π(σ2)

(13)

exp[-{(β-c)′T-1(β-c)+2b}/(2σ2)]

(14)

其中

σ2>0，a，b>0

用于識(shí)別先驗(yàn)認(rèn)知的參數(shù)就是NIG先驗(yàn)分布中的c、T、a和b等參數(shù)。如果在模型中使用了有信息先驗(yàn)分布，問題就會(huì)變得復(fù)雜，這主要體現(xiàn)在：一是需要指定或分派NIG先驗(yàn)分布中參數(shù)的數(shù)值；二是后驗(yàn)分布難以分析，因?yàn)樾枰獙?duì)參數(shù)ρ和σ2進(jìn)行積分以獲得參數(shù)β的后驗(yàn)表達(dá)式；同理，也需對(duì)β和ρ積分以得到σ2的后驗(yàn)分布。Hepple根據(jù)無信息先驗(yàn)分布得到的簡化表達(dá)式去尋求貝葉斯SAR模型的解，為解決這類問題提供了一個(gè)很好的思路[5]。Hepple基于一個(gè)無信息先驗(yàn)分布來替代上文中所提到的NIG先驗(yàn)分布，并假定ρ與β、σ2這些參數(shù)相互獨(dú)立，這樣就會(huì)得到一個(gè)簡化后的后驗(yàn)分布如式(15)，因此對(duì)于β，設(shè)置其先驗(yàn)為無信息先驗(yàn)，設(shè)定c=0并為β指定一個(gè)很大的先驗(yàn)方差T=IK·1010，且β向量參數(shù)間協(xié)方差為0；通過設(shè)定a=b=0，可以得到σ2的無信息先驗(yàn)；對(duì)于參數(shù)ρ，可以指定一個(gè)預(yù)先值，這意味著所有在可行范圍(1/λmin，1/λmax) 內(nèi)的結(jié)果都是可能的。

根據(jù)貝葉斯理論，模型參數(shù)的后驗(yàn)分布形式為：

p(β,σ,ρ|D)∝p(D|β,σ,ρ)π(β,σ)π(ρ)

(Ay-Xβ)]π(ρ)

(15)

根據(jù)這個(gè)簡化的后驗(yàn)分布，可將σ作為冗余參數(shù)處理并將其從式(15)中積分去掉，由此可得：

p(β,σ,ρ|D)∝|A|{(n-k)s2+

(16)

其中

c*=(X′X)-1X′Ay

s2=(Ay-Xc*)′(Ay-Xc*)/(n-k)

取ρ的條件概率，式(16)就變?yōu)橐粋€(gè)多維t分布，這樣就能對(duì)β積分，得到ρ的邊緣后驗(yàn)分布為：

(17)

只要已知式(17)中ρ的先驗(yàn)分布，就可按照單變量數(shù)值積分方法得到ρ的后驗(yàn)分布，據(jù)此做出統(tǒng)計(jì)推斷。其先驗(yàn)分布的決定形式一般有三種方法： 一是根據(jù)ρ在理論上的可行范圍(1/λmin，1/λmax)來決定；二是可以給參數(shù)ρ指定一個(gè)先驗(yàn)范圍(-1,1)，LeSage和Parent引入了ρ的另一種先驗(yàn)分布Be(α,β)，對(duì)均勻先驗(yàn)U(-1,1)是一個(gè)很好的替代[6]；三是考慮由于負(fù)的空間依賴對(duì)特定問題沒有多大意義，可將ρ的先驗(yàn)定義在區(qū)間(0,1)上。

如果ρ的先驗(yàn)分布一旦給定，就可做出關(guān)于ρ的簡化似然函數(shù)，按照Pace和Barry提出的方法，可通過使用一個(gè)ρ的取值在[ρmin，ρmax]之間的q×1階向量估計(jì)對(duì)數(shù)似然值[7]。只要提供出一個(gè)充分細(xì)的柵格q值并對(duì)其進(jìn)行向量化，則通過內(nèi)推法就可使插值點(diǎn)到達(dá)任何需要的精度。該向量可以通過運(yùn)用Simpson法則等簡單方法進(jìn)行單變量的數(shù)值積分。

2.MCMC仿真模擬

解析法是直接對(duì)參數(shù)后驗(yàn)分布的解析形式進(jìn)行統(tǒng)計(jì)分析，計(jì)算繁瑣，LeSage采用了另一種目前廣泛使用的方法——馬爾科夫鏈蒙特卡羅方法(MCMC)，對(duì)模型后驗(yàn)分布進(jìn)行推斷[1]。

MCMC抽樣算法是基于Markov鏈對(duì)變量進(jìn)行抽樣模擬，由于生成于同一條Markov鏈，MCMC抽樣算法產(chǎn)生的是相依樣本。當(dāng)構(gòu)造的Markov鏈?zhǔn)諗恐聊繕?biāo)分布，即達(dá)到穩(wěn)態(tài)，則此分布為所求的后驗(yàn)分布，然后根據(jù)穩(wěn)態(tài)后驗(yàn)分布中的樣本點(diǎn)計(jì)算Monte Carlo積分。因此，MCMC抽樣算法的關(guān)鍵是從基于穩(wěn)態(tài)后驗(yàn)分布的概率密度函數(shù)中進(jìn)行抽樣。MCMC可以通過檢驗(yàn)一個(gè)來自后驗(yàn)分布的大的隨機(jī)抽樣來實(shí)現(xiàn)，如果來自后驗(yàn)分布p(θ|D)的樣本足夠大，可以運(yùn)用核密度估計(jì)或者直方圖獲得概率密度的近似形式，而不需要精確的密度函數(shù)解析式，其中兩種最為廣泛使用的方法是Metropolis-Hastings(M-H)抽樣算法和Gibbs抽樣算法。

M-H算法依賴于建議分布，對(duì)于初始值θ0，在t+1時(shí)分布為f(θ|θt)，其中θt已知。候選點(diǎn)θ*從建議分布中抽樣獲得，并且：

如果接受概率為：

(18)

則接受θt+1=θ*；否則θt+1=θt。

c*=(X′X+T-1)-1(X′(In-ρ0W)y+T-1c)

T*=(X′X+T-1)-1

(19)

將抽樣參數(shù)向量記為β1并代替參數(shù)向量β0。

2)利用逆伽瑪分布IG(a*,b*)抽樣p(σ2|β1,ρ0)：

p(σ2|β1,ρ0)～I(xiàn)G(a*,b*)

a*=a+n/2

b*=b+(Ay-Xβ1)′(Ay-Xβ1)/2

A=In-ρ0W

(20)

3)利用接受概率為式(18)的M-H算法，將得到的值記為ρ1并返回到第一步。

從第一步到第三步的一個(gè)循環(huán)就完成了一次抽樣。經(jīng)過大量抽樣，再去除前面的退化期，從每次抽樣中收集參數(shù)樣本，就可利用抽樣結(jié)果進(jìn)行后驗(yàn)估計(jì)與推斷。Hastings證明了MCMC抽樣過程在排除最初的若干次抽樣后能達(dá)到穩(wěn)定狀態(tài)，一旦抽樣池達(dá)到穩(wěn)定狀態(tài)，將抽樣分布作為后驗(yàn)分布解釋就是合理的[8]。

貝葉斯空間高斯類模型的相關(guān)理論研究文獻(xiàn)一經(jīng)發(fā)表，就因其突出的優(yōu)點(diǎn)、比較廣泛的適應(yīng)性吸引了學(xué)術(shù)界的注意，越來越多的學(xué)者將其應(yīng)用到學(xué)術(shù)研究中，應(yīng)用領(lǐng)域不斷擴(kuò)展，研究內(nèi)容不斷豐富，也為其他類型的貝葉斯空間計(jì)量模型研究奠定了基礎(chǔ)。如Seya等為了研究日本從1989—2007年間的收入收斂性，構(gòu)建了一個(gè)空間SDM模型，探討了其理論背景，做出了收入的空間分布圖；模型的估計(jì)結(jié)果給出了日本在這一時(shí)期各個(gè)階段不同的收入差異和幾種不同的收斂性[9]。Yu等借鑒了過去研究能源消耗的時(shí)間序列和截面數(shù)據(jù)方法，考慮了能源消耗的時(shí)間-空間關(guān)系及面板空間數(shù)據(jù)的異質(zhì)性，針對(duì)中國30個(gè)省市從2007—2009年能源消耗數(shù)據(jù)，用貝葉斯MCMC抽樣方法得到后驗(yàn)仿真樣本，并用之計(jì)算樣本邊緣分布和后驗(yàn)分布的矩，研究了中國的能源消耗[10]。

(二)非高斯模型

y*～TMVN(μ,Ω)

(21)

p(β|ρ,y*)∝N(c*,T*)

c*=(X′X+T-1)-1(X′Sy*+T-1c)

T*=(X′X+T)-1

S=In-ρW

p(ρ|β,y*)∝|In-ρW|exp

(22)

在SAR模型中，y*的條件分布服從截尾正態(tài)分布。對(duì)于非截尾的n維多元正態(tài)分布，可以依次從n個(gè)條件一元正態(tài)分布中抽樣以獲得n個(gè)參數(shù)；但Geweke指出該方法不適用于截尾多元正態(tài)分布，否則處理SAR模型時(shí)會(huì)出現(xiàn)錯(cuò)誤[12]。

LeSage等使用了截尾多元正態(tài)分布的精度矩陣[3]，將該矩陣標(biāo)記為Ψ=Ω-1。由于z為多元正態(tài)分布，其中各個(gè)元素Zi基于其他所有元素Z-i的條件概率分布可以表示為具有條件均值和條件方差的一元分布，且計(jì)算簡便。在一定約束下，這些均值和方差的表達(dá)式可用于從一元截尾分布中進(jìn)行抽樣，由此就可用Gibbs抽樣從多元截尾正態(tài)分布中構(gòu)建一個(gè)樣本。

對(duì)于非截尾多元正態(tài)分布，Geweke使用分塊對(duì)稱矩陣的逆構(gòu)造E(zi|z-i)=γ-iz-i，其中γ-i=-Ψ-i/Ψi,i，Ψ-i是Ψ中第i行除該元素外的其他所有元素。而對(duì)于截尾正態(tài)分布，使用的條件分布應(yīng)為：

zi|z-i=γ-iz-i+hivi；hi=(Ψi,i)-1/2

(23)

vi～N(0,1)的抽樣服從以下截尾約束：

由此可以產(chǎn)生一個(gè)由zi(i=1,2,…,n)構(gòu)成的向量z，在對(duì)元素zi進(jìn)行抽樣時(shí)，可使用前期抽樣的結(jié)果z1,z2,…,zi-1。此外，在更新zi時(shí)，還會(huì)使用從前一次Gibbs抽樣中得到的zi+1,zi+2,…,zn。通過對(duì)含有n個(gè)觀測值的向量z進(jìn)行m次抽樣，可以得到y(tǒng)*=v+z(m)，其中y*則可以用于從模型參數(shù)β和ρ的條件分布中進(jìn)行抽樣。

在對(duì)截尾多元正態(tài)分布的抽樣中，設(shè)定MCMC抽樣器的迭代步驟m=10，每次MCMC抽樣僅需y*的樣本值，由于產(chǎn)生y*值時(shí)進(jìn)行了數(shù)千次抽樣，所得y*相對(duì)精確，但Gibbs抽樣中，僅需一步即可，即m=1，此時(shí)用由前一次MCMC取樣器產(chǎn)生的z值。m值的減少將顯著減少模型估計(jì)的時(shí)間，可以提高計(jì)算機(jī)運(yùn)行速度。同時(shí)，因?yàn)樵摲椒軌驅(qū)⒋髽颖緮?shù)據(jù)劃分為從一元條件分布中抽取的一系列樣本，可以不需要借助大容量計(jì)算機(jī)內(nèi)存就可以處理大樣本數(shù)據(jù)。上述過程還可以采用一些其他有效的計(jì)算方法，比如可通過使用精度矩陣Ψ，而不需要計(jì)算逆矩陣[(In-ρW)′(In-ρW)]-1，從而降低了對(duì)電腦內(nèi)存的要求。

將上述的貝葉斯MCMC方法進(jìn)行拓展，可以處理有序空間概率單位模型、空間Tobit模型和多項(xiàng)空間概率單位模型等多種非高斯模型，限于篇幅不一一贅述，具體可參見Qiang Zeng等的研究[13]。

非高斯貝葉斯空間計(jì)量模型的具體應(yīng)用可參見：Andrea Amaral等為了研究亞洲20世紀(jì)90年代的銀行危機(jī)，構(gòu)建了一個(gè)空間概率模型用于測度不同國家銀行系統(tǒng)間的傳染效應(yīng)，并在空間概率模型中考察了潛變量及其滯后項(xiàng)的空間效應(yīng)，同時(shí)也考察了誤差項(xiàng)的空間效應(yīng)，用Gibbs抽樣方法估計(jì)模型并分析了這種傳染效應(yīng)的產(chǎn)生機(jī)制[14]；Changjoo Kim等人建立了一個(gè)貝葉斯多元Tobit模型，模擬了家庭成員的交互行為是如何影響其他個(gè)體旅行行為的[15]。在他們所采用的數(shù)據(jù)集中，涉及到個(gè)人旅行數(shù)據(jù)存在大量刪失值，此時(shí)多元Tobit模型則可以使用截尾多元t條件后驗(yàn)分布對(duì)潛在的無法觀測的效用進(jìn)行分析。這些文獻(xiàn)涉及經(jīng)濟(jì)預(yù)測、經(jīng)濟(jì)結(jié)構(gòu)分析和經(jīng)濟(jì)行為模擬分析等多個(gè)領(lǐng)域，為實(shí)證分析研究者們提供了極好的借鑒。

(三)廣義線性模型

空間廣義線性模型中被解釋變量與其滯后項(xiàng)和解釋變量之間不再是呈線性關(guān)系，高斯模型和非高斯模型都可被視為廣義線性模型的特殊情形，因此這兩類模型都可用常規(guī)方法加以處理。對(duì)于其他的一些空間廣義線性模型，如Griffith研究了一個(gè)引入空間依賴的自邏輯斯蒂(autologistic)回歸模型[16]：

p(Yi=1|α,WiY)

(24)

其中α為一個(gè)能夠捕獲大尺度變化的參數(shù)，可依據(jù)Xi來進(jìn)行設(shè)置；ρ為空間自相關(guān)參數(shù)；Wi為相對(duì)于位置i的由wij構(gòu)成的行向量；行向量Wi構(gòu)成的空間權(quán)重矩陣W為對(duì)稱的。

該模型的空間自相關(guān)性能夠被Cliff等提出的聯(lián)合計(jì)數(shù)(join-count)統(tǒng)計(jì)量所測定[17]179-183。兩個(gè)1在地理位置相鄰記為BB，一個(gè)1和一個(gè)0地理相鄰記為BW，若兩個(gè)0地理相鄰記WW，把這些聯(lián)合計(jì)數(shù)統(tǒng)計(jì)量標(biāo)準(zhǔn)化就能得到其相應(yīng)的Z值。進(jìn)一步地，如果在模型(24)中考慮n×(K+1)階協(xié)變量矩陣X，矩陣X造成了地理二項(xiàng)數(shù)據(jù)的空間自相關(guān)，其均值響應(yīng)的設(shè)置將會(huì)使得模型的空間依賴參數(shù)值為0；由于變量缺失于均值響應(yīng)設(shè)置，因此空間自相關(guān)出現(xiàn)在殘差項(xiàng)中，可引入特征向量濾子加以處理。特征向量濾子是一種非參數(shù)方法，這種方法利用了空間自相關(guān)性的設(shè)定偏誤，即假設(shè)空間自相關(guān)性是由具有空間自相關(guān)的缺失外生變量帶來的；空間自相關(guān)潛在于協(xié)變量矩陣X，矩陣X造成了因變量Y的空間自相關(guān)性，其基本思想是把廣義線性回歸模型中固有的空間自相關(guān)性視為缺失變量效應(yīng)而將其剔除，從而利用一組空間代理變量來控制空間自相關(guān)性。代理變量為二項(xiàng)矩陣W的修正形式矩陣：

(25)

的特征向量，其中I為單位矩陣，1為n×1階1向量，n為區(qū)域個(gè)數(shù)。

利用空間權(quán)值矩陣可以計(jì)算出空間濾子的每一個(gè)特征向量E1，E2，…，En，這些特征向量代表了不同的空間模式，也相關(guān)于給定的空間自相關(guān)水平。因此，模型(24)可以被下述模型替代：

(26)

其中K表示n個(gè)特征向量的一個(gè)子集中所含特征向量的個(gè)數(shù)，即Ei,K為一個(gè)n×K階矩陣，其列是K個(gè)被選擇的特征向量。這樣，就把模型中的空間依賴效應(yīng)轉(zhuǎn)換為由Ei,Kβ所表示的大尺度變異項(xiàng)，使得參數(shù)ρ為0?？梢钥吹?，模型(24)和模型(26)間存在下述關(guān)系：

p(Y=1|WY)=p(Y=1|EΛETY)

=p(Y=1|Eδ)

(27)

其中Λ為一個(gè)特征值對(duì)角矩陣，δ為系數(shù)向量。這樣的處理方法關(guān)注了特征向量空間濾子而避免了處理協(xié)變量，大大簡化了所要研究的問題。

特征向量空間濾子方法把一個(gè)最小充分的特征向量集合作為代理變量施加于線性預(yù)測子集合，因此如果考慮到獨(dú)立性，就可以排除觀測值間的空間自相關(guān)性。這種方法比MCMC運(yùn)算更容易操作，參數(shù)估計(jì)比起擬極大似然估計(jì)更為有效，而且由于其對(duì)地理二項(xiàng)數(shù)據(jù)中隱含的自相關(guān)模式進(jìn)行了分解，這使得模型的空間自相關(guān)效應(yīng)能夠得到更準(zhǔn)確的解釋。

四、貝葉斯空間計(jì)量模型選擇理論及其應(yīng)用

空間計(jì)量模型選擇主要分為：空間和非空間模型的比較、不同空間權(quán)重結(jié)構(gòu)模型的比較和納入不同候選解釋變量模型的比較?；谪惾~斯統(tǒng)計(jì)理論，Zellner建立了進(jìn)行多個(gè)模型比較的方法[18]。該方法設(shè)定每個(gè)模型的先驗(yàn)概率、回歸參數(shù)的先驗(yàn)分布，計(jì)算后驗(yàn)?zāi)Ｐ透怕?，并用于推斷不同回歸模型與先驗(yàn)信息的一致性，其基本思想考慮m個(gè)空間回歸模型：M=M1,M2,…,Mm,對(duì)于這m個(gè)模型給予一個(gè)先驗(yàn)概率π(Mi)，根據(jù)樣本的條件分布與參數(shù)的先驗(yàn)分布，可以得到第i個(gè)模型的后驗(yàn)分布，據(jù)此則可計(jì)算模型的后驗(yàn)概率。

對(duì)于一個(gè)空間計(jì)量模型，可以把模型和參數(shù)的聯(lián)合后驗(yàn)分布表示為：

(28)

其中π(M)為模型的先驗(yàn)概率，π(η|M)為模型參數(shù)在給定模型條件下的條件概率，π(D|η,M)為樣本數(shù)據(jù)的條件分布。視η為冗余參數(shù)，則模型的邊緣條件分布可以表示為：

(29)

特別地，對(duì)于SAR模型，模型參數(shù)的后驗(yàn)分布可表示為：

(30)

e=(In-ρW)y-αιn-Xβ

模型參數(shù)中β服從條件正態(tài)先驗(yàn)N(β0，σ2(gX′X)-1)，σ2為倒伽瑪分布IG(a,b)，截距參數(shù)α為無信息先驗(yàn)，ρ為Be(d,d)先驗(yàn)。因此，式(28)可寫為：

其中

S(ρ)=e(ρ)′e(ρ)

C=gX′X

(31)

若僅對(duì)p(η|D)中的β和σ2積分，可消掉這兩個(gè)參數(shù)，可以得到：

(32)

其中

A=In-ρW

(33)

因此，在考慮模型比較的情況下，可采用貝葉斯理論的后驗(yàn)概率比的思想，對(duì)于有m個(gè)不同加權(quán)矩陣W(i),i=1,2,…,m的情形，根據(jù)似然函數(shù)和參數(shù)的先驗(yàn)分布，每個(gè)模型視為是不同的模型，則：

(34)

因此給定不同的模型(空間權(quán)重矩陣)，就可得到后驗(yàn)?zāi)Ｐ透怕剩?/p>

(35)

依賴不同模型的后驗(yàn)概率比為：

(36)

該方法的優(yōu)點(diǎn)在于比較模型的類別廣泛，可用于不同的空間權(quán)重矩陣的比較，而且不同模型設(shè)定比較或不同的解釋變量比較，皆可用后驗(yàn)概率比或貝葉斯因子來比較非嵌套模型。 但是，采用上述方法時(shí)需要注意的問題是，要求避免Lindley提出的一個(gè)悖論：當(dāng)比較不同數(shù)量參數(shù)的模型時(shí)，這些參數(shù)依賴擴(kuò)散先驗(yàn)分布，無論樣本數(shù)據(jù)信息如何，簡單的模型比復(fù)雜模型更好，即具有相同數(shù)量參數(shù)的兩個(gè)模型可以用擴(kuò)散先驗(yàn)來比較，但對(duì)包含參數(shù)個(gè)數(shù)變化的模型比較必須開發(fā)和采用先驗(yàn)策略。

對(duì)于包含不同解釋變量的替代線性回歸模型進(jìn)行貝葉斯模型平均，Madigan和York介紹了一種馬爾科夫鏈蒙特卡羅模型組成的(MC3)方法用于空間回歸模型[19]。對(duì)于有k個(gè)解釋變量的回歸模型，有2k種可能方式來選擇或排除回歸元。因此，模型平均需要納入足夠多的變量以克服遺漏變量偏誤和因包括冗余變量而降低模型的精度。

該方法的原理是建立一個(gè)策略隨機(jī)馬爾科夫鏈過程，用M代表模型目前的狀態(tài)，對(duì)于每個(gè)M∈M，定義一個(gè)鄰域nbd(M)，其中包括了模型M和其他模型的集合族，該集合族中的集合較之M具有更多或更少的邊緣。同時(shí)，定義一個(gè)轉(zhuǎn)移概率矩陣q，若對(duì)于一切M′?nbd(M)，則q(M→M′)=0；對(duì)一切M′∈nbd(M)，則q(M→M′)=常數(shù)。按照接受概率：

(37)

持續(xù)進(jìn)行抽樣，并得到建議模型M′。在貝葉斯方法中，將解釋變量設(shè)定的隨機(jī)性融入估計(jì)和推斷的解決方法是平均不同模型的設(shè)定，這與很多依賴單一模型設(shè)定的方法形成對(duì)比，采用不同的模型比較準(zhǔn)則得到最優(yōu)模型。程序確定了與特定解釋變量相關(guān)的模型，對(duì)每一個(gè)模型指定了一個(gè)后驗(yàn)?zāi)Ｐ透怕剩渲笜?biāo)可以直接計(jì)算，易于操作，結(jié)果更具有穩(wěn)健性，不受樣本量大小的限制。

有關(guān)貝葉斯空間計(jì)量模型選擇的其他論述還可參見陶長琪等的研究，他們比較了傳統(tǒng)的方法與MCMC方法，仿真研究顯示MCMC方法的檢驗(yàn)效度更高[20]。在模型選擇的具體應(yīng)用方面，曾召友等基于貝葉斯空間面板計(jì)量模型理論，研究了中國電信服務(wù)的外溢性，對(duì)空間面板計(jì)量模型簇中的嵌套和非嵌套情形進(jìn)行了整合處理，使得研究結(jié)果具有穩(wěn)健性[21]。

五、結(jié) 語

空間計(jì)量學(xué)從20世紀(jì)80年代誕生以來，受到了計(jì)量經(jīng)濟(jì)學(xué)家的廣泛關(guān)注，經(jīng)歷了從萌芽到成熟、從邊緣到主流的一個(gè)飛躍；眾多統(tǒng)計(jì)學(xué)、經(jīng)濟(jì)學(xué)、數(shù)學(xué)和計(jì)算科學(xué)研究學(xué)者為其方法體系的建立、發(fā)展和完善做出了大量的努力，其中貝葉斯方法的引入為空間計(jì)量模型的研究注入了新的活力。傳統(tǒng)的空間計(jì)量模型設(shè)定參數(shù)是確定的，在這種設(shè)定下建立的模型往往不符合現(xiàn)實(shí)經(jīng)濟(jì)結(jié)構(gòu)和行為方式，從而導(dǎo)致模型估計(jì)結(jié)果偏差很大，預(yù)測結(jié)果極為不準(zhǔn)確。如今，貝葉斯統(tǒng)計(jì)方法與計(jì)算科學(xué)的飛速發(fā)展可以很好地解決這一問題，這不得不讓經(jīng)濟(jì)學(xué)家們開始審視空間模型參數(shù)隨機(jī)性的問題。從這個(gè)意義上講，貝葉斯空間計(jì)量模型的發(fā)展與應(yīng)用無疑會(huì)對(duì)未來經(jīng)濟(jì)科學(xué)的研究范式革新與理論邊界擴(kuò)展產(chǎn)生影響。從現(xiàn)有文獻(xiàn)來看，貝葉斯空間計(jì)量模型對(duì)于一些常見模型已獲得一定有意義的成果，可將其運(yùn)用于經(jīng)濟(jì)研究各個(gè)領(lǐng)域，并且用其分析的結(jié)果較為準(zhǔn)確；從模型的靈活性來看，貝葉斯空間計(jì)量模型也優(yōu)于非貝葉斯模型。然而在實(shí)際問題運(yùn)用中，貝葉斯統(tǒng)計(jì)分析方法在模型推斷和模型選擇中的應(yīng)用還是有限制的，這種限制是由于參數(shù)數(shù)量增加帶來的計(jì)算困難，而不是傳統(tǒng)方法的自由度不足。

因此，貝葉斯空間計(jì)量模型的理論與應(yīng)用研究也面臨著一些巨大的挑戰(zhàn)。首先，空間計(jì)量模型大數(shù)據(jù)問題：目前的空間計(jì)量模型所采用的數(shù)據(jù)基本是人工記錄數(shù)據(jù)，但隨著科學(xué)技術(shù)的進(jìn)步，未來的空間計(jì)量模型會(huì)利用到自動(dòng)記錄數(shù)據(jù)。在自然或社會(huì)科學(xué)中的傳感器系統(tǒng)能夠從時(shí)間維度、空間維度和個(gè)體維度提供更大量的數(shù)據(jù)，目前的貝葉斯空間計(jì)量模型的理論方法則不適用于這種大數(shù)據(jù)情形。其次，隨著空間計(jì)量大數(shù)據(jù)時(shí)代的到來，與之相應(yīng)的就是對(duì)計(jì)算技術(shù)要求的提高。要能夠處理大數(shù)據(jù)中復(fù)雜的時(shí)空交互反應(yīng)，就必須要求有更高效、更快速、更準(zhǔn)確的新的計(jì)算技術(shù)，例如分布式計(jì)算技術(shù)、云計(jì)算和手持裝備的運(yùn)用。

中國學(xué)者對(duì)貝葉斯統(tǒng)計(jì)理論與方法的研究相對(duì)較晚，與西方國家相比尚有很大差距，這主要是因?yàn)樨惾~斯方法較之經(jīng)典統(tǒng)計(jì)方法在中國的應(yīng)用較少，該方法還不為大家所熟悉；并且其計(jì)算復(fù)雜，對(duì)于高維數(shù)據(jù)計(jì)算尤為困難，缺乏相應(yīng)的軟件支撐，這妨礙了貝葉斯方法在實(shí)際問題中的使用。目前而言，貝葉斯空間計(jì)量模型既有廣闊的發(fā)展空間，又面臨前沿問題的挑戰(zhàn)，這正是促進(jìn)其不斷完善、不斷進(jìn)步的動(dòng)力。

[1] LeSage J P.Bayesian Estimation of Spatial Autoregressive Models[J].International Regional Science Review，1997，1&2(20).

[2] LeSage J P，Pace R K.Introduction to Spatial Econometrics[M].New York：CRC Press，2009.

[3] LeSage J P，Kelley Pace，R，Lam N，et al.New Orleans Business Recovery in the Aftermath of Hurricane Katrina[J].Journal of the Royal Statistical Society.Series A，2011,174(4).

[4] Nelder J，Wedderburn R.Generalized Linear Models[J].Journal of the Royal Statistical Society.Series A，1972，135(3).

[5] Hepple L W.Bayesian Techniques in Spatial and Network Econometrics 1，Model Comparison and Posterior Odds[J].Environment and Planning A,1995，27(2).

[6] LeSage J P，Parent O.Bayesian Model Averaging for Spatial Econometrics Models[J].Geographical Analysis,2007,3(39).

[7] Pace R K，Barry Ronald.Quick Computation of Spatial Autoregressive Estimators[J].Geographic Annal，1997，29(3).

[8] Hastings W K.Monte Carlo Sampling Methods Using Markov Chains and Their Applications[J].Biometrika,1970，57(7).

[9] Seya Hajime，Tsutsumi Morito，Yamagata Yoshiki.Income Convergence in Japan：A Bayesian Spatial Durbin Model Approach[J].Economic Modelling,2012，29(1).

[10] Yu Yanpin，Li Xi， Qian Zhengming.A Study on China's Energy Consumption by Bayesian Spatial Econometric Model[J].Energy Procedia,2012(16).

[11] Albert J H，Chib S.Bayesian Analysis of Binary and Polychotomous Response Data[J].Journal of The American Statistical Association,1993,88(22).

[12] Geweke J.Efficient Simulation from the Multivariate Normal and Student-t Distributions Subject to Linear Constraints and the Evaluation of Constraint Probabilities[C].In Proceeding of 23rd Symposium on the Interface Between Computing Science and Statistics，Interface Foundation of North America Inc.Fairfax,1991.

[13] Qiang Zeng，Wei Huiying，Huang Helai.A Bayesian Spatial Random Parameters Tobit Model for Analyzing Crash Rates on Roadway Segments[J].Accident analysis and Prevention,2017，100(3).

[14] Andrea Amaral，Margarida Abreu，Victor Mendes.The Spatial Probit Model—An Application to the Study of Banking Crises at the End of the 1990's[J].Physica A：Statistics Mechanics and Its Applications,2014，415(415).

[15] Changjoo Kim，Olivier Parent.Modeling Individual Travel Behaviors Based on Intra-household Interactions[J].Regional Science and Urban Economics,2016，57(5).

[16] Griffith D A.A Spatial Filtering Specification for the Auto-logistic Model[J].Environment and Planning,2004，36(10).

[17] Cliff A D，Ord J.Spatial Processes[M].London：Pion,1981.

[18] Zellner A.Bayesian Econometrics[J].Econometrica,1985,53(2).

[19] Madigan D，York J.Bayesian Geographical Models for Discrete Data[J].International Statistical Review,1995，63(2).

[20] 陶長琪，楊海文.空間計(jì)量模型選擇及其模擬分析[J].統(tǒng)計(jì)研究,2014,31(8).

[21] 曾召友，龍志和，董大勇.基于Bayes理論的空間面板計(jì)量模型遴選框架——以中國電信服務(wù)外溢性分析為例[J].統(tǒng)計(jì)與信息論壇，2010(3).