Minitab在統(tǒng)計(jì)學(xué)專(zhuān)業(yè)實(shí)踐課程中的案例拓展

束 慧， 徐宏光

(南京郵電大學(xué) 理學(xué)院，南京210023)

1 引 言

統(tǒng)計(jì)學(xué)專(zhuān)業(yè)實(shí)踐課程是培養(yǎng)學(xué)生從“知識(shí)”到“應(yīng)用”的重要環(huán)節(jié)，在整個(gè)統(tǒng)計(jì)學(xué)專(zhuān)業(yè)培養(yǎng)方案中具有不可替代的地位.近年來(lái)，一些新的、側(cè)重理論學(xué)習(xí)與實(shí)踐操作相融合的教學(xué)理念在國(guó)內(nèi)高校中被廣泛推廣[1]，越來(lái)越多的高校統(tǒng)計(jì)學(xué)、數(shù)學(xué)方向的教師開(kāi)始設(shè)計(jì)實(shí)踐教學(xué)的新模型、新方法[2-3].作為一款商業(yè)統(tǒng)計(jì)軟件，Minitab在統(tǒng)計(jì)學(xué)專(zhuān)業(yè)的實(shí)踐教學(xué)中具有不可替代的地位，尤其是涉及到6σ管理等質(zhì)量控制課程的實(shí)踐時(shí)，其地位更加凸顯，圍繞Minitab設(shè)計(jì)的實(shí)踐教學(xué)案例陸續(xù)涌現(xiàn)[4].

Minitab是由美國(guó)賓夕法尼亞州立大學(xué)開(kāi)發(fā)的一款商用統(tǒng)計(jì)軟件，其特色功能為質(zhì)量控制模塊，被廣泛運(yùn)用于工業(yè)企業(yè)的產(chǎn)品質(zhì)量控制過(guò)程[5].以Minitab17版本為例，該軟件的質(zhì)量控制模塊由控制圖、質(zhì)量控制兩個(gè)子模塊組成.在控制圖模塊中，整合了子組變量、單值變量、屬性、時(shí)間加權(quán)、多屬性變量以及稀有事件控制圖等各類(lèi)控制圖的繪制功能.在運(yùn)用Minitab繪制各類(lèi)控制圖時(shí)，控制限的選擇是開(kāi)放給用戶(hù)的，提供了“標(biāo)準(zhǔn)差的倍數(shù)”或者指定“控制上、下限”的權(quán)限[5].這樣的功能設(shè)計(jì)，具有操作上的簡(jiǎn)便性，但同時(shí)也增加了用戶(hù)在繪制控制圖時(shí)的隨機(jī)性，尤其是繪制一些檢測(cè)生產(chǎn)過(guò)程中的中小變異的控制圖，如指數(shù)加權(quán)移動(dòng)平均(Exponentially Weighted Moving Average，EWMA)[6-9]，累計(jì)和CUSUM (Cumulative Sum)[10-11]等控制圖時(shí)，控制限設(shè)置存在一定的優(yōu)化空間.

本文從統(tǒng)計(jì)學(xué)專(zhuān)業(yè)實(shí)踐教學(xué)案例拓展的角度出發(fā)，對(duì)控制圖的重要優(yōu)化指標(biāo)進(jìn)行以Minitab為基礎(chǔ)的算法設(shè)計(jì)和計(jì)算，并根據(jù)運(yùn)算結(jié)果對(duì)控制圖的關(guān)鍵參數(shù)進(jìn)行優(yōu)化選擇.本文的研究既可以作為統(tǒng)計(jì)學(xué)專(zhuān)業(yè)實(shí)踐課程的案例拓展，也可以作為Minitab統(tǒng)計(jì)軟件質(zhì)量控制模塊功能優(yōu)化研究，可以為統(tǒng)計(jì)學(xué)專(zhuān)業(yè)實(shí)踐課程提供實(shí)踐內(nèi)容和實(shí)踐工具兩個(gè)方面的參考.

本文第二部分選取典型的控制圖，簡(jiǎn)要介紹其優(yōu)化指標(biāo)的計(jì)算方法，以及基于指標(biāo)的控制圖優(yōu)化原理；以理論方法為基礎(chǔ)，設(shè)計(jì)控制圖優(yōu)化指標(biāo)在Minitab上的實(shí)現(xiàn)算法，并提供一個(gè)案例計(jì)算結(jié)果；第三部分是本文的分析總結(jié)，基于文章的結(jié)果提供一些Minitab軟件在高校統(tǒng)計(jì)學(xué)專(zhuān)業(yè)實(shí)踐教學(xué)環(huán)節(jié)中的應(yīng)用拓展建議.

