寧方華，騫文成，屠震元，陳智峰

（1.浙江理工大學(xué) 機(jī)械工程學(xué)院，浙江 杭州 310018；2.杭州?？低晹?shù)字技術(shù)股份有限公司，浙江 杭州 310051）

0 引言

在生產(chǎn)加工過程中，為了保證產(chǎn)品質(zhì)量需要同時監(jiān)控每個方面的質(zhì)量指標(biāo)?，F(xiàn)有研究人員將傳統(tǒng)單變量質(zhì)量控制圖技術(shù)拓展到多元質(zhì)量控制圖技術(shù)，通過一個控制圖監(jiān)控生產(chǎn)過程中的多元質(zhì)量指標(biāo)，以提升監(jiān)控效率［1］。

近年來，HotellingT2 控制圖技術(shù)一經(jīng)提出便得到了廣泛應(yīng)用［1-9］。然而，隨著客戶需求個性化和品種多樣化，各類產(chǎn)品生產(chǎn)過程開始步入小批量生產(chǎn)模式，傳統(tǒng)HotellingT2 控制圖技術(shù)無法直接適用于小批量生，原因為該技術(shù)通過樣本均值μˉ和樣本協(xié)方差矩陣S的計算結(jié)果來判斷生產(chǎn)質(zhì)量，但小批量生產(chǎn)過程中同批次的生產(chǎn)數(shù)量一般不大，采集樣本數(shù)據(jù)較少，μˉ、S易受異常值影響，因此計算結(jié)果不夠準(zhǔn)確，甚至發(fā)生錯誤，導(dǎo)致傳統(tǒng)HotellingT2 控制圖技術(shù)在小批量生產(chǎn)過程中的監(jiān)控效果與實際不符。

為此，本文采用貝葉斯理論［10-11］改進(jìn)傳統(tǒng)HotellingT2控制圖，通過有效估計小批量生產(chǎn)過程在質(zhì)量數(shù)據(jù)樣本的μˉ、S，以提升其對異常值的識別和處理能力，給出了新的基于貝葉斯理論的HotellingT2 工序質(zhì)量控制圖方法。同時，結(jié)合生產(chǎn)實例的比較研究顯示，在小批量生產(chǎn)過程中，本文所提方法相較于傳統(tǒng)HotellingT2 控制圖方法的計算結(jié)果更合理，具有更好的監(jiān)控效果。

1 傳統(tǒng)HotellingT2質(zhì)量控制圖

1.1 基本原理

當(dāng)某個工序需要控制p個質(zhì)量指標(biāo)，且多元質(zhì)量指標(biāo)的隨機(jī)公差服從正態(tài)分布［12］，質(zhì)量指標(biāo)記為X=(X1，X2，…，Xp)，將當(dāng)前加工該工序的多個工件的多元質(zhì)量指標(biāo)測量值作為樣本信息，平均值向量=(μ1，μ2，…，μp)表示加工過程多元質(zhì)量指標(biāo)標(biāo)準(zhǔn)值，方差向量表示多元質(zhì)量指標(biāo)誤差的收斂程度，通過協(xié)方差求解公式進(jìn)行求解：

式中：E(XiXj)表示第i個質(zhì)量指標(biāo)和第j個質(zhì)量指標(biāo)乘積的期望值；E()Xi、E(Xj)分別表示第i個質(zhì)量指標(biāo)和第j個質(zhì)量指標(biāo)期望值。

p個質(zhì)量指標(biāo)的協(xié)方差矩陣S為：

式中：當(dāng)i≠j時，Sij表示第i個質(zhì)量指標(biāo)和第j個質(zhì)量指標(biāo)的協(xié)方差數(shù)值且Sij=Sji；當(dāng)i=j時，Sii為第i個質(zhì)量指標(biāo)的方差，即Sii=。由式（1）、式（2）可反映此時加工工序多元質(zhì)量指標(biāo)波動情況的綜合指標(biāo)。傳統(tǒng)HotellingT2 控制圖的瞬時統(tǒng)計量為：

式中：xi表示本工序第i個工件的p個質(zhì)量指標(biāo)的測量值；S-1為協(xié)方差矩陣的逆矩陣。經(jīng)文獻(xiàn)［1］驗證T2統(tǒng)計量服從具有參數(shù)p和n-p的F分布，即：

式中：n表示加工工件的數(shù)量；p為每個工件需要控制的質(zhì)量指標(biāo)數(shù)量。

若取虛報概率為α，則HotellingT2 控制圖的控制上限UCL和控制下限LCL可表示為：

