基于相依區(qū)間值數(shù)據(jù)的EWMA FIMSE控制圖

王美涵,陳培樂(lè),潘曉映,施建華,2,3,4*

（1.閩南師范大學(xué)數(shù)學(xué)與統(tǒng)計(jì)學(xué)院,福建漳州363000；2.福建省粒計(jì)算及其應(yīng)用重點(diǎn)實(shí)驗(yàn)室,福建漳州363000；3.福建省數(shù)據(jù)科學(xué)與統(tǒng)計(jì)重點(diǎn)實(shí)驗(yàn)室,福建漳州363000；4.數(shù)字福建氣象大數(shù)據(jù)研究所,福建漳州363000）

在傳統(tǒng)的統(tǒng)計(jì)過(guò)程控制（SPC）中，質(zhì)量特征的觀測(cè)值總是假設(shè)滿足相互獨(dú)立的要求，而當(dāng)質(zhì)量特征具有相依性時(shí)，傳統(tǒng)的指數(shù)加權(quán)移動(dòng)平均（EWMA）控制圖發(fā)揮了良好的效果.1959年，Roberts[1]首次提出EWMA控制圖，此后眾多專(zhuān)家學(xué)者對(duì)監(jiān)控過(guò)程均值變化[2]或監(jiān)控過(guò)程標(biāo)準(zhǔn)差變化[3]的EWMA 控制圖進(jìn)行研究.而在實(shí)際生產(chǎn)生活中，許多情形需要同時(shí)監(jiān)控過(guò)程均值和過(guò)程標(biāo)準(zhǔn)差的變化，因此一些學(xué)者提出了可以同時(shí)監(jiān)控均值和標(biāo)準(zhǔn)差的EWMA控制圖[4].

然而，隨著測(cè)量技術(shù)的發(fā)展，越來(lái)越多的數(shù)據(jù)以相依區(qū)間值的形式來(lái)描述或記錄.例如氣象部門(mén)需要每小時(shí)檢測(cè)一次空氣中顆粒物的濃度，顆粒物殘留在空氣中的時(shí)間較長(zhǎng)，這就導(dǎo)致了檢測(cè)值是相互依賴的，即后面所檢測(cè)的濃度還是會(huì)受到前面顆粒物滯留的影響.當(dāng)收集到的原始數(shù)據(jù)是相依區(qū)間值數(shù)據(jù)時(shí)，就不適合應(yīng)用傳統(tǒng)的EWMA控制圖來(lái)監(jiān)測(cè)質(zhì)量特征的變化.因此，如何設(shè)計(jì)基于相依區(qū)間值數(shù)據(jù)質(zhì)量特征的控制圖是一個(gè)重要而有趣的課題.

當(dāng)數(shù)據(jù)是區(qū)間值數(shù)據(jù)時(shí)，傳統(tǒng)的EWMA控制圖和EWMA MSE 控制圖通常是用對(duì)區(qū)間取平均值后的單一值來(lái)構(gòu)建控制圖，但是這樣的做法有可能會(huì)丟失一些重要的信息.在質(zhì)量特征具有相依性且是區(qū)間值的數(shù)據(jù)形態(tài)情況下，如何在保留相依區(qū)間值數(shù)據(jù)集的原始信息的同時(shí)，還能及時(shí)檢測(cè)質(zhì)量特征的變化，是質(zhì)量管理的重要課題.Yang[5]評(píng)估了基于相依區(qū)間值數(shù)據(jù)的模糊相對(duì)加權(quán)移動(dòng)平均（FRWMA）控制圖的監(jiān)控性能，發(fā)現(xiàn)該圖的監(jiān)控性能優(yōu)于傳統(tǒng)的移動(dòng)平均控制圖，可以更快地檢測(cè)到過(guò)程均值和方差的變化.

針對(duì)相依區(qū)間值數(shù)據(jù)，模糊區(qū)間均方誤差（FIMSE）作為均值偏移目標(biāo)值和過(guò)程均值方差變異的綜合度量，可以快速地檢測(cè)到區(qū)間值數(shù)據(jù)的變異，而指數(shù)加權(quán)移動(dòng)平均（EWMA）方法對(duì)全部歷史數(shù)據(jù)分配了逐步減少的權(quán)重，可以快速地檢測(cè)到相依數(shù)據(jù)小偏移的變異.首次將指數(shù)加權(quán)移動(dòng)平均（EWMA）方法與模糊區(qū)間均方誤差（FIMSE）統(tǒng)計(jì)量相結(jié)合，提出了基于相依區(qū)間值數(shù)據(jù)的指數(shù)加權(quán)移動(dòng)平均的模糊區(qū)間均方誤差（Exponentially Weighted Moving Average Fuzzy Interval Mean Square Error；EWMA FIMSE）控制圖來(lái)監(jiān)控過(guò)程質(zhì)量特征的變化，這樣可以在最大程度上保留區(qū)間數(shù)據(jù)所蘊(yùn)含的信息并提高監(jiān)控性能.

