黎小麗 李時敏,? 孔盈瑩

(1.廣東財經(jīng)大學(xué)統(tǒng)計與數(shù)學(xué)學(xué)院,廣州,510320；2.廣東財經(jīng)大學(xué)會計學(xué)院,廣州,510320)

1 引言

在大數(shù)據(jù)時代,數(shù)據(jù)特征表現(xiàn)的越來越復(fù)雜多樣.正態(tài)分布只是一種理論上的分布,實際數(shù)據(jù)往往呈現(xiàn)偏斜分布,并且至多是近似正態(tài)分布[1].但偏斜分布可能對統(tǒng)計量的穩(wěn)態(tài)性產(chǎn)生致命的影響.第一,由于普通的均值和方差不具備抵抗異常值的性能,從而使得實際樣本數(shù)據(jù)失去代表性,模型產(chǎn)生偏差,最終導(dǎo)致得到的結(jié)果不可靠.第二,異常值的出現(xiàn)會擴大或縮小控制圖的控制限間距,從而使得在統(tǒng)計過程控制(Statistic Process Control,簡稱SPC)中出現(xiàn)虛發(fā)報警或漏發(fā)報警的現(xiàn)象.這不僅增加SPC監(jiān)控成本,也不利于產(chǎn)品設(shè)計.以往學(xué)者常利用異常值診斷來直接剔除異常值,再采用經(jīng)典估計方法進行建模.但異常值往往隱含著重要的信息,刪除異常值得到的結(jié)論亦可導(dǎo)致決策發(fā)生錯誤.因此,有必要采用穩(wěn)健估計方法避免異常值產(chǎn)生的影響并建立正確的模型.

1953年,Box提出了“穩(wěn)健性”的定義,后來經(jīng)過研究得出了大量的穩(wěn)健估計方法.其中Huber,Hampel,Rousseeuw和Tukey等對穩(wěn)健估計進行了深入的研究.周江文[2]提出了等價權(quán)的概念,將M估計最小二乘化使得傳統(tǒng)最小二乘法具備了抗差能力,并提出了兩種有效的估計方案——IGG方案和IGGII方案.楊元喜[3]對等價權(quán)原理進行了擴充,提出了IGGIII方案.Huber[4]也給出了“穩(wěn)健估計”的具體定義.此后,穩(wěn)健估計方法大量應(yīng)用于各個領(lǐng)域.

郭亞帆[5]深入剖析了當(dāng)實際樣本數(shù)據(jù)并不服從正態(tài)分布時,幾種常用統(tǒng)計量并不具備抵御異常值的能力,并指出若注重于正態(tài)分布附近的較高效率使用切尾均值剔除較重尾分布的數(shù)據(jù)則選用Huber M統(tǒng)計量.脫立文等[6]則利用誤差分布律的概率分布函數(shù)來定義權(quán)函數(shù),并通過實例證明該方法的穩(wěn)健優(yōu)越性.Christophe C等[7]基于Holt-Winters預(yù)測方法對時間序列數(shù)據(jù)構(gòu)造的控制圖作了深入的研究,使用可靠的技術(shù)防止異常點影響控制圖的控制限,并構(gòu)造了穩(wěn)健的時間序列控制圖.王曉輝等[8]將常用的13種穩(wěn)健估計方法應(yīng)用于概率積分參數(shù)的選取,并表明其中三種穩(wěn)健估計方法在該應(yīng)用的顯著成效,其中包括IGGIII權(quán)函數(shù).

姜佃高等[9]通過仿真實驗也對常用的13種穩(wěn)健估計方法進行比較,并應(yīng)用于多元線性回歸中,得出其中兩種方法在多元線性回歸中具有相對更為有效的穩(wěn)健估計效果.蘇擁英等[10]采用Huber權(quán)函數(shù)對時間序列中不同的點施加不同權(quán)重,構(gòu)建穩(wěn)健自相關(guān)函數(shù),并實證分析得出該方法具有良好的穩(wěn)健性效果.李雄英等[11]基于Rousseeuw[12]提出的Fast-MCD方法建立穩(wěn)健主成分聚類算法,并結(jié)合模擬分析和實際含異常值的高維數(shù)據(jù)進行案例分析,表明穩(wěn)健主成分聚類算法可有效地抵御異常值產(chǎn)生的影響.

