李智航，馬岱君

(廣西師范大學(xué) 數(shù)學(xué)與統(tǒng)計(jì)學(xué)院，廣西 桂林 541004)

0 引言

變點(diǎn)問(wèn)題在統(tǒng)計(jì)學(xué)中是一個(gè)較熱門的課題.所謂變點(diǎn)，即是針對(duì)某一序列，在某點(diǎn)前后序列來(lái)自的總體是不同的，此點(diǎn)即為變點(diǎn).變點(diǎn)問(wèn)題廣泛適用于實(shí)際生活中的許多領(lǐng)域，如經(jīng)濟(jì)、地質(zhì)、生物醫(yī)學(xué)等等.自Page在文獻(xiàn)[1]中研究生產(chǎn)線上的生產(chǎn)穩(wěn)定性后，有越來(lái)越多的學(xué)者研究變點(diǎn)問(wèn)題.在變點(diǎn)問(wèn)題伊始，主要是運(yùn)用參數(shù)方法對(duì)變點(diǎn)進(jìn)行研究.陳希孺院士在文獻(xiàn)[2]中對(duì)參數(shù)變點(diǎn)問(wèn)題給出了系統(tǒng)的研究，Cs?rg?和Horváth在文獻(xiàn)[3]中對(duì)變點(diǎn)問(wèn)題的極限情況方面給出了較為系統(tǒng)的結(jié)論.實(shí)際中，很多時(shí)候模型的分布我們不甚了解，于是非參的方法逐漸引起學(xué)者的關(guān)注并形成以非參的方法處理變點(diǎn)模型為主的趨勢(shì)，其中便有不少學(xué)者利用經(jīng)驗(yàn)似然方法對(duì)變點(diǎn)問(wèn)題進(jìn)行研究.

Owen在文獻(xiàn)[4]提出經(jīng)驗(yàn)似然(一種非參數(shù)統(tǒng)計(jì)推斷方法)以及在文獻(xiàn)[5]進(jìn)行總結(jié)整理.經(jīng)驗(yàn)似然具有非參方法的穩(wěn)健性，相較于經(jīng)典的統(tǒng)計(jì)方法有著不少突出的優(yōu)勢(shì)，像Bartlctt糾偏性，域保持性，變換不變性以及其置信域形態(tài)由樣本決定等等.因此，越來(lái)越多學(xué)者對(duì)其應(yīng)用及理論進(jìn)行研究，其中有部分學(xué)者把經(jīng)驗(yàn)似然方法用到變點(diǎn)模型中，Einmahl和Mc Keague在文獻(xiàn)[6]中提及了用經(jīng)驗(yàn)似然的方法去檢驗(yàn)變點(diǎn)；而Liu，Zou和Zhang在文獻(xiàn)[7]運(yùn)用經(jīng)驗(yàn)似然的方法研究線性回歸模型系數(shù)變點(diǎn)問(wèn)題，證明了在沒(méi)有變點(diǎn)的情況下，他們構(gòu)建的統(tǒng)計(jì)量的漸近分布與Gumbel極值分布有關(guān)并且證明了變點(diǎn)存在時(shí)，在一定正則條件下該方法得到的變點(diǎn)估計(jì)具有相合性；Ning在文獻(xiàn)[8]也運(yùn)用經(jīng)驗(yàn)似然檢測(cè)帶有線性趨勢(shì)的均值變點(diǎn)；Liu和Qian在文獻(xiàn)[9]運(yùn)用經(jīng)驗(yàn)似然方法對(duì)分段簡(jiǎn)單線性模型進(jìn)行變點(diǎn)檢測(cè)，模擬出其檢驗(yàn)統(tǒng)計(jì)量于零假設(shè)下的極限分布為Gumbel極值分布以及變點(diǎn)估計(jì)的相合性，同時(shí)還給出了回歸系數(shù)的估計(jì)；李云霞和劉偉棠在文獻(xiàn)[10]利用經(jīng)驗(yàn)似然方法研究Logistic回歸模型的變點(diǎn)檢驗(yàn)，給出其經(jīng)驗(yàn)對(duì)數(shù)似然比統(tǒng)計(jì)量的極值分布，變點(diǎn)的估計(jì)及其相合性以及通過(guò)模擬與實(shí)例證明方法的可行性.Owen在文獻(xiàn)[4-5]提出用經(jīng)驗(yàn)歐氏似然代替經(jīng)驗(yàn)似然，由此減弱計(jì)算的復(fù)雜性，羅旭在文獻(xiàn)[11-12]證明了經(jīng)驗(yàn)歐氏似然與經(jīng)驗(yàn)似然具有類似的大樣本性質(zhì)，但是計(jì)算上較為簡(jiǎn)單，因此，本文將基于經(jīng)驗(yàn)歐氏似然方法對(duì)多元線性回歸模型中的系數(shù)變點(diǎn)進(jìn)行檢測(cè)，結(jié)合文獻(xiàn)[7]與文獻(xiàn)[9]構(gòu)建檢驗(yàn)函數(shù)的思想來(lái)構(gòu)造經(jīng)驗(yàn)歐氏似然比檢驗(yàn)函數(shù)并尋找其漸進(jìn)分布，同時(shí)討論由此得到的變點(diǎn)估計(jì)的相合性，最后通過(guò)數(shù)值模擬及實(shí)例進(jìn)一步對(duì)本文提出的方法進(jìn)行說(shuō)明.