1 EWMA FIMSE控制圖的設(shè)計(jì)

1959年，Roberts[1]首次提出EWMA控制圖，其中EWMA統(tǒng)計(jì)量形式如下：

其中，μ為總體均值，Xˉi為樣本均值，λ為平滑系數(shù)，0<λ≤1.與傳統(tǒng)的控制圖相比，它的優(yōu)勢(shì)在于賦予歷史數(shù)據(jù)以一定的權(quán)重，距離單位時(shí)間點(diǎn)越近的數(shù)據(jù)權(quán)重值越大，并且隨著距離單位時(shí)間點(diǎn)距離的增加，這種權(quán)重以指數(shù)的形式遞減.EWMA控制圖可以很好地應(yīng)用于相依單一值數(shù)據(jù)均值的監(jiān)測(cè).

2004年，Chen[4]等提出了一種新的EWMA控制圖，其統(tǒng)計(jì)量形式如下：

其中T?i=(n/σ2)Ti，Ti=(Xˉi?μ)2+(n?1)S2i/n,i=1,2,…,λ為平滑系數(shù)，0<λ≤1，與EWMA控制圖相比，它的優(yōu)勢(shì)在于可以同時(shí)監(jiān)測(cè)相依單一值數(shù)據(jù)的均值和方差.

除均值和方差以外，在質(zhì)量控制領(lǐng)域中，接近目標(biāo)值也是一個(gè)值得關(guān)注的問(wèn)題.1998年，Spiring 等[6]首次提出的MSE控制圖可以同時(shí)監(jiān)測(cè)均值和方差以及與目標(biāo)值的接近程度，其中MSE統(tǒng)計(jì)量形式如下：

其中T為目標(biāo)值，Xˉj為樣本均值，S2j為樣本方差.

借鑒MSE 統(tǒng)計(jì)量的構(gòu)造，將Chen 等提出統(tǒng)計(jì)量中的μ用目標(biāo)值T來(lái)替換，構(gòu)造EWMA MSE 統(tǒng)計(jì)量來(lái)設(shè)計(jì)EWMA MSE 控制圖，它的好處是可以同時(shí)監(jiān)測(cè)相依單一值數(shù)據(jù)的均值和方差以及與目標(biāo)值的接近程度.

1.1 EWMA MSE控制圖

設(shè)xi1,xi2,…,xin,i=1,2,…,是來(lái)自均值為μ，方差為σ2的正態(tài)分布N(μ,σ2)的隨機(jī)變量，其均方誤差統(tǒng)計(jì)量可表示為：

其中T為目標(biāo)值，令MSE*i=(n?1)MSEi/σ2，則MSE*i統(tǒng)計(jì)量服從自由度為n，參數(shù)為γ=n[(μ?T)/σ]2的非中心卡方分布，即MSE*i～，令μ?T=δ3σ，δ3>0，δ3為過(guò)程均值偏離目標(biāo)值的幅度，當(dāng)過(guò)程受控時(shí)，過(guò)程均值沒(méi)有偏離目標(biāo)值，即δ3=0時(shí)，MSE*i的均值和方差分別為：

構(gòu)造EWMA MSE統(tǒng)計(jì)量Qi如下：

其中Q0=n+n[(μ-T)/σ]2為初始值，λ為平滑系數(shù)，0<λ≤1，當(dāng)過(guò)程受控時(shí)，Qi的均值和方差分別為：

因?yàn)榻y(tǒng)計(jì)量Qi大于等于0，所以控制下限（LCL）為0，控制上限為：

其中L為控制界限系數(shù).

1.2 EWMA FIMSE控制圖