在SPC過程中,質(zhì)量特征值常表現(xiàn)出自相關(guān)現(xiàn)象,違背獨立性假設(shè).為了解決過程自相關(guān)情況下質(zhì)量控制問題,統(tǒng)計學(xué)家們提出了不同的處理方法,其中主要方法之一就是引入時間序列分析法.Alwan等[13]等提出如果在過程平穩(wěn)前提下自相關(guān)的樣本數(shù)據(jù)擬合準(zhǔn)確的時間序列模型,則得到的殘差序列相互獨立.因此在SPC中可使用殘差序列構(gòu)建殘差控制圖進行過程監(jiān)控.

孫靜[14]提出對于存在大量自相關(guān)現(xiàn)象的海量數(shù)據(jù),使用殘差控制圖替代傳統(tǒng)控制圖可得到更好的控制效果.王斌會等[15]通過運用蒙特卡洛模擬方法進一步研究質(zhì)量過程自相關(guān)條件下的殘差控制圖并與傳統(tǒng)控制圖比較,表明殘差控制圖能夠有效避免在受控狀態(tài)下虛發(fā)警報的現(xiàn)象.王志堅[16]則利用Hampel權(quán)函數(shù)來建立穩(wěn)健的ARMA模型進而得到獨立同分布的殘差序列,從而構(gòu)造出對異常值具有強抗干擾性的穩(wěn)健殘差控制圖.

本文在探究穩(wěn)健估計方法IGGIII權(quán)函數(shù)的最適調(diào)和參數(shù)的基礎(chǔ)上,構(gòu)建穩(wěn)健自協(xié)方差函數(shù),對含有異常值的時間序列進行穩(wěn)健變換后,進而建立穩(wěn)健的ARMA模型,從而得到獨立同分布的殘差序列,最后結(jié)合均值-方差穩(wěn)健估計量來構(gòu)造控制圖的控制限,構(gòu)造出穩(wěn)健的ARMA殘差控制圖.模擬和實證分析表明本文提出的基于IGGIII權(quán)函數(shù)構(gòu)造的穩(wěn)健ARMA殘差控制圖具有更好的抗差性.

2 原理與方法

2.1 傳統(tǒng)ARMA控制圖的不穩(wěn)健性

假設(shè){rt}是一個自相關(guān)的序列,當(dāng)它平穩(wěn)時,可通過如下的ARMA模型對其進行擬合:

當(dāng)隨機過程Xt=μ+εt,εt～N(0,σ2)時,μ=E(Xt),且μ,σ分別為過程的均值和標(biāo)準(zhǔn)差.傳統(tǒng)ARMA控制圖的控制限是采用3σ原則設(shè)置的,其控制上限(UCL)、中心限(CL)、控制下限(LCL)如下

然而,在實際中,很多數(shù)據(jù)往往呈現(xiàn)為非正態(tài)的有偏分布,并且往往會由于一些突發(fā)因素導(dǎo)致數(shù)據(jù)中出現(xiàn)一些極端值.如在金融時間序列數(shù)據(jù)中,金融危機、政策的頒布和在SPC監(jiān)控中機器發(fā)生故障、人為操作記錄錯誤等因素都可導(dǎo)致異常值的產(chǎn)生.若對這些含有異常值的數(shù)據(jù)采用傳統(tǒng)的方法建模、進行參數(shù)估計以及預(yù)測無疑會影響結(jié)果的準(zhǔn)確性.

