基于多重相关状态采样的多元EWMA控制图

刘欣怡　朱永忠

摘 要：为提高控制图的监测效率，提出了一种基于多重相关状态采样的多元EWMA控制图，并利用改进后的马尔可夫链方法计算控制图的平均运行长度。根据不同参数下控制图的平均运行长度，分析了控制图在失控和受控状态下的性能表现，并与其它多元EWMA控制图进行比较。模拟结果表明，该控制图具有良好的监测能力。最后用一组模拟数据来说明该方法的使用。

关键词：多重相关状态采样;多元EWMA控制图;平均运行长度;马尔可夫链

中图分类号：O212      文献标识码：A

A New Multivariate EWMA Control Chart Using

Multiple Dependent State Sampling

LIU Xin-yi， ZHU Yong-zhong

（School of Science， Hohai University， Nanjing， Jiangsu 211100， China）

Abstract：To improve the efficiency of control chart， a new multivariate EWMA control chart has been developed using multiple dependent state sampling. The parameter of the proposed control chart as well as the average run length has been evaluated for different process settings using an improved Markov chain model. The proposed control chart has been evaluated using the average run length of the in-control and out-of-control process. The comparison of performance of the proposed and existing control chart has been studied. It has been shown that the proposed control chart using MDS sampling is more sensitive than the existing control chart. A simulated data is used to illustrate the procedure of the proposed control chart.

Key words：multiple dependent state sampling; multivariate EWMA control chart; average run length; Markov chain model

在产品质量控制和管理中，统计过程控制（Statistical Process Control， SPC）起着很重要的作用，其最重要的方法就是绘制统计过程控制图。現今，很多工业过程需要同时监测多个相关的质量特性，因此很多研究者着眼于多元控制图的研究。Hotelling[1]首先提出了用于监测较大偏移的T2控制图，随后Crosier[2]和Lowry等[3]分别于1988年和1992年提出了多元累计和（Multivariate Cumulation Sum， MCUSUM）和多元指数加权移动平均（Multivariate Exponentially Weighted Moving Average， MEWMA）控制图。MacGregor和Kourti[4]，Lowry和Montgomery[5]及Montgomery[6]详细介绍了已有的几种多元控制图。

为了进一步提高控制图的灵敏度，很多学者提出了利用各种采样方法作为运行准则的控制图。Wortham和Baker[7]于1976年首先提出了多重相关状态（Multiple Dependent State， MDS）采样方法。Aslam等[8]于2015年将MDS采样和传统的计数控制图相结合，提出了一种基于MDS的计数控制图。Aldosari等[9]将MDS采样方法运用到了变化抽样区间的控制图中。2016年，Aslam等[10]提出了一种基于MDS的用于监测服从指数分布的质量特性的控制图。Raza等[11]研究了当质量特性服从多元泊松分布时，采用MDS方法对控制图的影响。最近，Khan等[12]提出了一种基于MDS的自适应控制图。以上这些控制图都表明了利用MDS采样方法会提高控制图的灵敏度。但该方法还没有应用于质量特性服从多元正态分布的情形。因此将MDS采样方法拓展到质量特性服从正态分布的问题中，利用EWMA统计量构造一个新的控制图（MDS-MEWMA图），通过计算该控制图的平均运行长度（Average Run Length，ARL），并对比其他MEWMA控制图的ARL，说明该方法的有效性。

结构安排如下：第1节介绍MDS采样方法，并建立基于MDS的MEWMA控制图。第2节通过马尔科夫链方法计算MDS-MEWMA图的ARL，并通过和其他控制图的比较来分析其性能表现。第3节通过模拟来具体说明该控制图的应用。最后第4节对文章进行了总结。

1 基于MDS的MEWMA控制图

1.1 多重相关状态采样方法的描述

多重相关状态采样方法最早由Wortham和Baker[7]提出，用于判断某一批次生产的产品是否合格，其核心思想为延迟决策。设MDS（r，b，m）表示一个判断准则，其中r（r≥0）表示无条件接受的最大缺陷数;b（b>0）表示当有附加条件时，可额外接受的缺陷的最大数量;m（m>0）表示当有附加条件时，延迟决策所需要的未来批次的数量。具体过程如下：

第一步：对于批次k，从当前批次中选择n个单位的随机样本，并确定缺陷数c。

第二步：若c≤r，则接受该批次，认为该批次生产的产品为合格品;

若r+1≤c≤r+b，则推迟决策直到得到接下来m个批次的样本，并确定每个批次样本的缺陷数c（c1，c2，…cm）。若对j=1，2，…，m，都有cj≤r，则接受当前批次，认为该批次生产的产品为合格品;若有任何一个批次的产品被判定为不合格品，则拒绝当前批次，认为当前批次生产的产品不合格;

若c>r，则拒绝当前批次，认为该批次生产的产品不合格。

第三步：令k=k+1，并返回Step 1。

该方法将当前批次的信息和未来批次的信息结合，使决策的结果更加准确和可靠。

