基于IGGIII權函數的穩(wěn)健ARMA殘差控制圖構建

黎小麗 李時敏,? 孔盈瑩

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

1 引言

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

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

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

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

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

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

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

2 原理與方法

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

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

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

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

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

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

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

2.2 基于IGGIII權函數的樣本自相關函數

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

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

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

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

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

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

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

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

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

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

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

記xt為原始序列,通過權函數對其加權后的序列記為

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

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

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

3 模擬分析

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

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

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

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

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

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

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

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

圖3 與圖5是基于傳統(tǒng)方法建立的控制圖,其控制限是基于3σ原則構造的.對含有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)健技術改進后的方法建立的控制圖.對含有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)健技術改進后構建的穩(wěn)健ARMA殘差控制圖受異常值的影響較小,對異常值有更好的抗干擾性和抗差性,檢測效果更好.

4 實證檢驗

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

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

圖7 AMZN對數收益率時序圖

圖8 AMZN對數收益率箱線圖

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

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

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

表4 AMZN.rtn數據的正態(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)模型的參數估計結果

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

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

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

由圖10可知,由于數據中存在大量異常值導致控制圖的控制上下限之間的距離被拉大,傳統(tǒng)控制圖出現了漏發(fā)警報的現象.而由圖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型殘差控制圖識別出的異常值位置異同

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

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

5 小結

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