注意到，當(dāng)數(shù)據(jù)為區(qū)間值數(shù)據(jù)時(shí)，上述的EWMA MSE 控制圖只能將區(qū)間的最大值最小值取平均值，將區(qū)間值轉(zhuǎn)化為單一值來(lái)計(jì)算，這樣的方法有可能會(huì)損失大量重要的信息.為了在最大程度上保留區(qū)間值數(shù)據(jù)所蘊(yùn)含的信息，首次提出EWMA FIMSE控制圖，它可以同時(shí)監(jiān)測(cè)相依區(qū)間值數(shù)據(jù)的均值和方差以及與目標(biāo)值的接近程度.

假設(shè)有一組來(lái)自均值為μ，方差為σ2的正態(tài)分布N(μ,σ2)的相依區(qū)間值隨機(jī)變量序列(xi1,yi1),(xi2,yi2),…,(xin,yin),i=1,2,…,且滿足關(guān)系式xin

令FIMSE1i=min(MSE1i,MSE2i)，F(xiàn)IMSE2i=max(MSE1i,MSE2i)，定義模糊區(qū)間均方誤差統(tǒng)計(jì)量FIMSEi為：

令FIMSE*i=(n?1)FIMSEi/σ2，F(xiàn)IMSE*i服從參數(shù)為γ=n[(μ?T)/σ]2=nδ23，自由度為n的非中心卡方分布.

當(dāng)過(guò)程受控時(shí)，F(xiàn)IMSE*i的均值和方差分別為：

定義EWMA FIMSE統(tǒng)計(jì)量FQi如下：

則

其中FQ0=[n+n[(μ?T)/σ]2，n+n[(μ?T)/σ]2]為初始值，λ為平滑系數(shù)，0<λ≤1，當(dāng)過(guò)程受控時(shí)，F(xiàn)Qi的均值和方差分別為：

Wu等[7]介紹關(guān)于模糊區(qū)間變量的假設(shè)檢驗(yàn)方法，本文借鑒該方法構(gòu)建FQi的假設(shè)檢驗(yàn)問(wèn)題，設(shè)

因?yàn)榻y(tǒng)計(jì)量FQi大于等于0，且越趨于0 代表與目標(biāo)值的差異越小，所以設(shè)控制下限（LCL）為0，在α的顯著水平下，F(xiàn)Qi落入拒絕域的概率函數(shù)為：

當(dāng)過(guò)程受控時(shí)，控制上限為：

其中L為控制界限系數(shù).

2 模擬分析

本節(jié)將使用蒙特卡羅模擬法來(lái)生成模糊區(qū)間值的時(shí)間序列數(shù)據(jù)，這些模糊區(qū)間值數(shù)據(jù)用于評(píng)估EW‐MA MSE控制圖的平均運(yùn)行長(zhǎng)度（ARL）值，ARL指第一個(gè)質(zhì)量特征落在控制界限之外時(shí)所期望的樣本子組的數(shù)量，受控狀態(tài)下用ARL0表示，失控狀態(tài)下用ARL1表示，ARL1越小說(shuō)明控制圖的統(tǒng)計(jì)性能越好.并通過(guò)均方誤差值將EWMA FIMSE 控制圖的性能與傳統(tǒng)的EWMA MSE 控制圖、EWMA控制圖進(jìn)行比較.

計(jì)算ARL[8]的步驟如下：

1）生成一組服從均值為μ，標(biāo)準(zhǔn)差為σ的正態(tài)分布的模糊區(qū)間樣本觀測(cè)值.

2）利用生成的模糊區(qū)間樣本觀測(cè)值分別采用式(1)、式(3)、式(17)計(jì)算統(tǒng)計(jì)量EWMA，EWMA MSE，EWMA FIMSE.

3）利用式(21)、式(22)計(jì)算當(dāng)觀測(cè)統(tǒng)計(jì)值大于控制下限（LCL）且小于控制上限（UCL）時(shí)的ARL0，當(dāng)ARL0的值近似370（固定α=0.0027，即ARL0=1/α=370）時(shí)，得到控制界限系數(shù)L的值.

4）生成一組服從均值為μ+δ1σ，δ1>0，標(biāo)準(zhǔn)差為δ2σ，δ2>1 的正態(tài)分布的模糊區(qū)間樣本觀測(cè)值，計(jì)算失控狀態(tài)下，各觀測(cè)統(tǒng)計(jì)值小于控制下限（LCL）和大于控制上限（UCL）時(shí)的ARL1.