本文安排如下：第一節(jié)介紹多元線性回歸系數(shù)單變點(diǎn)模型，給出經(jīng)驗(yàn)歐氏似然比檢驗(yàn)函數(shù)；第二節(jié)給出主要理論結(jié)果及證明；第三節(jié)為數(shù)值模擬結(jié)果；第四節(jié)進(jìn)行實(shí)例分析；第五節(jié)對(duì)本文進(jìn)行總結(jié).

1 模型與方法

考慮如下線性回歸系數(shù)單變點(diǎn)模型：

其中Xi∈Rd均為來(lái)自總體分布為G的獨(dú)立隨機(jī)向量，α，β(β≠α)∈Rd為未知參數(shù)；k0為未知變點(diǎn)，εi為獨(dú)立同分布隨機(jī)變量，且滿足Eεi=0=σ2< ∞，E(Xiεi)=0.當(dāng)沒(méi)有變點(diǎn)時(shí)，顯然有k0≥n成立；當(dāng)變點(diǎn)存在時(shí)則有1≤k0

對(duì)于回歸模型(1)，當(dāng)沒(méi)有變點(diǎn)時(shí)，記E(Xi(Yi-α))=0，而其回歸系數(shù)的最小二乘估計(jì)如下：

其中X=(X1，X2，…，Xn)T，Y=(Y1，Y2，…，Yn)T，此時(shí)其擬合殘差^Yi=XiT^α對(duì)所有i應(yīng)滿足E(Xi(Yi-^Yi))=0，記Zi)=Xi(Yi-，于是對(duì)假設(shè)檢驗(yàn)問(wèn)題(2)和每個(gè)固定的k，我們可以構(gòu)建如下兩樣本經(jīng)驗(yàn)歐氏似然比檢驗(yàn)函數(shù)：

其中pi≥0分別為對(duì)應(yīng)的概率質(zhì)量，=1，由拉格朗日乘子法容易求得(具體過(guò)程可參考文獻(xiàn)[12])：

從(3)式來(lái)看，對(duì)整個(gè)序列{Zi而言，若不存在變點(diǎn)，都應(yīng)該是收斂到0的，故有較大的概率使得(3)式較小；相反，若存在變點(diǎn)顯然不收斂到0，故較大的概率使得(3)式較大；由于k是未知的，所以我們選取如下的檢驗(yàn)函數(shù)來(lái)對(duì)假設(shè)檢驗(yàn)問(wèn)題(2)做檢驗(yàn)：

但是當(dāng)k或n-k較小時(shí)，例如k

當(dāng)Mn足夠大，大于檢驗(yàn)的臨界值時(shí)，我們拒絕原假設(shè)，認(rèn)為變點(diǎn)存在，此時(shí)給出變點(diǎn)位置的估計(jì)如下：

關(guān)于Mn在零假設(shè)下的極限分布，以及由(5)式給出的變點(diǎn)位置估計(jì)的極限性質(zhì)，我們?cè)谙乱还?jié)中給出.

2 主要理論結(jié)果

第一個(gè)定理給出檢驗(yàn)統(tǒng)計(jì)量Mn在原假設(shè)成立時(shí)的極限分布.

定理1 若存在δ∈ (0，1)使得E‖Zi(α)‖2(1+δ)=E‖Xiεi‖2(1+δ)< ∞，E‖Xi‖4< ∞ 且E(XXT)=∑1，∑1正定，當(dāng)H0成立時(shí)，有

證明：當(dāng)原假設(shè)成立時(shí)，由條件E‖Zi(α)‖2(1+δ)=E‖Xiεi‖2(1+δ)< ∞ 及重對(duì)數(shù)率有

由(6)式及重對(duì)數(shù)率我們有

