邵思雨, 龍 琪, 張 潔
(長春工業(yè)大學 數(shù)學與統(tǒng)計學院, 吉林 長春 130012)
時間序列數(shù)據(jù)是指把具有某一個相同統(tǒng)計指標的數(shù)值依據(jù)它發(fā)生的日期先后順序排列而成的數(shù)列。而取值為非負整數(shù)的時間序列數(shù)據(jù)作為一種計數(shù)數(shù)據(jù),其研究范圍更為廣泛。因為這類時間序列數(shù)據(jù)只能取整值,而傳統(tǒng)實數(shù)范圍上的時間序列數(shù)據(jù)能取到小數(shù),所以這兩種數(shù)據(jù)的研究存在較大差別。一般情況下,依據(jù)非負整數(shù)值時間序列數(shù)據(jù)有無上界,可以將此數(shù)據(jù)劃分成有上限的數(shù)據(jù)(0,1,…,n)和無上限的數(shù)據(jù)(0,1,…)。例如:某一地區(qū)一周七天下雨的天數(shù)數(shù)據(jù),以及某一地區(qū)患新型肺炎的人數(shù)等。這些年,研究者為了擬合這類數(shù)據(jù)建立了不少模型。迄今為止,最常見的兩種建模方法分別為基于潛過程構(gòu)造的狀態(tài)空間模型[1],以及基于稀疏算子提出的整數(shù)值自回歸模型[2-6]。
針對取值范圍有上限的整數(shù)值時間序列數(shù)據(jù),例如長春市每天有多少個城區(qū)發(fā)生交通事故數(shù)據(jù)、全國每個月有多少城市發(fā)生地震數(shù)據(jù)等。擬合這類時間序列數(shù)據(jù)較為常見的方法是McKenzie E[7]提出的First-order binomial autoregre-ssive(BAR(1))模型。Weiβ C H[8]研究BAR(1)模型的相應性質(zhì),并將其推廣到高階情形。此后,越來越多BAR模型的衍生模型涌現(xiàn)[9-11]。然而,二項稀疏算子是基于獨立同分布的伯努利隨機變量序列提出的,所以,如果將BAR(1)模型應用到個體之間具有相依結(jié)構(gòu)的數(shù)據(jù)中,可能會使該模型的擬合效果不好。基于此,Ristic M M等[12]提出一個新的算子——廣義二項稀疏算子,它的提出為具有相依結(jié)構(gòu)的整數(shù)值時間序列數(shù)據(jù)的研究提供了參考價值。近期,Kang Y等[13]和康堯[14]將廣義二項稀疏算子應用到二項自回歸模型中,得出一類一階廣義二項自回歸(GBAR(1))模型,并且對該模型進行了相應的統(tǒng)計推斷。由于此算子放松了對二項稀疏算子的獨立性假設,從而對于包含相依結(jié)構(gòu)的整數(shù)值時間序列數(shù)據(jù)的研究具有較大意義,并且能夠有效地解決古典BAR(1)模型解釋實際背景問題時的部分局限性。
對于上述模型的統(tǒng)計推斷研究,模型相應的參數(shù)估計問題是研究的重點內(nèi)容之一。除了常用的Yule-Walker 估計、條件最小二乘估計、條件最大似然估計和經(jīng)驗似然估計以外,貝葉斯估計也是一種對模型參數(shù)進行統(tǒng)計推斷的方法,不少學者對貝葉斯方法進行了研究。Lazar N A等[15]針對貝葉斯經(jīng)驗似然檢驗了后驗推斷的有效性。張哲等[16]利用貝葉斯方法研究了First-order integer-valued autoregressive(INAR(1))模型參數(shù)的估計問題。Xi R等[17]針對spike-and-slab先驗的經(jīng)驗似然的貝葉斯分位數(shù)回歸,研究了其中非參數(shù)貝葉斯變量選擇等問題。朱復康等[18]利用矩方法和貝葉斯方法研究Integer-valued generalized autoregressive conditional heteroscedastic(INGARCH (1,1))模型的參數(shù)估計問題。近期部分研究詳見文獻[19-20] 。
基于此,考慮針對GBAR(1)模型,研究相應的貝葉斯統(tǒng)計推斷問題,并將該模型應用于實際數(shù)據(jù)的分析研究中。通過貝葉斯估計方法得到有效的擬合結(jié)果,為具有相依結(jié)構(gòu)的實際問題研究提供理論支撐。
Kang Y等[13]將廣義二項稀疏算子 “ °θ”引入到BAR (1)模型中,提出一類一階廣義二項自回歸(GBAR (1))模型,并且給出模型的Yule-Walker估計、條件最小二乘估計、條件最大似然估計。文中通過進一步研究給出了GBAR (1)模型的貝葉斯估計。首先給出GBAR (1)模型的定義及相關性質(zhì)。
定義1如果{Xt}是有上限的非負整數(shù)值時間序列數(shù)據(jù),并且{Xt}t∈N滿足方程:
Xt=α°θΧt-1+β°θ(n-Χt-1),t≥1,
(1)
其中n∈N為數(shù)據(jù)上限,α、θ∈(0,1),各稀疏運算之間相互獨立,則稱{Xt}為廣義二項自回歸模型,記作GBAR(1)模型。
“α°θ”是廣義二項稀疏算子, 其定義為:
Ui=(1-Vi)Wi+ViZ,
式中:{Wi}i∈N,{Vi}i∈N----兩列獨立同分布的Bernoulli(α)和Bernoulli(θ)隨機變量序列;
Z----一個Bernoulli(α)隨機變量。
{Wi}i∈N,{Vi}i∈N和Z之間是相互獨立的,α、θ∈(0,1)。
對于給定的n∈N,GBAR(1)模型的轉(zhuǎn)移概率為:
P(Xt=k|Xt-1=l)=
(2)
在轉(zhuǎn)移概率已知的條件下進一步考慮GBAR(1)模型的貝葉斯估計。假設X1,X2,…,XT是來自GBAR(1)模型的T個樣本觀測值,α,β,θ為待估參數(shù)。由式(2)可得樣本的聯(lián)合似然函數(shù)為
P(x|α,β,θ)=L(x|α,β,θ)=
(3)
取參數(shù)α,β,θ的先驗分布分別為Beta(e,f)、Beta(p,q)和Beta(h,g),這里我們分別記作π(α),π(β),π(θ),即
(4)
(5)
(6)
由此得出聯(lián)合密度函數(shù)
P(x,α,β,θ)=L(x|α,β,θ)×π(α)×
π(β)×π(θ)。
(7)
在給定聯(lián)合似然函數(shù)以及先驗分布后,由貝葉斯原理可知,參數(shù)α,β,θ的聯(lián)合后驗為
(8)
(9)
(10)
(11)
對GBAR(1)模型進行數(shù)值模擬,通過參數(shù)估計的偏差和均方誤差來比較Conditional Least Squares (CLS) 估計、Conditional Maximum Likelihood (CML) 估計和Bayesian (BYE)估計的估計效果。根據(jù)參數(shù)α,β,θ的取值空間,選取以下四組不同的參數(shù)組合進行數(shù)值模擬:
1)n=5,α=0.25,β=0.75,θ=0.65;
2)n=5,α=0.40,β=0.60,θ=0.55;
3)n=8,α=0.40,β=0.20,θ=0.45;
4)n=8,α=0.80,β=0.20,θ=0.35。
選取參數(shù)α,β,θ的先驗分布為貝塔分布。令貝塔分布的均值等于參數(shù)真值的方法來選取超參數(shù),并令幾組參數(shù)組合的先驗分布分別為:
1)α~Beta(1,3),β~Beta(6,2),θ~Beta(2,1);
2)α~Beta(2,3),β~Beta(6,4),θ~Beta(2,2);
3)α~Beta(2,3),β~Beta(1,4),θ~Beta(1,1);
4)α~Beta(4,1),β~Beta(1,4),θ~Beta(1,2)。
鑒于在進行貝葉斯估計時,需要計算多重積分,文中采用樣本平均值法求解數(shù)值積分。此方法是利用樣本的均值來估計所求的定積分
若Y服從(c,d)上的均勻分布,則有
得出定積分
I=E(g(Y))×(d-c)。
因此,當序列{Yt,t=1,2,…,N}獨立同分布于(c,d)上均勻分布,則
Gi=g(Yi),i=1,2,…,N
同為一列獨立同分布的隨機變量序列,根據(jù)強大數(shù)定律,則有
所以
是I的強相合無偏估計,如此求解定積分的方法稱為樣本平均值法。
運用此法求積分
能將積分求解問題轉(zhuǎn)化為求g(Y)期望的問題。在數(shù)值模擬時,樣本平均值法對應的樣本量N取值為105。
模擬考慮在觀測數(shù)據(jù)樣本量T分別為20,50,100 和200的情況下進行,最終結(jié)果是在R軟件環(huán)境下取 500 次估計結(jié)果的平均值。四個參數(shù)組合的樣本路徑和ACF如圖 1所示。
從圖1可以看出, 序列為平穩(wěn)時間序列,并且存在短期自相關性。
(a) 樣本路徑 (b) ACF
GBAR(1)模型在四組參數(shù)下CLS 估計、CML估計和BYE估計的模擬結(jié)果分別見表1~表4。
表1 第1組參數(shù)的估計結(jié)果
表2 第2組參數(shù)的估計結(jié)果
表3 第3組參數(shù)的估計結(jié)果
表4 第4組參數(shù)的估計結(jié)果
為了展現(xiàn)三種估計效果的好壞,我們用偏差 (bias) 和均方誤差(MSE)兩種度量指標。通過表中模擬結(jié)果可以看出,CLS、CML和BYE三種估計的偏差和均方誤差都比較小,這表明以上三種估計方法的效果都比較合適。并且隨著樣本量的增大, 估計量的偏差和均方誤差都在減小,估計效果越來越好。其次,在樣本量相對比較小時,貝葉斯估計的偏差和均方誤差比條件最大似然估計的略小一些,而當樣本量增加后,貝葉斯估計的效果和條件最大似然估計的效果相差也不大,且貝葉斯估計的變化幅度相對較小, 后續(xù)我們通過對數(shù)據(jù)污染情況下參數(shù)估計的比較更是驗證了貝葉斯估計的穩(wěn)健性。
為了對三種估計方法做進一步的比較,接下來考慮三種估計方法對帶污染的GBAR(1)模型的估計效果。
定義2若一個隨機過程{Mt}t∈N滿足方程
Mt=(1-δt)Xt+δtξt,t∈N,
式中:{Xt}t∈N----一個沒有被污染的GBAR(1) 過程;
{ξt}t∈N,{δt}t∈N----一列i.i.d.的隨機變量序列,分別服從B(n,p)和B(1,r);
r----污染比例,該模型為帶污染的GBAR(1)模型。
樣本量為100和200時,污染比例r=0.3,0.5時GBAR(1)模型的三種估計方法模擬結(jié)果分別見表5~表6。
從表中可知,對于有污染的數(shù)據(jù),貝葉斯估計的效果要優(yōu)于CML和CLS估計。由此可知,當數(shù)據(jù)存在污染時,貝葉斯估計方法更為可靠。進而說明,貝葉斯估計更加穩(wěn)健,估計結(jié)果合理有效。
實際生活中,有上限范圍的整數(shù)值時間序列數(shù)據(jù)是很常見的,例如:某一地區(qū)一周七天下雨天數(shù)數(shù)據(jù)、全國每個月有多少城市發(fā)生地震數(shù)據(jù)等。文中從天氣網(wǎng)中查詢了2020年1月1日至2021年12月31日湖南省長沙市的歷史天氣,然后統(tǒng)計其每周下雨天數(shù)的數(shù)據(jù)情況,通過整理得到了有上限范圍的時間序列數(shù)據(jù),其中上限范圍n=7,樣本量為105,數(shù)據(jù)的相應統(tǒng)計信息見表7 。
表7 基本統(tǒng)計量
周降雨天數(shù)數(shù)據(jù)的樣本路徑和自相關圖如圖2所示。
(a) 樣本路徑 (b) ACF
從圖2(a)可以看出,數(shù)據(jù)沒有向上或向下的趨勢,而是在 3 附近波動,說明數(shù)據(jù)是平穩(wěn)的。觀察圖2(b)可以看出,自相關系數(shù)在2階后衰減到了標準差以內(nèi),也能確定序列是平穩(wěn)時間序列。并且序列具有自相關性,因此用上限n=7的GBAR(1)模型來刻畫這組數(shù)據(jù)是合理的。
考慮將GBAR(1)模型應用于該數(shù)據(jù),并與文獻[7] 提出的BAR(1)模型、文獻[21]提出的衍生模型、BAR(1)、(IZ-BAR(1))模型和(ZIB-AR(1))模型進行比較。其中,GBAR(1)模型的貝葉斯估計中超參數(shù)設置為e=5,f=5,p=2,q=3,h=6,g=4,對應樣本平均值法的N取106。
通過將GBAR(1)模型與各種其他模型參數(shù)的赤池信息準則(AIC)和貝葉斯信息準則(BIC)進行對比,用來說明GBAR(1)模型擬合的優(yōu)良性。
一方面,無論是GBAR(1)模型的貝葉斯估計還是條件最大似然估計,其AIC、BIC值都小于BAR(1)模型、IZ-BAR(1)模型和ZIB-AR(1)模型,從而得出GBAR(1)模型的擬合效果優(yōu)于BAR(1)模型以及其衍生模型;另一方面,GBAR(1)模型的貝葉斯估計對于樣本量接近100的數(shù)據(jù)估計結(jié)果和AIC、BIC值與條件最大似然估計的相差很小,說明應用貝葉斯方法估計模型參數(shù)也能得到一個比較好的擬合效果,GBAR(1) 模型與其他模型AIC、BIC比較見表8。
表8 GBAR(1) 模型與其他模型AIC、BIC比較
進一步對GBAR(1)模型的皮爾遜殘差進行分析,其中皮爾遜殘差公式為
根據(jù)GBAR(1)模型的條件期望和條件方差公式,我們分別計算了GBAR(1)模型在貝葉斯和條件最大似然兩種估計方法下殘差的均值和方差。認為擬合殘差的均值和方差越接近0和1,模型的擬合效果越好。在本節(jié)實例數(shù)據(jù)中,所擬合的貝葉斯方法下GBAR(1)模型殘差的均值和方差分別為0.025 6和1.293 3。該結(jié)果也再次說明GBAR(1)模型擬合這組數(shù)據(jù)是合適的。
為了進一步研究具有相依結(jié)構(gòu)的廣義二項自回歸模型的參數(shù)估計問題,充分利用樣本信息、總體信息,以及通過歷史經(jīng)驗或者規(guī)律得到的參數(shù)的先驗信息,求得GBAR(1)模型參數(shù)的貝葉斯估計表示式。證明了對于GBAR(1)模型在小樣本情況下,貝葉斯估計的效果要優(yōu)于條件最大似然估計和條件最小二乘估計。通過對帶污染的GBAR(1)模型進行數(shù)值模擬,驗證了貝葉斯估計方法的穩(wěn)健性。最后將GBAR(1)模型應用于長沙市下雨天數(shù)據(jù)集進行擬合,與幾類二項自回歸模型進行比較,得到最小的AIC與BIC值。表明GBAR(1)模型的擬合效果優(yōu)于BAR(1)模型以及其衍生模型。