表1列出了n=5，δ3=0.0 和1.0，λ=0.1，0.2 和0.3，EWMA控制圖，EWMA MSE 控制圖和EWMA FIMSE 控制圖的ARL0近似等于370 時(shí)的控制界限系數(shù)L 的值.隨著δ3的增加，EWMA FIMSE 控制圖和EWMA MSE 控制圖對(duì)應(yīng)的L 值增加，EWMA控制圖對(duì)應(yīng)的L 值不變.隨著λ的增加，EWMA FIMSE 控制圖對(duì)應(yīng)的L值減少，EWMA控制圖對(duì)應(yīng)的L值增加.說(shuō)明EWMA FIMSE控制圖和EWMA MSE控制圖的控制界限系數(shù)L 的值的選取與參數(shù)δ3和λ的選取有關(guān)，EWMA控制圖的控制界限系數(shù)L 的值的選取只與參數(shù)λ的選取有關(guān).

表1 基于EWMA 、EWMA MSE、EWMA FIMSE三個(gè)控制圖的ARL0Tab.1 ARL0 based on EWMA ,EWMA MSE and EWMA FIMSE control charts

表1 基于EWMA 、EWMA MSE、EWMA FIMSE三個(gè)控制圖的ARL0Tab.1 ARL0 based on EWMA ,EWMA MSE and EWMA FIMSE control charts

δ3 λ EWMA images/BZ_7_401_2895_442_2945.pngL 2.702 2.860 2.925 2.702 2.860 2.925 ARL0 369.9 370.1 370.4 369.9 370.1 370.4 EWMA MSE L 2.731 3.168 3.445 11.482 9.997 9.407 ARL0 370.1 370.8 370.4 370.1 369.9 370.0 0.0 1.0 0.1 0.2 0.3 0.1 0.2 0.3 EWMA FIMSE L 5.130 4.849 4.779 19.908 16.040 14.345 ARL0 370.3 370.8 369.8 370.8 370.5 369.8

表2列出了λ=0.1 時(shí)EWMA控制圖、EWMA MSE 控制圖和EWMA FIMSE 控制圖的ARL1.顯然，EWMA FIMSE 控制圖的ARL1都小于EWMA MSE 控制圖和EWMA控制圖的ARL1，例如：當(dāng)δ3=0.0，δ1=0.5，δ2=1.1 時(shí)，EWMA FIMSE 控制圖的ARL1=11.83，EWMA MSE 控制圖的ARL1=18.03，EWMA控制圖的ARL1=26.28，這表明EWMA FIMSE 控制圖比EWMA MSE 控制圖和EWMA控制圖更敏感，特別是當(dāng)δ1和δ2取較小值時(shí)，EWMA FIMSE控制圖的監(jiān)測(cè)效果更好.

表2 基于EWMA 、EWMA MSE、EWMA FIMSE三個(gè)控制圖的ARL1Tab.2 ARL1 based on EWMA ,EWMA MSE and EWMA FIMSE control charts

表2 基于EWMA 、EWMA MSE、EWMA FIMSE三個(gè)控制圖的ARL1Tab.2 ARL1 based on EWMA ,EWMA MSE and EWMA FIMSE control charts

δ3 δ1 δ2 0.1 0.0 0.5 1.0 0.1 1.0 0.5 1.0 1.1 1.3 1.5 1.1 1.3 1.5 1.1 1.3 1.5 1.1 1.3 1.5 1.1 1.3 1.5 1.1 1.3 1.5 EWMA images/BZ_7_401_2895_442_2945.png154.25 76.36 46.12 26.28 23.17 20.56 9.683 9.550 9.396 154.25 76.36 46.12 26.28 23.17 20.56 9.683 9.550 9.396 EWMA MSE 44.81 11.76 6.610 18.03 8.514 5.546 6.747 4.893 3.887 40.42 16.59 10.22 10.17 7.598 6.079 4.795 4.154 3.707 EWMA FIMSE 34.51 10.15 6.194 11.83 7.006 5.110 5.168 4.280 3.680 31.97 14.68 9.552 9.705 6.526 5.916 4.446 3.816 3.328

3 實(shí)證分析

作為應(yīng)用，搜集了漳州市2019年3月每天5 個(gè)時(shí)間段內(nèi)一氧化碳（CO）的最高濃度和最低濃度的共155 個(gè)區(qū)間數(shù)據(jù)（https://quotsoft.net/air/），將其按日期分成31 組進(jìn)行分析，試圖檢測(cè)出這個(gè)時(shí)間段的CO濃度異常的時(shí)間點(diǎn).為此，分別構(gòu)建EWMA控制圖，EWMA MSE 控制圖和EWMA FIMSE 控制圖.利用S-W檢驗(yàn)數(shù)據(jù)是否服從正態(tài)分布，算得p=0.111 2>0.05，說(shuō)明數(shù)據(jù)服從正態(tài)分布.