若監(jiān)控過程中T2統(tǒng)計量一直處于UCL、LCL 之間，證明加工過程處于受控狀態(tài)，若統(tǒng)計量超出范圍則證明出現(xiàn)異常，加工過程也不再處于受控狀態(tài)，控制圖將發(fā)出警報。根據(jù)這些特點，管理人員可進(jìn)一步分析生產(chǎn)不穩(wěn)定的原因，從而作出相應(yīng)調(diào)整或改進(jìn)生產(chǎn)過程的控制指標(biāo)，使加工過程回到受控狀態(tài)，提升生產(chǎn)合格率和穩(wěn)定性。

1.2 局限性

傳統(tǒng)HotellingT2 控制圖原理反映了加工工序的整體穩(wěn)定性，不僅監(jiān)控每個指標(biāo)在控制限制內(nèi)的分布或波動情況，還監(jiān)控每個指標(biāo)間的相關(guān)關(guān)系是否異常。然而，由式（3）可見(xi-)反映了變量取值的波動，由于(xi-)前后各乘一次，T2值為二次型，協(xié)方差為正定矩陣S，因此T2值為非負(fù)，控制圖只有上控制線，零為自然下限。

關(guān)系異常指每個質(zhì)量指標(biāo)無法保持穩(wěn)態(tài)下相互波動性的關(guān)系。由于在小批量加工過程中同批次加工樣本數(shù)據(jù)較少，因此μˉ、S對異常值非常敏感，容易導(dǎo)致計算結(jié)果會產(chǎn)生較大偏差，甚至可能導(dǎo)致監(jiān)控失效。

本文構(gòu)造了一個例子說明在小批量加工過程中，異常值對控制圖參數(shù)的影響。隨機(jī)生成一組雙變量數(shù)據(jù)，稱為數(shù)據(jù)組1（即不存在異常值的數(shù)據(jù)，如表1 所示），將6 號正常值60.5 和17.5 替換為異常值41.25 和20.75，稱為數(shù)據(jù)組2（即存在離群值的數(shù)據(jù)）。采用傳統(tǒng)方法計算數(shù)據(jù)組1、2的均值和協(xié)方差，如表2所示。

Table 1 Data of 10 bivariate variables表1 10個雙變量數(shù)據(jù)

Table 2 Comparison of mean and covariance matrix between data group 1 and data group 2表2 數(shù)據(jù)組1和數(shù)據(jù)組2的均值和協(xié)方差矩陣比較

由表2 可知，當(dāng)小批量數(shù)據(jù)中存在一個異常值時，數(shù)據(jù)組1、2 的μˉ、S存在顯著差異，將導(dǎo)致控制圖監(jiān)控失效，證明了傳統(tǒng)HotellingT2 控制圖在批量生產(chǎn)過程中易受異常值影響。因此，有必要研究一種能滿足小批量加工過程的工序質(zhì)量控制圖，提升計算控制圖參數(shù)的準(zhǔn)確性和合理程度。

2 基于貝葉斯理論的HotellingT2 工序質(zhì)量控制圖

貝葉斯理論是基于先驗信息、樣本信息等的一種統(tǒng)計推斷［13-14］。首先依據(jù)各種信息推導(dǎo)先驗分布；然后綜合樣本信息，根據(jù)貝葉斯公式得到后驗分布；最后根據(jù)后驗分布推斷未知參數(shù)。該方法充分利用歷史、樣本等數(shù)據(jù)，可極大提升控制圖參數(shù)的計算準(zhǔn)確性和合理程度。

2.1 質(zhì)量數(shù)據(jù)標(biāo)準(zhǔn)化

本文收集過去一定時間段內(nèi)滿足要求的不同批次歷史質(zhì)量數(shù)據(jù)進(jìn)行分析，然而不同批次質(zhì)量標(biāo)準(zhǔn)均不同，無法直接進(jìn)行參數(shù)求解，需要對質(zhì)量數(shù)據(jù)進(jìn)行標(biāo)準(zhǔn)化處理，使各批次質(zhì)量數(shù)據(jù)滿足分布要求。

數(shù)據(jù)標(biāo)準(zhǔn)化處理的常用方法包括標(biāo)準(zhǔn)化、歸一化和中心化。其中，標(biāo)準(zhǔn)化和中心化適用于標(biāo)準(zhǔn)值為中間數(shù)的數(shù)據(jù)標(biāo)準(zhǔn)化處理；歸一化適用標(biāo)準(zhǔn)值為最小值的數(shù)據(jù)標(biāo)準(zhǔn)化；相對公差法［15］處理后的質(zhì)量數(shù)據(jù)可直接反映相對標(biāo)準(zhǔn)值的誤差大小，數(shù)據(jù)處理方法如公式（7）所示。

