蔣 捷，鄭月晨，周 浩，張慧增

(杭州師范大學理學院，浙江 杭州 311121)

門限自回歸(threshold autoregressive, TAR)模型是典型的非線性參數模型,可分為沖量門限自回歸(momentum-TAR, MTAR)模型和自激勵門限自回歸(self-exciting TAR, SETAR)模型[1].在經濟學研究中，很多與經濟相關的時間序列大部分都具有非線性特征,如股價指數和匯率等.在處理非線性問題時,建立傳統(tǒng)的線性時間序列模型不能充分地提取數據信息,從而也顯示出了線性模型的局限性.SETAR模型是基于不同的門限變量把整個時間序列分為若干個小段(體制),在每個體制上用AR(autoregressive)模型去線性逼近,把全局空間上的非線性時間序列模型轉化成子空間上的線性模型.雖然二體制的TAR模型形式簡單,最小二乘法也可以估計出模型參數,但繁瑣的過程和估計精度很難應用到具體問題中.對于TAR模型,蒙特卡洛模擬(Markov Chain Monte Carlo,MCMC)方法不僅可以相對快速地估計出模型參數,還可以提高估計的精度,實施起來也比較便捷.本文運用了Gibbs抽樣,在機器學習中,很多情況無法確定一個概率分布的具體密度函數,只是知道樣本其中一個屬性在其他所有屬性下的條件概率,而Gibbs抽樣正是解決此類問題的方法之一.TAR模型的參數估計方法主要分為以下幾種：一是最小二乘法結合赤池信息量(akaike information criterion, AIC),但此方法計算量大且繁瑣,很難推廣并且實現;二是Tasy提出的F統(tǒng)計量以及t統(tǒng)計量來選擇各參數的值;三是利用MCMC方法.本文利用MCMC算法與Gibbs抽樣,可避免繁瑣的迭代以及復雜的似然比檢驗,對模型的選擇和參數的估計實施更為快捷有效.

1 門限自回歸模型及其參數的貝葉斯估計

1.1 門限自回歸模型

若y={yt,t=1,2,…}為一時間序列,滿足以下結構的時間序列模型稱為門限自回歸模型,記為TAR(k,p1,p2,…,pk):

(1)

下面以二體制(即k=2)的門限自回歸模型TAR(2,p1,p2)來說明TAR模型的參數估計，即

(2)

設?p∈N*,使得0≤p1,p2,k=2,d≤p,{y1,y2,…,yp}為時間序列{yt}的前p個觀測值.把{yp+1,yp+2,…,yn}向前移動d個單位,令πi為移動后新時間序列的第i個觀測值，即yπ1

f(Y|Φi,σ2,r,d;i=1,2)=

(3)

為了解決當數據維數增加時計算和書寫變復雜的問題,把式(3)寫成矩陣相乘的形式.令

則式(3)可寫為

(4)

1.2 基于貝葉斯分析下TAR模型的參數估計

貝葉斯學派與頻率學派的主要差異在于:前者利用收集到的先驗信息,形成先驗分布,在給定樣本分布和先驗分布之后,運用貝葉斯公式計算待估參數的后驗分布.在求得后驗分布后,模擬該分布生成的隨機數[4],若把隨機數看成樣本，則用樣本均值和樣本方差來估計與評價待估參數.

設(Φi,σ2,r,d)的先驗分布為π(Φi,σ2,r,d),由條件概率公式可知

π(Φi,σ2,r,d)=π(Φi|σ2,r,d)·π(σ2|r,d)·π(r|d)·π(d).

(5)

由式(4)和(5)可求得(Φi,σ2,r,d)的后驗分布

(6)

設θ是總體分布的參數向量,π(θ)是參數向量θ的先驗分布,若后驗密度函數與π(θ)有相同的密度函數形式,則π(θ)稱之為參數向量θ的共軛先驗分布.根據共軛先驗分布的特點,因為Φ1,Φ2相互獨立,可以得到Φi的條件后驗分布

(7)

σ2的條件后驗分布：

(8)

r的條件后驗分布：

p(r|Y,Φi,σ2,d)∝exp(-s2/2σ2)I(a

(9)

其中I表示示性函數.

后驗分布包含了對參數的信息，如均值、方差，但很多情況無法確定一個概率分布的具體密度函數，或密度函數太復雜，只是知道樣本其中一個屬性在其他所有屬性下的條件概率，此時上述后驗分布正好滿足這樣的要求，所以對每個參數的估計可以通過Gibbs抽樣來實現.

2 Gibbs抽樣

Gibbs抽樣是一種特殊的MCMC算法[5]，主要用于狀態(tài)空間為乘積空間的情況.往往在待估參數維數較高的情況下，貝葉斯統(tǒng)計推斷更傾向于求出總體對某參數的條件分布，Gibbs抽樣方法適用于條件分布比邊緣分布更容易求得的情形.該方法每次只更新狀態(tài)的一個分量，這樣就能把高維抽樣轉化為低維抽樣的問題.對于TAR模型的Gibbs抽樣步驟如下:

3 數值模擬

作為上述方法的檢驗,超參數的選取為:Φ1～N(0,0.2),Φ2～N(0,0.2),σ2～IG(1.2,0.5),r～U[p5,p95][6],d服從1,2,3上的離散均勻分布;其中p5,p95分別表示數據的第5和第95個百分位數,εt～N(0,σ2),以TAR(2,1,1)模型為例:

yt=[1+0.5yt-1I(yt-d<0.6)]+[1-0.9yt-1I(yt-d>0.6)]+εt

(10)

生成500個模型(10)的數據,數據時序圖如圖1所示.

圖1 生成樣本的時序圖Fig.1 Sequence diagrams of sample generated

表1 被估參數模擬結果Tab.1 Simulation results of parameters

利用上述生成的數據,用Gibbs抽樣方案進行迭代,被估參數的后驗均值和后驗方差如表1所示.表1記錄了模型(10)各個參數的真實值、后驗均值和標準誤.

4 結論

根據貝葉斯統(tǒng)計推斷,通過先驗分布(共軛先驗分布),對TAR模型中各參數的后驗分布進行分析推導,同時利用MCMC算法與Gibbs抽樣解決了參數從高維到低維的合理性轉變,提高了估計精度.最后的數值模擬得到的結果也與真實值偏差較小,說明該方法是有效的.