首先，為了估計(jì)過(guò)程受控時(shí)的樣本均值和標(biāo)準(zhǔn)偏差，先計(jì)算搜集的31 組區(qū)間值樣本數(shù)據(jù)每個(gè)時(shí)間段CO 的最高濃度和最低濃度的均值，繪制傳統(tǒng)的Xˉ控制圖、S控制圖，剔除失控的樣本組.再用受控的樣本組計(jì)算出過(guò)程受控時(shí)的平均值和標(biāo)準(zhǔn)偏差，計(jì)算全部的樣本組的平均值和標(biāo)準(zhǔn)偏差，從而估計(jì)出偏移值然后分別計(jì)算EWMA、EWMA MSE、EWMA FIMSE 統(tǒng)計(jì)量的值，并分別繪制EWMA控制圖、EWMA MSE控制圖、EWMA FIMSE控制圖.

為了計(jì)算方便，假設(shè)目標(biāo)值T=0.672 1（T值可視實(shí)際需要而定），剔除Xˉ控制圖和S控制圖中失控的樣本組26，27，28，30，31，即圖1、圖2中上下界外的的樣本組.算得過(guò)程受控時(shí)的樣本均值和標(biāo)準(zhǔn)偏差分別為=0.790 3和=0.118 2，估計(jì)的偏移值分別為

圖1 控制圖Fig.1 control chart

圖2 S控制圖Fig.2 S control chart

如圖3-圖5所示，在EWMA控制圖、EWMA MSE控制圖和EWMA FIMSE控制圖中，最先被監(jiān)測(cè)到出界的樣本組分別為第22號(hào)、第23號(hào)和第26號(hào)，即分別在3月22日、3月23日和3月26日發(fā)現(xiàn)空氣中CO濃度的異常.這表明，使用EWMA FIMSE 控制圖來(lái)監(jiān)控漳州市CO 濃度這一“質(zhì)量特征”的效果要優(yōu)于EWMA MSE 控制圖和EWMA控制圖.因此，當(dāng)質(zhì)量特征滿足模糊相依區(qū)間值數(shù)據(jù)時(shí)，EWMA FIMSE控制圖可以很好地應(yīng)用于工業(yè)、經(jīng)濟(jì)、環(huán)境和管理領(lǐng)域的質(zhì)量控制.

圖3 EWMA 控制圖Fig.3 EWMA control chart

圖4 EWMA MSE控制圖Fig.4 EWMA MSE control chart

圖5 EWMA FIMSE控制圖Fig.5 EWMA FIMSE control chart

據(jù)調(diào)查，2019年3月下旬正是春節(jié)假期結(jié)束之時(shí)，大多數(shù)企事業(yè)單位、工廠等在這段時(shí)間復(fù)工、復(fù)產(chǎn)，這些因素導(dǎo)致車(chē)流量急劇增加，汽車(chē)尾氣的排放導(dǎo)致了CO 濃度的升高，數(shù)據(jù)分析結(jié)果從CO 的濃度異常變化角度很大程度地反映出這一社會(huì)現(xiàn)象，因此本文提出的方法有助于政府有關(guān)部門(mén)采取相關(guān)措施，實(shí)時(shí)監(jiān)測(cè)這樣的社會(huì)活動(dòng)變化，做好預(yù)案.

4 結(jié)論

針對(duì)具有相依性的區(qū)間值數(shù)據(jù)，提出了一種基于相依區(qū)間值數(shù)據(jù)的指數(shù)加權(quán)移動(dòng)平均模糊區(qū)間均方誤差（EWMA FIMSE）控制圖.通過(guò)模擬和實(shí)證分析表明，所提出的EWMA FIMSE 控制圖的監(jiān)控性能優(yōu)于傳統(tǒng)的EWMA MSE 控制圖和EWMAXˉ控制圖，可以更快地監(jiān)測(cè)到過(guò)程均值和方差發(fā)生較小的變化時(shí)的失控現(xiàn)象.除了在氣象監(jiān)測(cè)上的應(yīng)用，也可以應(yīng)用在經(jīng)濟(jì)、工程等領(lǐng)域.