2 實(shí)踐案例構(gòu)建——基于Minitab 的控制圖優(yōu)化

2.1 控制圖優(yōu)化理論模型

(1)

在進(jìn)行質(zhì)量檢驗(yàn)過(guò)程中，如果抽樣的時(shí)間間隔及樣本容量n固定不變，則稱(chēng)相應(yīng)的控制圖為固定抽樣率(Fixed Sampling Rate，F(xiàn)SR)控制圖.對(duì)于一個(gè)FSR控制圖，定義從檢測(cè)開(kāi)始到它發(fā)出生產(chǎn)變異的報(bào)警標(biāo)識(shí)為止，抽取的平均樣本組數(shù)為平均運(yùn)行長(zhǎng)度(Average Run Length， ARL).當(dāng)過(guò)程可控時(shí)，質(zhì)量控制圖的報(bào)警就屬于誤報(bào)，可控的ARL越大越好；當(dāng)過(guò)程失控時(shí)，質(zhì)量控制圖應(yīng)盡早報(bào)警，失控的ARL越小越好.

2.2 控制圖優(yōu)化指標(biāo)ARL的計(jì)算方法

控制圖ARL的計(jì)算方法多樣，研究較豐富，運(yùn)用有限狀態(tài)Markov鏈的平均首達(dá)時(shí)間計(jì)算ARL是比較常見(jiàn)的[6，7，12，13].這部分內(nèi)容結(jié)合統(tǒng)計(jì)學(xué)專(zhuān)業(yè)實(shí)踐環(huán)節(jié)教學(xué)要求，可以設(shè)置為對(duì)學(xué)生的統(tǒng)計(jì)學(xué)理論基礎(chǔ)的考察項(xiàng)目.

在上述前提下，控制上限可以表示為

UCL=L·σs=h，

(2)

相應(yīng)的控制下限

LCL=-h，

(3)

其中，L為控制限參數(shù)，是優(yōu)化控制圖時(shí)的重要選擇參數(shù).

ARL計(jì)算的基本依據(jù)是有限狀態(tài)Markov鏈從非吸收態(tài)首次轉(zhuǎn)移到吸收態(tài)的平均時(shí)間，因此，需要對(duì)控制圖界限進(jìn)行Markov鏈狀態(tài)空間轉(zhuǎn)換.將控制圖的上下控制限之間的區(qū)間分成k(奇數(shù))個(gè)子區(qū)間，小區(qū)間的寬度d

(4)

Markov鏈的狀態(tài)空間可以用{Ij}j=1，2，…，k表示，且設(shè)Ij為小區(qū)間的中心點(diǎn)

(5)

當(dāng)EWMA檢驗(yàn)統(tǒng)計(jì)量St落在第j個(gè)區(qū)間內(nèi)時(shí)，即

就定義統(tǒng)計(jì)量St處于狀態(tài)Ij.EWMA控制圖是雙邊的，將大于上控制限h的區(qū)域與小于下控制限-h的區(qū)域這兩個(gè)吸收域合并為一個(gè).在上述設(shè)置下，該Markov鏈的一步轉(zhuǎn)移概率pij為

pij=P{St+1∈Ij|St∈Ii}

phj=0，j=1，2，…，k;

phh=1.

(6)

檢測(cè)統(tǒng)計(jì)量St在m時(shí)刻的一步轉(zhuǎn)移概率陣為

(7)

矩陣P的最后一列表示從狀態(tài)Ii出發(fā)一步轉(zhuǎn)移到吸收狀態(tài)Ih的概率，因此可記

以fih表示從狀態(tài)Ii出發(fā)第一次轉(zhuǎn)移到吸收狀態(tài)Ih所需要的步數(shù)(首達(dá)時(shí)間)，fh=(f1h，…，f(h-1)h)T，則當(dāng)檢測(cè)統(tǒng)計(jì)量的初始態(tài)為Ii時(shí)的ARL值即為E(fh)的第i個(gè)分量，特別地，初始狀態(tài)是從中心點(diǎn)出發(fā)，則為中間位置的分量.根據(jù)有限狀態(tài)馬氏鏈的平均首達(dá)時(shí)間計(jì)算方法，有