由上兩式和條件E‖Xi‖4<∞以及Marcinkewicz-Zygmund強(qiáng)大數(shù)率我們有

類似地有

注意到最小二乘估計(jì)的性質(zhì)，k Z1k+(n-k)Z2k=0，于是有

又因?yàn)?/p>

其中K=min(1/10，δ/2(1+2δ))，于是類似文獻(xiàn)[3]定理1.3.1證明中便可得到定理1，與其有區(qū)別的是，我們的結(jié)論由文獻(xiàn)[3]定理A.3.4導(dǎo)出，而非推論A.3.1，定理證畢.

定理1告訴我們Mn的極限分布為極值分布，給定顯著性水平，利用此分布我們可以給出假設(shè)檢驗(yàn)問(wèn)題(2)的漸近拒絕域，下面的定理則是有關(guān)由(5)式給出的變點(diǎn)位置估計(jì)的相合性.記k0為真實(shí)變點(diǎn)，

定理2若E‖Xiεi‖2<∞，E‖Xi‖4<∞，E(XiXTi)是非退化的，存在τ0∈(0，1)使得|τn0-τ0|=O(1)，則當(dāng)H1成立時(shí)，有|-τn0|=op(1)，且(log log n)

證明：由條件E‖Xiεi‖2< ∞，|τn0-τ0|=O(1)容易有

不妨先設(shè)k≤k0，注意k≥[n]，此時(shí)有

又由(7)式，條件E‖Xiεi‖2<∞，E‖Xi‖4<∞及大數(shù)定律，我們有如下式子成立

定理2告訴我們?cè)谝欢l件下，由(5)式給出的變點(diǎn)位置的估計(jì)是具有相合性的，同時(shí)，在備擇假設(shè)成立時(shí)相對(duì)于原假設(shè)成立時(shí)，統(tǒng)計(jì)量Mn具有更高階的發(fā)散速度的，由此可知，用Mn來(lái)對(duì)假設(shè)檢驗(yàn)問(wèn)題(2)做檢驗(yàn)是合適的.

3 數(shù)值模擬

本節(jié)給出數(shù)值模擬的結(jié)果，為方便與文獻(xiàn)[7]做比較，我們采用與文[7]相同的模擬條件.考慮如下簡(jiǎn)單線性回歸模型：

取n=300，顯著性水平α=0.05，xi服從標(biāo)準(zhǔn)正態(tài)分布，對(duì)于εi我們考慮四種分布，分別是εi=N(0，1)，εi=ζ(1)-1，εi=(2)-1(χ2(4)-4)，εi)-1t(4)，其中N為正態(tài)分布，ζ為指數(shù)分布，χ2為卡方，t為t分布；k0=[0.25n]，[0.33n]，[0.5n]，每種情況下重復(fù)模擬次數(shù)為M=10000次.為了比較的公平，對(duì)于每種情況我們通過(guò)模擬得到Mn在零假設(shè)下的經(jīng)驗(yàn)分布的上0.05分位數(shù)作為臨界值，然后計(jì)算其檢驗(yàn)功效.此外，為了考察變點(diǎn)位置估計(jì)的效果，在每種εi的分布及k0取值下，n分別取300，400，600，每種情況同樣重復(fù)模擬M次，然后計(jì)算估計(jì)值離真實(shí)值的平均絕對(duì)值距離D，即

模擬結(jié)果見(jiàn)表1和表2.

表1 檢驗(yàn)功效比較Tab.1 Comparisons of tests powers

表2 變點(diǎn)估計(jì)效果表Tab.2 Effect of change point estimation

表1為3種方法的功效比較表，其中EELR為經(jīng)驗(yàn)歐氏似然比方法，即本文提出的方法，ELR為文[7]中的經(jīng)驗(yàn)似然方法，UI為文獻(xiàn)[3]中的Union-intersection方法，而ELR和UI的模擬結(jié)果均來(lái)源于文獻(xiàn)[7]，三種方法的臨界值均用零假設(shè)下經(jīng)驗(yàn)分布的0.95分位數(shù)確定.從表1可以看出，在干擾項(xiàng)εi是正態(tài)分布時(shí)，EELR與ELR和UI方法有著差不多的表現(xiàn)，在干擾項(xiàng)是其余分布時(shí)，EELR比UI有著更好的表現(xiàn)，但是在指數(shù)分布和卡方分布時(shí)，變點(diǎn)位置k0為75，100時(shí)，EELR方法比ELR方法表現(xiàn)較差，但是k0=150時(shí)，EELR與ELR有差不多的表現(xiàn)，而在干擾項(xiàng)是t分布時(shí)，EELR有著比ELR更好的表現(xiàn).由此可以看出對(duì)于相對(duì)厚尾的分布，EELR的檢驗(yàn)功效會(huì)更好.

