基于Bootstrap方法的厚尾AR(p)序列均值變點(diǎn)檢驗

喬 瑞，楊云鋒，金 浩

(西安科技大學(xué) 理學(xué)院，陜西 西安 710054)

0 引 言

變點(diǎn)檢驗對于決策者們具有重要的現(xiàn)實(shí)意義,若忽略變點(diǎn)的存在,則會導(dǎo)致錯誤的建模,從而做出錯誤的決策,造成一定的風(fēng)險.為了規(guī)避風(fēng)險,應(yīng)用統(tǒng)計方法對變點(diǎn)進(jìn)行統(tǒng)計推斷就顯得尤為重要.在變點(diǎn)檢驗的各種問題中,均值變點(diǎn)檢驗在變點(diǎn)問題的分析中占據(jù)了非常重要的地位.Gardner[1]討論了方差為1時的獨(dú)立高斯序列均值變點(diǎn)檢驗.隨后Sen等[2]修正了Gardner[1]的檢驗統(tǒng)計量,在方差未知時探討了高斯序列的均值變點(diǎn)問題.Chen等[3]基于比值型監(jiān)測統(tǒng)計量監(jiān)測了長記憶時間序列的均值變點(diǎn)和方差變點(diǎn),發(fā)現(xiàn)此檢驗在越靠近監(jiān)測開始時間時具有更好的有限樣本性能.Lebarbier[4]基于懲罰最小二乘估計方法對高斯序列的多均值變點(diǎn)進(jìn)行了估計.韓四兒等[5]在累積和統(tǒng)計量的基礎(chǔ)上,提出擬似然比方法對多均值變點(diǎn)進(jìn)行檢驗.齊培艷等[6]研究了非參數(shù)回歸模型基于小波系數(shù)的均值變點(diǎn)的在線監(jiān)測問題.與傳統(tǒng)的Kolmogorov-Smirnov檢驗以及“滑窗”檢驗方法相比,自正則的K-S檢驗可以避免窗寬參數(shù)的選取,因此潘婉彬等[7]采用基于自正則的K-S方法對羊群行為的均值變點(diǎn)進(jìn)行了檢驗.

以上的均值變點(diǎn)檢驗問題的模型的新息過程都是高斯序列,要求方差存在.但隨著越來越多的高頻數(shù)據(jù)的出現(xiàn),例如環(huán)境、工業(yè)等數(shù)據(jù)都呈現(xiàn)出尖峰厚尾特征,它們的尾部存在很多異常值,導(dǎo)致大部分信息滯留在尾部,常規(guī)的高斯序列已經(jīng)無法對其進(jìn)行描述,而厚尾序列恰好可以刻畫這類數(shù)據(jù)的特征.因此,對厚尾序列變點(diǎn)的統(tǒng)計推斷是必不可少的.這一內(nèi)容吸引了一眾統(tǒng)計學(xué)家和計量經(jīng)濟(jì)學(xué)家對其進(jìn)行研究.劉艦東等[8]利用累積和統(tǒng)計量研究了方差無窮的ARCH序列均值變點(diǎn)的檢驗問題.針對檢驗統(tǒng)計量依賴未知尾指數(shù)的現(xiàn)象,Jin等[9]采用Bootstrap子抽樣方法檢測了厚尾序列的均值變點(diǎn).Wang等[10]通過利用Bootstrap抽樣分布來逼近檢驗統(tǒng)計量的極限分布,以實(shí)現(xiàn)厚尾序列均值變點(diǎn)的檢驗和估計.Zhou等[11]基于加權(quán)最小二乘估計法構(gòu)造加權(quán)累積和統(tǒng)計量推斷了厚尾序列的均值變點(diǎn),其檢驗統(tǒng)計量的極限分布不依賴于未知的厚尾指數(shù).更多關(guān)于厚尾序列變點(diǎn)的研究,可參見文獻(xiàn)[12-14].

基于Ratio統(tǒng)計量將Jin等[9]對獨(dú)立同分布的厚尾序列均值變點(diǎn)檢驗拓展到一般的AR(p)模型.由于檢驗功效易受變點(diǎn)位置的影響,通過顛倒統(tǒng)計量重新構(gòu)造出新的Ratio統(tǒng)計量來彌補(bǔ)此缺陷.另外,考慮到統(tǒng)計量的極限分布依賴于未知的尾指數(shù),利用Bootstrap抽樣方法來逼近統(tǒng)計量的極限分布從而實(shí)現(xiàn)厚尾序列均值變點(diǎn)的檢驗.