1.2 基于MDS的MEWMA控制图

假设某一生产过程需要监测的质量特性共有p（p≥2）个，由这些质量特性组成的随机向量X=（X1，X2，…，Xp）′服从均值为μ=（μ1，μ2，…，μp）′，协方差阵为Σ的二元正态分布，即X～Np（μ，Σ）。

定义第i个批次的MEWMA统计量Zi为：Zi=R（Xi-μ）+（I-R）Zi-1，其中Z0=（0，0，…，0）′，R=diag（r1，r2，…，rp）′，（0≤rj≤1，1≤j≤p）为平滑系数矩阵;I为单位阵。MEWMA控制图可以针对不同的质量特性Xi选择不同的平滑系数。当生产过程出现了较大的变动时，可以选择较小的rj，而当出现较小变动时，则需要选择较大的平滑系数。通常选择r1=r2=…=rp，在这种情况下，MEWMA统计量可简化为：

Zi=r（Xi-μ）+（1-r）Zi-1（1）

其协方差矩阵为ΣZi=r1-（1-r）2i2-rΣ。当生产的批次足够多，即i足够大时，协方差矩阵可近似为：

ΣZi=r2-rΣ（2）

将统计量T2i=Z'iΣ-1ZiZi繪制于MDS-MEWMA图中，并通过MDS采样方法判断过程是否受控。Lowry等[3]证明了MDS-MEWMA图的ARL值只与由均值向量μ和协方差阵Σ组成的非中心参数λ=（μ'Σ-1μ）1/2有关。

绘制基于MDS采样方法的MEWMA控制图的步骤如下：

第一步：选定第i个生产过程得到的样本，测量其要监测的质量特性Xi，并计算其MEWMA统计量Zi。

第二步：若Zi≤h1，则判定该生产过程处于受控状态;若Zi≥h2，则为失控状态。否则，执行第三步。

第三步：若h1≤Zi≤h2，则推迟决策直到得到接下来的k个生产过程的样本：若这k个生产过程都为受控状态，即Zj≤h1，j=1，2，…，k，则认为当前生产过程处于受控状态;否则认为该生产过程为失控状态。

其中，k称为MDS参数，由不同的生产过程决定。当k=0时，MDS-MEWMA控制图就是由Lowry等提出的传统的MEWMA控制图。

传统的MEWMA控制图只有上控制限，而MDS-MEWMA控制图在上控制限的基础上又增加了一条控制限，两条控制限分别为内控制限（h1）和外控制限（h2），可根据控制图的平均运行长度ARL来确定。通常控制图有两类ARL，即受控的ARL（ARL0）和失控的ARL（ARL1）。ARL的常用计算方法主要有马尔可夫链方法，积分方程法和随机模拟法。基于已有研究方法，采用马尔可夫链方法来计算。

在计算得到状态空间后，可通过转移概率矩阵计算出ARL，但该控制图的转移概率矩阵具有较强的稀疏性，在此基础上进行计算会很困难。Woodall和Reynolds[13]针对转移概率矩阵具有较强的稀疏性的情况，提出了一种新的计算ARL的方法。基于此方法可计算出MDS-MEWMA控制图的ARL。

为了确定该控制图的状态空间，定义W′=（W1，W2，…，Wk+1），X′=（X1，X2，…，Xk+1），其中

Wj=1，Zj∈（h1，h2）0，其他     ，

Xj=Wj，∑jh=1Wh

此时，X′=（X1，X2，…，Xk+1）表示此控制图的一个状态。得到控制图的状态空间后，可以很容易地计算出ARL。

以k=3的控制图为例，将整个空间划分为三个区域：（0，h1），（h1，h2），（h2，+SymboleB@），则该控制图共有7个状态，如表1所示。设pj（j=1，2，3）表示监测统计量落在第j个区域内的概率，则可得到该控制图运行长度的概率分布。

对于初始状态为i的控制图，可通过如下公式来计算其ARL：

ARL∑nh=1hP（Ni=h）+λ^P（Ni=n）×n1-λ^+1（1-λ^）2

其中λ^=1-∑nh=1P（Ni=h）1-∑n-1h=1P（Ni=h），n=ARL0，P（Ni=h）表示初始状态为i，经过h步到达吸收态的概率。

2 MDS-MEWMA控制圖的性能表现

通过平均运行长度来分析控制图的性能表现，将其与已有的MEWMA控制图进行比较，以度量他们的监测效率。设定控制图受控状态下的ARL为200，计算并比较ARL1。

计算p=2，k=1，2，3，r的取值在0.05到3之间时MDS-MEWMA控制图的ARL1值，结果见表2，表3，表4。分析后可知，当固定MDS参数k的值时，随着平滑系数r的减小，ARL的值在减小。如表2所示，当固定λ=0.25，r=0.2时，控制图的ARL为65.89，而r=0.8时该控制图的ARL值为105.64。当其他参数不变，MDS参数k的值增加时，该控制图的ARL也会减小。表2中k=1，λ=0.5，r=0.4时，该控制图ARL的值为38.42，而表3中其他条件相同，k=2时，该控制图的ARL值为32.42。