表2中的數(shù)據(jù)為對(duì)應(yīng)情況下的D，而括號(hào)里面的則是相應(yīng)的標(biāo)準(zhǔn)差.由表2可以看出，無(wú)論干擾項(xiàng)服從什么分布，變點(diǎn)位置如何選取，D及其相應(yīng)的標(biāo)準(zhǔn)差都表現(xiàn)出隨著樣本量增大而變??；對(duì)于相同的變點(diǎn)位置及樣本量，無(wú)論干擾項(xiàng)是什么分布，D值都比較接近，由此可以認(rèn)為干擾項(xiàng)的分布對(duì)變點(diǎn)位置估計(jì)的效果影響不大，綜合來(lái)看，本文提出的方法對(duì)于變點(diǎn)位置的估計(jì)有著較好的效果，也驗(yàn)證了前一節(jié)提到的，該估計(jì)具有相合性的說(shuō)法.

4 實(shí)例分析

我們選用的實(shí)例是美國(guó)黃石國(guó)家公園噴泉，該例子被文獻(xiàn)[7]與文獻(xiàn)[13]等研究過(guò)，美國(guó)黃石國(guó)家公園噴泉是間歇性噴泉，每隔一段時(shí)間會(huì)噴發(fā)一次，文獻(xiàn)[14]將1980年10月美國(guó)黃石國(guó)家公園噴泉數(shù)據(jù)畫成散點(diǎn)圖，其中橫軸x為噴發(fā)持續(xù)時(shí)間；縱軸y為噴發(fā)時(shí)間間隔.我們將該數(shù)據(jù)(數(shù)據(jù)來(lái)源于文獻(xiàn)[15])同樣畫出散點(diǎn)圖，即圖1，其中x軸單位為秒，y軸單位為分鐘，一共有270組數(shù)據(jù).由圖1可知，該噴泉噴發(fā)的持續(xù)時(shí)間大致分為兩部分；同時(shí)可以用線性回歸模型對(duì)持續(xù)時(shí)間與時(shí)間間隔進(jìn)行擬合，其中持續(xù)時(shí)間為自變量，時(shí)間間隔為因變量.文獻(xiàn)[13]利用參數(shù)似然方法對(duì)該數(shù)據(jù)進(jìn)行研究，得到在沒(méi)有回歸系數(shù)變點(diǎn)的原假設(shè)下p值為0.17的結(jié)論，文獻(xiàn)[7]則利用截?cái)嗟慕?jīng)驗(yàn)似然方法得到原假設(shè)下的p值為0.11，而利用本文提出的方法計(jì)算得出統(tǒng)計(jì)量Mn=7.32，利用定理1計(jì)算出p值為0.22，給定顯著性水平α=0.05，則可認(rèn)為變點(diǎn)不存在，可以用簡(jiǎn)單線性回歸模型對(duì)其進(jìn)行擬合，這與文獻(xiàn)[7，13]的研究結(jié)果一致.

圖1 美國(guó)黃石國(guó)家公園噴泉噴發(fā)時(shí)間間隔與持續(xù)時(shí)間散點(diǎn)圖Fig.1 The scatter plot of geyser eruption interval and duration in Yellowstone National Park,USA

5 結(jié)論

本文基于線性回歸系數(shù)單變點(diǎn)模型的特點(diǎn)，結(jié)合文[7，9]的思想構(gòu)建了截?cái)嗟慕?jīng)驗(yàn)歐氏似然比檢驗(yàn)統(tǒng)計(jì)量，并證明了零假設(shè)成立下該檢驗(yàn)統(tǒng)計(jì)量與文[7]中構(gòu)建的檢驗(yàn)統(tǒng)計(jì)量有著一樣的極限分布；變點(diǎn)存在時(shí)，給出了變點(diǎn)位置的估計(jì)并證明了在一定條件下，該估計(jì)具有相合性；此外，還證明了該檢驗(yàn)統(tǒng)計(jì)量趨于無(wú)窮的速度在有變點(diǎn)時(shí)比沒(méi)有變點(diǎn)時(shí)更快，然后通過(guò)數(shù)值模擬分析說(shuō)明了本文提出的方法的有效性，最后將其運(yùn)用至美國(guó)黃石國(guó)家公園噴泉的例子上，得到與文[7，13]差不多的結(jié)論.