2.2 先驗信息提取

將標(biāo)準(zhǔn)化處理后的歷史質(zhì)量數(shù)據(jù)作為貝葉斯參數(shù)估計的先驗信息，即抽樣前相同受控條件下的一些歷史質(zhì)量信息。一般情況下，先驗信息主要來源于相同工序、相同設(shè)備，不同批次的歷史加工質(zhì)量數(shù)據(jù)。設(shè)有m個批次的加工質(zhì)量信息，且每個批次需要控制p個質(zhì)量指標(biāo)，經(jīng)過標(biāo)準(zhǔn)化處理后的各批次歷史質(zhì)量數(shù)據(jù)服從正態(tài)分布，i=1，2，…，m；j=1，2，…，p，則每個批次不同質(zhì)量指標(biāo)的均值μij和方差為：

式中：g為第i批歷史質(zhì)量數(shù)據(jù)的樣本容量。

由于小批量生產(chǎn)的條件多變，與當(dāng)前批次更接近的歷史批次更能反映當(dāng)前生產(chǎn)狀況，因此根據(jù)加工序列對各批次歷史質(zhì)量數(shù)據(jù)進(jìn)行指數(shù)加權(quán)。由于在不同批次加工過程中的數(shù)據(jù)采集量不盡相同，因此采用變樣本容量獲得歷史質(zhì)量數(shù)據(jù)不同批次的組間均值和組間方差。

式中：k為加權(quán)系數(shù)，通常取(N+14)/14；N為各批次生產(chǎn)的平均天數(shù)；m為歷史質(zhì)量數(shù)據(jù)總批次。

綜上，本文將先驗信息得到的分布作為工序質(zhì)量控制的先驗分布，即先前不同批次受控狀態(tài)下該工序加工過程的多元質(zhì)量指標(biāo)分布狀況，為合理進(jìn)行貝葉斯參數(shù)估計提供信息基礎(chǔ)。

2.3 貝葉斯理論參數(shù)估計

由于現(xiàn)有樣本Y 的第j個質(zhì)量指標(biāo)均值μj與方差間存在相互影響，需要以先驗分布的組間均值和組間方差對現(xiàn)有樣本Y 的均值μj和方差進(jìn)行估計。根據(jù)共軛貝葉斯理論［13］，設(shè)方差的先驗分布為逆伽馬（逆Gamma）分布，記為IGa(α，β)，設(shè)μj的先驗分布為N(a，b2)，故其共軛先驗分布具有的乘積形式為π()。其中：

式中：a為先驗分布的均值，b2為先驗分布的方差；α為先驗?zāi)尜ゑR分布的形狀參數(shù)，β為先驗?zāi)尜ゑR分布的尺度參數(shù)。

2.3.1 方差估計值

依據(jù)逆Gamma 分布的相關(guān)性質(zhì)，先驗分布～I(xiàn)Ga(α，β)，則數(shù)學(xué)期望E和方差Var分別為：

聯(lián)立式（13）、式（14）可得：

在獲取現(xiàn)有樣本Y 后，依據(jù)共軛函數(shù)分布特性可得先驗分布形式與后驗分布形式一致。故可設(shè)：

式中：α'為后驗?zāi)尜ゑR分布|Y 的形狀參數(shù)；β'為后驗?zāi)尜ゑR分布|Y 的尺度參數(shù)；n為現(xiàn)有樣本Y的容量；為現(xiàn)有樣本Y中第j個質(zhì)量指標(biāo)的數(shù)據(jù)均值；k0為先驗信息在實時批次中所占的分量，在小批量生產(chǎn)過程中一般取常數(shù)，k0=30。

因此，可得的貝葉斯估計值為：

2.3.2 均值估計值

依據(jù)正態(tài)分布的相關(guān)性質(zhì)，先驗分布μj|～N(a，b2)的數(shù)學(xué)期望為：

獲取現(xiàn)有樣本Y 后，依據(jù)共軛貝葉斯的相關(guān)性質(zhì)可知：

式中：a'為后驗分布μj|，Y的均值；b'2為后驗分布μj|，Y的方差。

(μj|)聯(lián)合先驗密度函數(shù)π(μj|)為：

μj的邊緣后驗密度函數(shù)π(μj|Y)為：

本文自由度為(2α+n)的t分布，其中：

均值μj的貝葉斯估計值μjB為：