為了說明異常值對傳統(tǒng)控制圖的影響,我們首先生成100個標(biāo)準(zhǔn)正態(tài)分布的隨機數(shù),使用傳統(tǒng)控制圖檢測出樣第46號數(shù)據(jù)不在控制上下限之間(圖1).然后將第25號和第26號數(shù)據(jù)分別換成異常值3.5和4(圖2),此時第46號數(shù)據(jù)卻沒有超出控制上下限,出現(xiàn)了“遮蔽效應(yīng)”.這是由于替換的異常值使得原始序列中的均值和標(biāo)準(zhǔn)差都發(fā)生了變化,從而導(dǎo)致控制限的間距被拉大,出現(xiàn)了漏報情況.上述現(xiàn)象表明傳統(tǒng)控制圖易受異常值的影響.由此可見,數(shù)據(jù)中異常值的出現(xiàn)可能會導(dǎo)致過程監(jiān)控中產(chǎn)生漏發(fā)報警的現(xiàn)象.因此,有必要對傳統(tǒng)控制圖進行穩(wěn)健改進.

圖1 傳統(tǒng)控制圖檢測出46號異常值

圖2 傳統(tǒng)控制圖未檢測出46號異常值

2.2 基于IGGIII權(quán)函數(shù)的樣本自相關(guān)函數(shù)

Box等[17]指出異常值對時間序列模型的識別、參數(shù)估計、診斷檢驗甚至預(yù)測都會產(chǎn)生重大影響.異常值的出現(xiàn)導(dǎo)致殘差序列不服從獨立同分布,從而使傳統(tǒng)殘差控制圖對異常值的檢測性能失效.這不僅不利于在SPC中有效控制產(chǎn)品質(zhì)量,還給企業(yè)造成在監(jiān)控過程中經(jīng)濟成本的提升.

ARMA模型是擬合時間序列常用的模型,矩估計(也稱Yule-Walker估計)是對ARMA模型常用的參數(shù)估計方法之一,其基本思想是用樣本自協(xié)方差函數(shù)W(k)來估計模型的參數(shù).設(shè)時間序列為{Xt},則其自相關(guān)函數(shù)為

當(dāng)時間序列中存在異常值時,自協(xié)方差函數(shù)是不穩(wěn)健的,這會導(dǎo)致參數(shù)估計值偏離實際值.由于樣本自協(xié)方差函數(shù)W(k)是對ARMA模型的參數(shù)進行矩估計的核心,因此,為了避免異常值的影響,本文采用楊元喜[13]提出的IGGIII權(quán)函數(shù)對樣本自協(xié)方差函數(shù)W(k)進行穩(wěn)健改進.IGGIII權(quán)函數(shù)的表達式如下：

其中,εt表示第t個數(shù)據(jù)對應(yīng)的殘差值,σ表示標(biāo)準(zhǔn)差,b,c是調(diào)節(jié)系數(shù),分別取值1～1.5,2.5～3.通過模擬分析可知權(quán)函數(shù)對于調(diào)節(jié)系數(shù)的確定也是敏感的.郭立志[18]等通過探究含有粗差的基坑位移變化量數(shù)據(jù),得到IGGIII權(quán)函數(shù)的調(diào)節(jié)系數(shù)均在0.928～1.136范圍時,穩(wěn)健估計的抗差效果更為明顯.我們通過不斷選取IGGIII權(quán)函數(shù)的調(diào)節(jié)系數(shù)作數(shù)值模擬(模擬達1000次以上),發(fā)現(xiàn)當(dāng)調(diào)節(jié)系數(shù)選取為b=1,c=3時,構(gòu)造的控制圖的控制限更穩(wěn)健,間距沒有發(fā)生很大變化,并且能夠把數(shù)據(jù)中的所有異常值檢測出.

根據(jù)以往學(xué)者的研究,發(fā)現(xiàn)不同的權(quán)函數(shù)的選取會產(chǎn)生不一樣的穩(wěn)健估計效果.權(quán)函數(shù)對于越小的殘差絕對值賦予的權(quán)重越大,對于越大的殘差絕對值則賦予的權(quán)重越小.從幾何角度來解釋,權(quán)函數(shù)壓縮了數(shù)據(jù)中潛在的異常值.采用上述IGGIII權(quán)函數(shù)對自協(xié)方差函數(shù)進行改進,得到穩(wěn)健的樣本自協(xié)方差函數(shù)：

利用它可以構(gòu)建出穩(wěn)健的ARMA模型.

2.3 穩(wěn)健ARMA殘差控制圖的構(gòu)建