ARL=E(fh)=(I-R)-11，

(8)

2.3 基于Minitab的ARL計(jì)算及控制圖優(yōu)化

根據(jù)上述理論模型，控制圖優(yōu)化的基本思路是：給定變異系數(shù)δ，狀態(tài)空間劃分k，優(yōu)化參數(shù)：平滑系數(shù)λ、控制限L，使得控制圖在受控態(tài)下(變異系數(shù)δ=0)的ARL較大，而失控狀態(tài)(變異系數(shù)δ≠0)下的ARL較小，同時(shí)分析狀態(tài)空間劃分對(duì)相應(yīng)的ARL優(yōu)化的影響.

借助Minitab的數(shù)據(jù)整理、數(shù)據(jù)(矩陣)計(jì)算功能，以Minitab基于列變量的運(yùn)算規(guī)則為語(yǔ)句設(shè)計(jì)特點(diǎn)，在Minitab宏命令“GMACRO”功能框架下，定義ARL計(jì)算過(guò)程腳本，設(shè)計(jì)基于Minitab的控制圖優(yōu)化指標(biāo)ARL計(jì)算及反饋優(yōu)化機(jī)制(圖1).

圖1 基于Minitab的EWMA控制圖優(yōu)化指標(biāo)計(jì)算及反饋

Step1 設(shè)置不同的控制限L，偏移系數(shù)δ，狀態(tài)空間k，平滑系數(shù)λ值，計(jì)算相應(yīng)的LCL，d;

Step2 在計(jì)算轉(zhuǎn)移概率矩陣時(shí)，借助Minitab的“Mesh”命令，可以產(chǎn)生計(jì)算轉(zhuǎn)移概率矩陣的網(wǎng)格化指標(biāo)(i，j)，替代一般編程語(yǔ)言中的循環(huán)結(jié)構(gòu);

Step3 按列存儲(chǔ)(i，j，LCL，L，λ，d)，其中的(LCL，L，λ，d)利用Minitab的產(chǎn)生模板數(shù)據(jù)，和對(duì)應(yīng)的網(wǎng)格數(shù)據(jù)對(duì)一致， 以(i，j，LCL，L，λ，d)值為輸入，計(jì)算相應(yīng)的轉(zhuǎn)移概率值pij;

Step4 利用Minitab的拆分列“Unstack”、列轉(zhuǎn)置“Transpose”，以及數(shù)據(jù)復(fù)制“Copy”功能，將列變量存儲(chǔ)為轉(zhuǎn)移概率矩陣R；

Step5 計(jì)算ARL值，利用Minitab的“Diagonal”產(chǎn)生單位矩陣，“Subtract”矩陣運(yùn)算，“Invert”矩陣求逆，“Multiply”矩陣乘積等宏命令；

Step6 (反饋機(jī)制)記錄不同的偏移系數(shù)δ、狀態(tài)空間劃分k，控制限L，平滑系數(shù)λ下的ARL計(jì)算結(jié)果，分組擬合ARL與k，λ值的關(guān)系曲線(xiàn)，從中篩選出符合要求的優(yōu)化參數(shù)設(shè)置.具體，根據(jù)EWMA控制圖檢測(cè)中小偏移為主的特征，設(shè)置偏移系數(shù)δ=0，δ=1，分別擬合出ARL與(k，λ)的等值線(xiàn)圖(圖2)； ARL與k的關(guān)系曲線(xiàn)圖(圖3)；ARL與λ的關(guān)系曲線(xiàn)圖(圖4).以上圖形都是用Minitab作圖功能繪制.

(k，λ)組合對(duì)ARL的影響是綜合性的：δ=0時(shí)，不同控制限制L的設(shè)置下，ARL取值隨著λ取值的變化，始終呈現(xiàn)“梯度”特征，且在平滑系數(shù)λ取值較低的階段，ARL均較長(zhǎng).這樣的結(jié)果說(shuō)明了過(guò)程受控時(shí)，若平滑系數(shù)λ取值較低(0.1附近)，EWMA控制圖可以減少“誤報(bào)”(圖2).