根據(jù)貝葉斯理論進(jìn)行參數(shù)估計得到均值估計向量μB=(μ1B，μ2B，…，μpB)和方差估計向量，然后結(jié)合式（1）、式（13）—式（20）得到p個質(zhì)量指標(biāo)的協(xié)方差矩陣SB=SP×P，以充分利用先驗信息、樣本信息等，在先驗分布的基礎(chǔ)上估計后驗分布參數(shù)，降低小批量加工過程中因樣本量小帶來的控制圖參數(shù)計算誤差，從而提升控制圖監(jiān)控的準(zhǔn)確性。

2.4 基于貝葉斯理論的HotellingT2工序質(zhì)量控制圖統(tǒng)計量

基于貝葉斯理論的HotellingT2 工序質(zhì)量控制圖設(shè)計理念是利用貝葉斯理論對小批量生產(chǎn)的質(zhì)量數(shù)據(jù)進(jìn)行參數(shù)估計，得到貝葉斯估計均值向量μB和協(xié)方差矩陣SB，再結(jié)合傳統(tǒng)方法設(shè)計多元工序質(zhì)量控制圖，以監(jiān)控加工過程中多元質(zhì)量指標(biāo)的波動狀態(tài)。

具體的，根據(jù)式（1）、式（7）—式（27）得到貝葉斯均值估計向量μB和估計協(xié)方差矩陣SB，結(jié)合傳統(tǒng)式（3）構(gòu)造出反映多元質(zhì)量指標(biāo)波動狀況的綜合指標(biāo)，基于貝葉斯理論的HotellingT2工序質(zhì)量控制圖統(tǒng)計量為：

式中：為貝葉斯估計協(xié)方差矩陣SB的逆矩陣為當(dāng)前加工過程中第i個工件的p個質(zhì)量指標(biāo)的測量值標(biāo)準(zhǔn)化處理后組成的向量；μB為p個質(zhì)量指標(biāo)貝葉斯估計均值向量。

由式（4）—式（6）可知，T2控制圖統(tǒng)計量服從具有參數(shù)p和m-p的F 分布，取一定的虛報概率α即可確定控制UCL 和LCL。當(dāng)T2(i) ≥0 時只需根據(jù)式（5）確定其控制上限UCL，控制下限LCL為零。

綜上，使用基于貝葉斯理論的HotellingT2 工序質(zhì)量控制圖對小批量加工過程的多元質(zhì)量指標(biāo)進(jìn)行監(jiān)控的原理為：若加工過程的T2統(tǒng)計量處于UCL 和LCL 內(nèi)，則說明加工過程處于受控狀態(tài)，否則過程處于不受控狀態(tài)。本文所提方法能從歷史數(shù)據(jù)中提取當(dāng)前產(chǎn)品的質(zhì)量信息，充分利用先驗信息、樣本信息等，在先驗分布的基礎(chǔ)上估計后驗分布參數(shù)，降低小批量加工過程中因樣本量小而帶來的統(tǒng)計誤差，進(jìn)一步提升控制圖監(jiān)控的準(zhǔn)確性。

3 實例分析

為了本文方法的有效性和優(yōu)越性，計算分析某發(fā)動機(jī)凸輪軸支承軸直徑、各支撐軸的同軸度和支撐軸圓度，如表3所示。

Table 3 Measurement of the diameter of the engine camshaft support shaft，coaxiality of each support shaft，and roundness of the support shaft表3 發(fā)動機(jī)凸輪軸支承軸直徑、各支撐軸的同軸度和支撐軸圓度測量值

表3 中x1、x2、x3分別表示直徑、同軸度和圓度；、為標(biāo)準(zhǔn)化處理后的數(shù)據(jù)，可比較傳統(tǒng)方法和基于貝葉斯理論的HotellingT2 工序質(zhì)量控制圖對加工過程的質(zhì)量監(jiān)控能力差異。

凸輪軸支承軸直徑標(biāo)準(zhǔn)值為30.01 mm，上限為30.023 5 mm，下限為30.000 5 mm；各支撐軸的同軸度標(biāo)準(zhǔn)值為0，上限為0.05 mm，下限為0；支撐軸圓度標(biāo)準(zhǔn)值為0，上限為0.4 μm，下限為0，樣本中8、17、20 號為異常多元數(shù)據(jù)。表4 為最近3 批不同型號產(chǎn)品加工過程該工序多元質(zhì)量指標(biāo)測量值標(biāo)準(zhǔn)化處理后的分布情況。其中，質(zhì)量指標(biāo)的均值表示穩(wěn)定狀態(tài)的控制平均值，協(xié)方差矩陣反映穩(wěn)定狀態(tài)多元質(zhì)量指標(biāo)波動收斂程度。