由上文知,傳統(tǒng)控制圖的控制限由樣本均值μ和標(biāo)準(zhǔn)差σ確定,而這兩個參數(shù)的估計均對異常值敏感,是導(dǎo)致傳統(tǒng)控制圖不穩(wěn)健的因素之一,因此有必要利用樣本均值μ和標(biāo)準(zhǔn)差σ的穩(wěn)健估計量對控制限進行穩(wěn)健構(gòu)造.

中位數(shù)、切尾均值、縮尾均值與加權(quán)三均值是常用的均值穩(wěn)健估計量.四分位間距、縮尾標(biāo)準(zhǔn)差、絕對離差中位數(shù)與平均絕對離差是常用的標(biāo)準(zhǔn)差估計量.本文通過仿真實驗,對于含有異常值的序列,確定用加權(quán)三均值與平均絕對離差分別對均值μ和標(biāo)準(zhǔn)差σ進行估計而得到的控制圖的控制限具有更好的穩(wěn)健性.對于一個序列x1,x2,···,xN其加權(quán)三均值與平均絕對離差為

其中QL,QU,MED分別表示次序統(tǒng)計量中的樣本下四分位數(shù),上四分位數(shù)和中位數(shù).中位數(shù)由于不受極端值的影響,因此具有穩(wěn)健性.

下面我們給出采用穩(wěn)健IGGIII權(quán)函數(shù)構(gòu)造穩(wěn)健ARMA模型,進而構(gòu)建穩(wěn)健ARMA殘差控制圖的步驟.

記xt為原始序列,通過權(quán)函數(shù)對其加權(quán)后的序列記為

1)根據(jù)xt的中位數(shù),計算出殘差序列,進而得到殘差絕對值序列；

2)根據(jù)殘差絕對值序列的范圍,賦予原序列xt相應(yīng)的權(quán)重w(εt),從而得到加權(quán)的序列

4)根據(jù)穩(wěn)健ARMA模型得到殘差序列,并根據(jù)穩(wěn)健的均值和標(biāo)準(zhǔn)差統(tǒng)計量構(gòu)造控制圖的控制限,從而構(gòu)造出穩(wěn)健ARMA殘差控制圖.

3 模擬分析

為檢驗基于IGGIII權(quán)函數(shù)穩(wěn)健改進后的效果,本節(jié)通過含有異常值的樣本數(shù)據(jù)來說明穩(wěn)健ARMA殘差控制圖比傳統(tǒng)ARMA殘差控制圖的檢測效果更優(yōu).

首先,隨機模擬生成一組樣本量為500的AR(1)且的時間序列數(shù)據(jù),再從這500個數(shù)據(jù)中隨機抽取,分別構(gòu)造污染率為4%和8%的位置污染分布：

其中ξ表述污染率,μ0和σ2分別為AR(1)模型的均值和方差,μ=4.5或?4.2.

然后,對受污染的序列分別進行建模,擬合結(jié)果見表1.

從表1可知,當(dāng)數(shù)據(jù)序列中含有8%的異常值時,通過傳統(tǒng)方法建立的模型,其顯著性較低,而通過穩(wěn)健方法得到的模型仍然保持著較好的顯著性,這說明通過穩(wěn)健改進的ARMA模型對異常值有良好的耐受性.

表1 傳統(tǒng)ARMA(1,1)模型和穩(wěn)健ARMA(1,1)模型的參數(shù)估計結(jié)果比較

此外,當(dāng)序列中含有4%的異常值時,擬合得到的傳統(tǒng)ARMA(1,1)模型的AIC值為1767,穩(wěn)健ARMA(1,1)模型的AIC值為1016；而當(dāng)序列中含有8%的異常值時,擬合得到的傳統(tǒng)ARMA(1,1)模型的AIC值為2068,穩(wěn)健ARMA(1,1)模型的AIC值為1288.通過比較可知,使用穩(wěn)健技術(shù)建立的模型更優(yōu),更有效.