图1给出了k=1时不同的平滑系数下，MDS-MEWMA控制图的ARL的比较。从图中可以看出该控制图在r=0.2时，其监控效果相对较好。因此，在实际生产过程中，可以选择较大的k值和较小的r值来获得更好地监控效果。

为了进一步研究MDS-MEWMA控制图的表现性能，不将同参数下MDS-MEWMA控制图和已有MEWMA控制图的ARL进行比较。设定每种控制图具有相同的受控ARL，并比较失控ARL，如图2所示。图2中k=0的曲线代表了传统的MEWMA控制图，从图中可以看出基于MDS采样方法的MEWMA控制图在ARL方面表现更好，对于较小变动的监测更有效。

3 数据案例

下面用一组模拟数据来说明MDS-MEWMA控制图的使用。首先利用MATLAB生成50组服从二元正态分布的随机数组，其中前35组为处于受控过程的均值为μ=（0，0），协方差阵为Σ=1001的质量特性，后15组为偏移程度为0.5的失控状态下的质量特性。

取k=3，r=0.2，绘制MDS-MEWMA控制图，结果见图3。该控制图绘制过程与传统的MEWMA控制图做法一样，不同之处在于对于失控状态的判断，故将两种控制图绘制在一张图中，以便比较。图3中有三条控制限，分别为h，h1和h2。其中h为传统的MEWMA控制图的控制限，h1和h2分别为MDS-MEWMA控制图的内控制限和外控制限。从图中可以看出，当使用MEWMA控制图对该过程进行监控时，失控信号出现在第41个点，而当使用MDS-MEWMA控制图时，在第38个点处就已经发出了失控信号，比传统控制图提前了3个点。因此MDS-MEWMA控制图可以更有效的对较小的过程偏移进行监测。

4 结 论

为了充分利用未来信息，更灵敏地监测较小的过程偏移，讨论了基于MDS采样方法的一种新的MEWMA控制图。运用马尔可夫链方法估计了该控制图的参数，计算平均运行长度ARL，并和传统的MEWMA控制图进行比较。研究表明MDS-MEWMA控制图的监控效率优于一般的MEWMA控制图。当过程失控时，MDS-MEWMA控制图能更快的发出失控信号，提高生产效率。最后利用一个模拟数据来具体说明该控制图的应用，同时也表明MDS-MEWMA控制图的监控效果优于传统的MEWMA控制图。未来的研究可着眼于对基于MDS采样方法的控制图的优化方面，对其中的参数进行调整，以便达到更好地监测效果。

参考文献

[1] HOTELLING H. Multivariate quality control illustrated by air testing of sample bombsights[J]. Techniques of Statistical Analysis， 1947：11-184.

[2] CROSIER R B. Multivariate generalizations of cumulative sum quality-control schemes[J]. Technometrics， 1988， 30（3）：291-303.

[3] LOWRY C A， WOODALL W H. A multivariate exponentially weighted moving average control chart[J]. Technometrics， 1992， 34（1）： 46-53.

[4] MACGREGOR J F， KOURTI T. Statistical process control of multivariate processes[J]. Control Eng Practice， 1995， 3（3）：403-414.

[5] LOWRY C A， MONTGOMERY D C. A review of multivariate control charts[J]. A I I E Transactions， 1995， 27（6）：800-810.

[6] MONTGOMERY D C. Introduction to statistical quality control[M]. New York： John Wiley &Sons， 2005.

[7] WORTHAM A W， BAKER R C. Multiple deferred state sampling inspection[J]. International Journal of Production Research， 1976， 14（6）：719-731.

[8] ASLAM M， NAZIR A， JUN C H. A new attribute control chart using multiple dependent state sampling[J]. Transactions of the Institute of Measurement and Control， 2015， 37（4）：569-576.

[9] ALDOSARI M S， ASLAM M， JUN C H. A new attribute control chart using multiple dependent state repetitive sampling[J]. IEEE Access， 2017， PP（1）：1-1.

[10]ASLAM M， AHMAD L， JUN C H， et al. A control chart for an exponential distribution using multiple dependent state sampling[J]. Quality & Quantity， 2016， 32（8）：2803-2812.

[11]ALI R M， MUHAMMAD A. Design of control charts for multivariate poisson distribution using generalized multiple dependent state sampling[J]. Quality Technology and Quantitative Management， 2019， 16：629-650.

[12]KHAN N， ASLAM M. Design of an EWMA adaptive control chart using MDS sampling[J]. Journal of Statistics and Management Systems， 2019， 22（3）：535-555.

[13]WOODALL W H， REYNOLDS M R. A discrete Markov chain representation of the sequential probability ratio test[J]. Communications in Statistics-Sequential Analysis. 1983， 2（1）：27-44.