δ=1時(shí)， ARL取值隨著λ取值的變化雖然仍呈現(xiàn)“梯度”特征，但在不同控制限制L的設(shè)置下，ARL與平滑系數(shù)λ取值的變化并沒(méi)有呈現(xiàn)相似的變化趨勢(shì).L=3，平滑系數(shù)λ取值較低區(qū)間[0.1，0.4]，ARL較短；L=2，平滑系數(shù)λ取值落在[0.4，0.6]區(qū)間，ARL較短.這樣的結(jié)果說(shuō)明了過(guò)程非受控時(shí)， EWMA控制圖可以盡早“報(bào)警”，但需要將控制限和平滑系數(shù)進(jìn)行動(dòng)態(tài)的調(diào)整(圖2).

圖2 控制圖受控、失控時(shí)，ARL與(k，λ)的等值線(xiàn)圖

固定(k，λ)中的一個(gè)，ARL與另一參數(shù)間的關(guān)系則更加特征鮮明：

δ=0時(shí)，不同控制限制L的設(shè)置下，ARL受狀態(tài)空間k的影響較小，僅僅在平滑系數(shù)λ取值較低(0.1)時(shí)，出現(xiàn)一定的影響作用(圖3(a)); δ=1時(shí)，不同控制限制L的設(shè)置下，ARL受狀態(tài)空間劃分k的影響均較小(圖3(b)).

(a) 受控時(shí)，ARL與k的關(guān)系擬合圖 (b) 非受控時(shí)，ARL與k的關(guān)系擬合圖圖3 ARL與k的關(guān)系擬合圖

δ=0時(shí)，不同控制限制L的設(shè)置下，ARL受平滑系數(shù)λ設(shè)置的影響較大，和等值線(xiàn)圖顯示的結(jié)果一致，平滑系數(shù)λ取值較低(0.1)時(shí)，ARL均較長(zhǎng)(圖4(a)).

(a) 受控時(shí)，ARL與λ的關(guān)系擬合圖 (b) 非受控時(shí)，ARL與λ的關(guān)系擬合圖圖4 ARL與λ的關(guān)系擬合圖

δ=1時(shí)，不同控制限制L的設(shè)置下，ARL隨著平滑系數(shù)λ的取值不同，呈現(xiàn)差異化的變化趨勢(shì)：控制限較大，ARL隨著λ的取值變大而變長(zhǎng)；控制限中等，則ARL隨著λ的取值變大而先變短后又變長(zhǎng)；控制限較小時(shí)，ARL隨著λ的取值變大而變短(圖4(b)).

綜合以上仿真分析結(jié)果，遵循控制圖優(yōu)化的基本原理，受控時(shí)，EWMA控制圖應(yīng)選取較小的平滑系數(shù)值；而當(dāng)過(guò)程非受控時(shí)，應(yīng)結(jié)合控制限L的寬度，相應(yīng)的平滑系數(shù)值也要作出變化：L寬，則λ小；L中等，則λ中等；L窄，則λ大.這樣的取值關(guān)系，可以使得受控時(shí)的ARL較長(zhǎng)，非受控時(shí)ARL較短，實(shí)現(xiàn)控制圖優(yōu)化的目標(biāo).

3 結(jié) 論

本文以統(tǒng)計(jì)軟件Minitab為基礎(chǔ)，設(shè)計(jì)了一類(lèi)統(tǒng)計(jì)學(xué)專(zhuān)業(yè)實(shí)踐教學(xué)中的拓展案例.該案例以質(zhì)量控制領(lǐng)域的控制圖優(yōu)化理論為基礎(chǔ)，結(jié)合Minitab的數(shù)據(jù)處理特征、功能模塊以及宏命令語(yǔ)言，實(shí)現(xiàn)了控制圖的仿真優(yōu)化.該案例完整、高效的實(shí)現(xiàn)過(guò)程說(shuō)明，功能性統(tǒng)計(jì)軟件在高校統(tǒng)計(jì)學(xué)專(zhuān)業(yè)實(shí)踐教學(xué)中具有非常大的潛力，這類(lèi)軟件的教學(xué)推廣，是對(duì)廣泛流行的SAS、SPSS、R以及Stata等主流統(tǒng)計(jì)軟件的有益補(bǔ)充.

致謝作者非常感謝相關(guān)文獻(xiàn)對(duì)本文的啟發(fā)以及審稿專(zhuān)家提出的寶貴意見(jiàn).