對上述具有不同污染率的數(shù)據(jù),分別建立傳統(tǒng)ARMA(1,1)模型和穩(wěn)健ARMA(1,1)模型然后,構(gòu)建相應(yīng)的ARMA殘差控制圖,結(jié)果如圖3?圖6 所示.

圖3 與圖5是基于傳統(tǒng)方法建立的控制圖,其控制限是基于3σ原則構(gòu)造的.對含有4%(20個)異常值的序列,由圖3及表2可知,控制圖只識別出13個異常值,漏報率為35%.對含有8%(40個)異常值的序列,由圖5及表3可知,傳統(tǒng)ARMA型控制圖只識別出28個異常值,漏報率為30%.這說明傳統(tǒng)控制圖的監(jiān)控性能不佳,主要原因是異常值拉大了上下控制限之間的距離.

圖3 含4%異常值時建立的傳統(tǒng)ARMA殘差控制圖

圖5 含8%異常值時建立的傳統(tǒng)ARMA殘差控制圖

圖4 與圖6是基于穩(wěn)健技術(shù)改進后的方法建立的控制圖.對含有4%(20個)異常值的序列,由圖6及表2可知,控制圖識別出20個異常值,漏報率為0%.對含有8%(40個)異常值的序列,由圖6及表3可知,穩(wěn)健ARMA型控制圖識別出38個異常值,漏報率為5%.

表2 傳統(tǒng)和穩(wěn)健ARMA殘差控制圖在序列含有4%異常值下的監(jiān)控情況

表3 傳統(tǒng)和穩(wěn)健ARMA殘差控制圖在序列含有8%異常值下的監(jiān)控情況

圖4 含4%異常值時建立的穩(wěn)健ARMA殘差控制圖

圖6 含8%異常值時建立的穩(wěn)健ARMA殘差控制圖

綜上,通過穩(wěn)健技術(shù)改進后構(gòu)建的穩(wěn)健ARMA殘差控制圖受異常值的影響較小,對異常值有更好的抗干擾性和抗差性,檢測效果更好.

4 實證檢驗

在證券投資分析中,證券數(shù)據(jù)往往表現(xiàn)出波動聚集性.本節(jié)以證券市場中的實例來檢驗穩(wěn)健時間序列建模方法對傳統(tǒng)殘差控制圖的改進效果.我們選取亞馬遜(AMZN)公司從2013年3月25日至2020年3月25日股票收盤價的日對數(shù)收益率(數(shù)據(jù)來源：雅虎財經(jīng)網(wǎng)站),共1763個觀測數(shù)據(jù).對數(shù)收益率定義如下：

其中Pt為第t天的股票收盤價,Rt為第t天的對數(shù)收益率.首先對AMZN的日對數(shù)收益率進行探索性數(shù)據(jù)分析.由圖7可以看到股票收盤價日對數(shù)收益率的數(shù)據(jù)有波動聚集性的現(xiàn)象,且有大量數(shù)據(jù)顯著偏離均值,故初步判斷數(shù)據(jù)中存在異常值.再由收益率箱線圖(圖8)可知,數(shù)據(jù)偏離正態(tài)分布.

圖7 AMZN對數(shù)收益率時序圖

圖8 AMZN對數(shù)收益率箱線圖

事實上,表4的正態(tài)性檢驗結(jié)果拒絕了AMZN.rtn數(shù)據(jù)服從正態(tài)分布的原假設(shè),但是從表5的平穩(wěn)性檢驗結(jié)果知,AMZN.rtn數(shù)據(jù)序列是平穩(wěn)的.此外,通過圖9中的自相關(guān)圖(ACF)和偏自相關(guān)圖(PACF)可判斷序列AMZN.rtn數(shù)據(jù)序列存在相關(guān)性.我們對其分別建立傳統(tǒng)和穩(wěn)健的時間序列模型如下(擬合結(jié)果見表6)：

表5 序列AMZN.rtn的平穩(wěn)性檢驗結(jié)果

圖9 序列AMZN.rtn的ACF和PACF圖

表4 AMZN.rtn數(shù)據(jù)的正態(tài)性檢驗