1 模型與假設(shè)

考慮新息過程為AR(p)過程的均值變點(diǎn)模型:

yt=μ+Δ·I{t>k*}+εt

(1)

εt=ρ1εt-1+…+ρpεt-p+ηt,t=1,…,T

(2)

這里考慮假設(shè)檢驗問題:

H0:Δ=0

H1:在未知的時刻k*,Δ≠0

下面給出本文所需的假設(shè)和引理:

假設(shè)1假設(shè)特征多項式ρ(z)=1-ρ1z-…-ρpzp的根都在單位圓之外.

假設(shè)2{ηt}∈D(κ),這里D(κ)為尾指數(shù)κ∈(1,2]的穩(wěn)定吸收域,且Eηt=0.

引理1如果{ηt}是獨(dú)立同分布序列,且滿足假設(shè)2,則:

其中L1(r)和L2(r)分別是在[0,1]上的κ-Lévy過程和κ/2-Lévy過程.Kokoszka等[15]進(jìn)一步證明aT可以簡化為aT=T-κ(T).L1(·)是一個穩(wěn)定過程,其可以表達(dá)為:

這里{Ut}是獨(dú)立同分布在區(qū)間[0,1]上的均勻分布隨機(jī)變量,{δt}是獨(dú)立同分布的隨機(jī)變量,滿足P(δt=1)=1-q,P(δt=-1)=q.Γ1,Γ2,…是具有勒貝格測度的泊松過程的到達(dá)時間,且{Ut,δt,Γt}相互獨(dú)立.

經(jīng)典的累積和統(tǒng)計量已被廣泛用來檢測均值變點(diǎn),但累積和統(tǒng)計量需要通過序列的長期方差來標(biāo)準(zhǔn)化,特別是在新息過程是相依的情形下.Breiman[16]和Antoch等[17]指出找到具有良好性質(zhì)的長期方差估計量是極其困難的,即使在獨(dú)立情形下,長期方差估計量也常常表現(xiàn)不佳,在備擇假設(shè)下顯得尤為突出.因此,Horváth等[18]提出了Ratio統(tǒng)計量:

(3)

定義新的Ratio檢驗統(tǒng)計量如下:

2 主要結(jié)論

下面給出統(tǒng)計量在原假設(shè)下的極限分布及在備擇假設(shè)下的一致性.

定理1若隨機(jī)序列{yt}由(1)～(2)生成,{εt}是厚尾AR(p)過程,且假設(shè)1～2成立,在原假設(shè)H0下,當(dāng)T→∞時,有:

(4)

則式(4)可寫為矩陣形式:

G=Jρ+ξ

因此,ρ的最小二乘估計為:

(5)

將G=Jρ+ξ代入式(5),則:

(6)

(7)

令：

聯(lián)立式(6)和(7),則:

(8)

對于t=p+1,…,k,有：

則：

(9)

同理可證得:

(10)

定理1表明修正的Ratio統(tǒng)計量R*在原假設(shè)下的極限分布是Lévy過程的泛函.下面給出備擇假設(shè)下統(tǒng)計量的一致性結(jié)論.

定理2若隨機(jī)序列{yt}由(1)～(2)生成,{εt}是厚尾AR(p)過程,且假設(shè)1～2成立,在備擇假設(shè)H1下,當(dāng)T→∞時：

證明當(dāng)k*≤k時,基于樣本yt,t=1,…,k,得到殘差:

(11)

(12)

將式(11)代入式(2),有:

其中:

(13)

因式(13)的第二項不是隨機(jī)變量,再次利用BN分解,得:

=Op(T)

(14)

令：

根據(jù)定理1的證明,結(jié)合式(12)和式(14),有:

(15)

因此,當(dāng)t=p+1,…,k時：

則:

令i=[Ts],k=[Tr],k*=[Tr*],當(dāng)k*≤k時,有:

=3-1r*2(r-r*)2r-1Δ2(1-ρ1-…-ρp)2

和

=4-1r*2(r-r*)2r-1Δ2(1-ρ1-…-ρp)2.

從而容易推出:

由于樣本yt,t=k+1,…,T不包含均值變點(diǎn),則統(tǒng)計量分母的極限分布與其在原假設(shè)下的極限分布相一致.因此,(i)得證.類似的,(ii)同理可證.則定理2證畢.

定理2證明了Ratio統(tǒng)計量在備擇假設(shè)下的一致性.不難發(fā)現(xiàn),統(tǒng)計量的發(fā)散速度隨著厚尾指數(shù)κ,樣本容量T、跳躍幅度Δ的增大而增大.特別地,若p=1,則發(fā)散速度隨著回歸系數(shù)ρ1的減小而增大.

統(tǒng)計量的極限分布依賴于未知的尾指數(shù)κ,Mandelbrot[20]提出了一個粗略的方法來估計厚尾指數(shù)κ,但精確性不令人滿意.為了避免干擾參數(shù)κ的估計,利用Bootstrap抽樣方法來逼近原假設(shè)下統(tǒng)計量的極限分布以獲取精確的臨界值.具體步驟如下:

選擇一個通用的m使得其在任何情況下都是最優(yōu)的難以實(shí)現(xiàn),但Mcmurry等[21]提供了控制經(jīng)驗水平和經(jīng)驗勢的一個理想選擇m=[4T/ln(T)].因此,在數(shù)值模擬中這個子樣本m也將繼續(xù)沿用.

3 數(shù)值模擬

本節(jié)通過蒙特卡洛數(shù)值模擬驗證基于Bootstrap方法的Ratio檢驗的有效性.考慮一階自回歸單均值變點(diǎn)模型:yt=μ+Δ·I{t>k*}+εt,εt=ρ1εt-1+ηt,其中{ηt}是零期望獨(dú)立同分布厚尾序列.先利用Bootstrap抽樣方法確定統(tǒng)計量的臨界值.不失一般性,設(shè)定顯著性水平α=0.05,厚尾指數(shù)κ={1.1,…,2},循環(huán)次數(shù)為B=3 000,樣本容量為T=2 000,均值μ=0,回歸系數(shù)ρ1={-0.5,0,0.5}.模擬結(jié)果均通過Matlab軟件實(shí)現(xiàn).

表1 基于Bootstrap的Ratio統(tǒng)計量臨界值,T=2 000

為了進(jìn)一步說明Ratio檢驗統(tǒng)計量的檢驗功效,設(shè)定跳躍幅度Δ={0,2,4},變點(diǎn)時刻r*={0.3,0.7},樣本容量T={300,500,1 000},并用 3 000 次試驗中拒絕原假設(shè)的百分?jǐn)?shù)作為經(jīng)驗勢函數(shù)值.圖1所有的藍(lán)線表示回歸系數(shù)ρ1=-0.5的統(tǒng)計量的拒絕率,綠線表示ρ1=0時的拒絕率,紅線表示ρ1=0.5時的拒絕率,實(shí)線表示R,虛線表示R*,橫坐標(biāo)為厚尾指數(shù)κ,縱坐標(biāo)為拒絕率.

圖1表示原假設(shè)下原統(tǒng)計量R和重新構(gòu)造的統(tǒng)計量R*的拒絕率.由圖中可以看出,當(dāng)樣本容量增大時,拒絕率越來越接近顯著性水平0.05,且在其附近的波動越來越小;厚尾指數(shù)、回歸系數(shù)的變化幾乎不影響Ratio檢驗統(tǒng)計量的拒絕率.這說明基于Bootstrap方法的Ratio檢驗統(tǒng)計量很好地控制了經(jīng)驗水平.

(a)樣本T=300 (b)樣本T=500 (c)樣本T=1 000圖1 原假設(shè)下Ratio檢驗的拒絕率Fig.1 The rejection rate of Ratio test under the null hypothesis

圖2～圖3表示備擇假設(shè)下跳躍幅度Δ=2,變點(diǎn)時刻r*={0.3,0.7}的原統(tǒng)計量R和重新構(gòu)造的統(tǒng)計量R*的拒絕率.隨著樣本容量的增大,統(tǒng)計量的拒絕率增大,這與定理2的結(jié)論相吻合,即樣本容量越大統(tǒng)計量發(fā)散程度越大.當(dāng)厚尾指數(shù)κ減小時,統(tǒng)計量的拒絕率減小.這是由于隨著厚尾指數(shù)κ越小,序列包含的異常值會增多,導(dǎo)致臨界值越大,出現(xiàn)難以拒絕原假設(shè)的現(xiàn)象，從而使得拒絕率減小.此外,κ越小,備擇假設(shè)下統(tǒng)計量的發(fā)散速度越慢,這說明了κ越小,拒絕率越低的合理性.當(dāng)回歸系數(shù)ρ1=-0.5時,拒絕率最大,回歸系數(shù)ρ1=0次之,回歸系數(shù)ρ1=0.5時,拒絕率最小.這可以從定理2的結(jié)論中得出,ρ1越小,統(tǒng)計量越發(fā)散,拒絕率越高.