Table 4 Standardized distribution information of multiple quality indicators of historical batch processes表4 歷史批次工序多元質(zhì)量指標(biāo)標(biāo)準(zhǔn)化分布信息

3.1 均值與方差

本文分別使用傳統(tǒng)方法和貝葉斯估計方法計算凸輪軸加工過程中各個指標(biāo)的均值與方差，結(jié)果如表5 所示。由此可知，在異常值影響下本文方法的均值為0.51、0.46，更接近穩(wěn)定歷史批次穩(wěn)定狀態(tài)均值，第2、3 質(zhì)量指標(biāo)的方差相較于傳統(tǒng)方法更小，證明了基于貝葉斯理論的統(tǒng)計方法較好地削弱了異常值對計算結(jié)果的影響。

Table 5 Comparison of calculation results of mean and variance between two methods表5 兩種方法均值和方差計算結(jié)果比較

3.2 協(xié)方差矩陣

兩種方法計算得到凸輪軸直徑、同軸度和圓度的協(xié)方差矩陣結(jié)果如表6 所示。由此可知，本文所提方法的第1、2 質(zhì)量指標(biāo)、第2、3 質(zhì)量指標(biāo)的協(xié)方差和第2 質(zhì)量指標(biāo)的方差遠(yuǎn)小于傳統(tǒng)方法，整體上由貝葉斯參數(shù)估計方法計算得到的協(xié)方差矩陣更小，證明了基于貝葉斯理論的統(tǒng)計方法具有較好的抗干擾能力。

Table 6 Comparison of covariance matrix calculation results between two methods表6 兩種方法協(xié)方差矩陣計算結(jié)果比較

3.3 質(zhì)量指標(biāo)統(tǒng)計量

根據(jù)式（3）、式（28）分別計算傳統(tǒng)方法和基于貝葉斯理論的HotellingT2工序質(zhì)量控制圖的統(tǒng)計量，如圖1所示。由此可見，傳統(tǒng)方法只監(jiān)測出第17 號異常多元質(zhì)量數(shù)據(jù)，存在漏報現(xiàn)象，證明傳統(tǒng)方法對異常值較為敏感，易受到異常值影響，導(dǎo)致監(jiān)控效果與實際情況存在偏差。基于貝葉斯理論的HotellingT2 工序質(zhì)量控制圖能識別出所有異常值，即第8、17、20 號異常多元質(zhì)量數(shù)據(jù)，并對異常情況發(fā)出警報，且中位數(shù)（2.44）小于傳統(tǒng)方法中位數(shù)（2.69），證明基于貝葉斯理論的HotellingT2 工序質(zhì)量控制圖收斂程度更好，對異常值影響具有良好的抵御能力，能準(zhǔn)確監(jiān)測失控狀態(tài)。

Fig.1 Comparison of quality indicator statistics圖1 質(zhì)量指標(biāo)統(tǒng)計量比較

綜上，在小批量加工過程中使用傳統(tǒng)方法的監(jiān)控效果并不理想，基于貝葉斯理論的HotellingT2 工序質(zhì)量控制圖的監(jiān)控效果明顯更可靠、準(zhǔn)確，能及時發(fā)現(xiàn)加工過程異常情況，并能提供較好的解決方法。

4 結(jié)語

傳統(tǒng)HotellingT2 控制圖作為加工過程質(zhì)量監(jiān)控的重要工具之一，監(jiān)控效果可作為評判加工過程是否處于受控狀態(tài)的重要依據(jù)。然而，在小批量加工過程中，由于樣本數(shù)據(jù)較小，傳統(tǒng)方法易受異常值影響，導(dǎo)致監(jiān)控準(zhǔn)確度下降，從而誤判加工過程的受控狀態(tài)，進(jìn)而給管理人員帶來錯誤的信息。

為此，本文結(jié)合貝葉斯理論與傳統(tǒng)方法，構(gòu)建基于貝葉斯理論的HotellingT2 工序質(zhì)量控制圖，以降低或消除小批量加工過程中異常值對監(jiān)控效果的影響，從而使監(jiān)控效果更符合實際，并通過發(fā)動機(jī)凸輪軸加工工序質(zhì)量控制實例分析。實例表明，當(dāng)數(shù)據(jù)中存在異常值時，傳統(tǒng)方法易出現(xiàn)漏發(fā)報警情況，而基于貝葉斯理論的HotellingT2 工序質(zhì)量控制圖可檢測所有異常狀況并發(fā)出警報，具有較好的穩(wěn)定性和準(zhǔn)確性，能更好地應(yīng)用于小批量加工過程中的工序質(zhì)量控制。