由表6知,傳統(tǒng)ARMA(1,1)模型的AIC值為?8990,穩(wěn)健ARMA(1,1)模型的AIC值為?11884,進一步說明穩(wěn)健ARMA(1,1)模型優(yōu)于傳統(tǒng)ARMA(1,1)模型.

表6 傳統(tǒng)和穩(wěn)健ARMA(2.3)模型的參數(shù)估計結(jié)果

接著,根據(jù)擬合得到的模型分別構(gòu)建ARMA殘差控制圖,如圖10,圖11所示.

圖10 傳統(tǒng)ARMA殘差控制圖

圖11 穩(wěn)健ARMA殘差控制圖

由圖10可知,由于數(shù)據(jù)中存在大量異常值導(dǎo)致控制圖的控制上下限之間的距離被拉大,傳統(tǒng)控制圖出現(xiàn)了漏發(fā)警報的現(xiàn)象.而由圖11可知,穩(wěn)健ARMA殘差控制圖能夠有效地檢測出SPC中的異常值,這表明穩(wěn)健ARMA殘差控制圖的檢測效果更優(yōu),對異常值有較好的抗差性和抗干擾性.

進一步,由表7可知,傳統(tǒng)ARMA殘差控制圖只檢測出48個異常值,而穩(wěn)健ARMA殘差控制圖檢測出85個異常值,且傳統(tǒng)ARMA控制圖檢測出來的異常值,都能被穩(wěn)健ARMA殘差控制圖檢測出,說明穩(wěn)健ARMA殘差控制圖的檢測效果確實優(yōu)于傳統(tǒng)ARMA殘差控制圖.

表7 傳統(tǒng)和穩(wěn)健ARMA型殘差控制圖識別出的異常值位置異同

現(xiàn)在,我們再來看看兩種不同的殘差控制圖在幾個具體的異常值點上的表現(xiàn).首先,通過查詢相關(guān)資料可知,在2020年3月11日亞馬遜股票開盤價為1875.85美元,收盤價為1820.86美元,股票價格下跌3.75%,振幅3.69%.這個異常值(第1753樣本點)出現(xiàn)的主要原因是受油價暴跌以及新冠肺炎疫情的影響,美股在該日之前近一個月出現(xiàn)多次熔斷,股市暴跌,導(dǎo)致亞馬遜的股票也受到影響.對此異常值點,傳統(tǒng)控制圖和穩(wěn)健控制圖都能檢測出來.其次,在2015年8月20日亞馬遜因為美國《紐約時報》公開報道批評該公司的員工管理模式,給公司帶來一定的負(fù)面影響,從而導(dǎo)致股價明顯下跌(第607樣本點),對此異常值傳統(tǒng)ARMA殘差控制圖沒有檢測出,而穩(wěn)健ARMA殘差控制圖則檢測出了這個異常值點.最后,在2018年3月27日亞馬遜股價因中美貿(mào)易摩擦也是大幅下跌(第1261樣本點),這一異常值點也被穩(wěn)健ARMA殘差控制圖檢測出,而沒有被傳統(tǒng)ARMA殘差控制圖檢測出.

可見,本文基于改進的IGGIII權(quán)函數(shù)對時序數(shù)據(jù)進行處理,可有效減少異常值對模型的影響,從而使構(gòu)造的穩(wěn)健ARMA殘差控制圖對異常值具有更好的耐受性.

5 小結(jié)

本文針對平穩(wěn)時間序列,首先說明了異常值的存在對傳統(tǒng)建模方法以及構(gòu)造控制圖產(chǎn)生的影響,然后采用IGGIII權(quán)函數(shù)構(gòu)建穩(wěn)健自協(xié)方差函數(shù),對含有異常值的序列進行穩(wěn)健變換,利用變換后的數(shù)據(jù)建立穩(wěn)健的ARAM模型.最后,利用均值和標(biāo)準(zhǔn)差的穩(wěn)健估計計算出控制圖的控制限,進而構(gòu)造出穩(wěn)健ARMA殘差控制圖.模擬和實證表明基于本文所構(gòu)建的穩(wěn)健ARMA殘差控制圖比傳統(tǒng)ARMA殘差控制圖具有更好的抗差性.