圖2表明當(dāng)r*=0.3時,原統(tǒng)計量R比重新構(gòu)造的統(tǒng)計量R*的拒絕率略高,這說明若變點(diǎn)時刻位于樣本前半段,原統(tǒng)計量R較容易檢測到變點(diǎn),但兩者之間的差異不大.但由圖3可以看出,若變點(diǎn)時刻位于樣本后半段r*=0.7時,重新構(gòu)造的統(tǒng)計量R*的拒絕率遠(yuǎn)高于原統(tǒng)計量R的拒絕率,且兩者之間差異愈加明顯.這恰好印證了本文所提檢驗統(tǒng)計量的優(yōu)勢.根本原因是相比原統(tǒng)計量的檢驗依賴于變點(diǎn)位置,重新構(gòu)造的Ratio檢驗統(tǒng)計量的經(jīng)驗勢不會因變點(diǎn)時刻位于樣本后半段而大幅度降低,從而減弱檢驗功效.

(a)樣本T=300 (b)樣本T=500 (c)樣本T=1 000圖2 備擇假設(shè)下Ratio檢驗的拒絕率,r*=0.3,Δ=2Fig.2 The rejection rate of Ratio test under the alternative hypothesis,r*=0.3,Δ=2

(a)樣本T=300 (b)樣本T=500 (c)樣本T=1 000圖3 備擇假設(shè)下Ratio檢驗拒絕率,r*=0.7,Δ=2Fig.3 Ratio test rejection rate under the alternative hypothesis,r*=0.7,Δ=2

圖4～圖5表示備擇假設(shè)下跳躍幅度Δ=4,變點(diǎn)時刻r*={0.3,0.7}的原統(tǒng)計量R和重新構(gòu)造的統(tǒng)計量R*的拒絕率.正如所期望的,統(tǒng)計量的檢驗功效依賴于跳躍幅度Δ,當(dāng)跳躍幅度增大時,統(tǒng)計量的拒絕率出現(xiàn)了明顯的增大.與Δ=2的情況類似,隨著樣本容量、厚尾指數(shù)的增大,統(tǒng)計量的拒絕率增大,但回歸系數(shù)越大統(tǒng)計量的拒絕率越小.總之,相比于原統(tǒng)計量R對變點(diǎn)在樣本后半段檢驗不顯著的缺陷,重新構(gòu)造的統(tǒng)計量R*的檢驗顯著性不受變點(diǎn)位置的影響.這表明基于Bootstrap方法的Ratio檢驗為檢測厚尾相依序列均值變點(diǎn)提供了一種行之有效的工具.

(a)樣本T=300 (b)樣本T=500 (c)樣本T=1 000圖4 備擇假設(shè)下Ratio檢驗拒絕率,r*=0.3,Δ=4 Fig.4 Ratio test rejection rate under the alternative hypothesis,r*=0.3,Δ=4

(a)樣本T=300 (b)樣本T=500 (c)樣本T=1 000圖5 備擇假設(shè)下Ratio檢驗拒絕率,r*=0.7,Δ=4 Fig.5 Ratio test rejection rate under the alternative hypothesis,r*=0.7,Δ=4

4 結(jié) 論

本文基于Bootstrap方法檢測了厚尾自回歸過程的均值變點(diǎn).為提升檢驗功效,通過顛倒統(tǒng)計量重新構(gòu)造出新的Ratio檢驗統(tǒng)計量,在原假設(shè)下推導(dǎo)了檢驗統(tǒng)計量的極限分布為Lévy過程的泛函,并給出了備擇假設(shè)下的一致性.為了避免厚尾指數(shù)的估計,Bootstrap抽樣方法被用來確定更精確的臨界值.最后由蒙特卡洛數(shù)值模擬驗證了基于Bootstrap方法的Ratio檢驗的有效性和可行性.