基于高斯混合模型的非参数控制图

Nonparametric Control Chart Based on Gaussian Mixture Model

摘要: 在工业生产过程中，由于设备参数、原料特性以及工艺条件的异质性，过程数据往往呈现出复杂的多模态分布特征。这类数据结构超出传统统计过程控制方法所依赖的单模态或同分布假设的适用范围。多模态数据复杂的多峰特性，使得基于参数化分布假设的传统建模方法难以准确刻画数据的真实分布，进而影响过程监控。针对这种情况，本文提出了一种数据驱动的监控方法。首先，引入高斯混合模型对多模态过程数据进行概率密度建模，克服了传统方法在分布形式上的先验假设局限；其次，基于该模型构建负对数似然统计量，并结合指数加权移动平均策略，实现对过程异常的检测。数值仿真与实例验证了本文所提控制图的可行性和实用价值，为多模态过程的质量监控提供了一种有效的新思路。

Abstract: In industrial production processes, process data often exhibit complex multimodal distribution characteristics due to heterogeneity in equipment parameters, raw material properties, and process conditions. Such data structures exceed the applicability of traditional statistical process control methods, which typically rely on unimodal or identically distributed assumptions. The intrinsic multi-peak nature of multimodal data makes it difficult for conventional parametric modeling approaches to accurately capture the true data distribution, thereby compromising process monitoring performance. To address this issue, this paper proposes a data-driven monitoring method. First, a Gaussian mixture model is introduced to model the probability density of multimodal process data, overcoming the limitations of prior distributional assumptions inherent in traditional methods. Based on this model, a negative log-likelihood statistic is constructed and integrated with an exponentially weighted moving average strategy to detect process anomalies. Numerical simulations and real-world case studies demonstrate the feasibility and practical value of the proposed control chart, offering an effective new approach for quality monitoring in multimodal processes.

1. 引言

统计过程控制(SPC)作为一种过程稳定性监控技术，已在工业制造领域得到广泛应用。传统SPC控制图(例如Shewhart $\overline{X}$ 图)的理论基础建立在一个关键假设之上，受控过程数据服从单模态概率分布，即其概率密度函数(Probability Density Function, PDF)呈现单一峰值特征[1]。然而，实际工业过程往往表现出显著的时变特性，包含动态和非线性等非平稳行为特征[2]。当生产过程涉及多种操作模式切换时，采集到的过程数据可能呈现多模态分布特性，即其PDF中存在多个峰值。在本文中，这种由多峰分布所生成的过程数据将被称为多模态数据。

多模态过程数据在复杂工业系统中都有出现，由于这类数据不满足单模态和同分布的基本假设，因此难以用某一具体的参数分布来准确建模。多模态过程数据的分布因存在多个局部最大值，其密度曲线结构比单模态分布更为复杂。这使得传统SPC方法在处理多模态数据时可能表现不佳。在设计监控、监督和控制系统时，必须充分考虑过程数据的特性。然而，专门针对多模态数据设计的控制图在文献中有限。Quiñones-Grueiro等人[2]对多模态过程的研究进行了总结，部分学者提出先通过转换函数将多模态特征去除，进而使用常规的监控方法解决问题。例如Deng等人[3]通过采用服从单峰分布的局部概率密度估计，将多模态数据转化为概率密度。然而，对于复杂系统而言，转换函数很难准确确定，并且数据转换的计算成本可能较高，为每个观测定义标准化模型是具有挑战性的。

另一种可行的选择是非参数或无分布假设的控制图，这类方法虽不依赖于某一具体的参数分布形式，但是，这些非参数控制图通常依赖于过程是单模态的假设，即数据是源自相同分布，因此不适用于存在多模态过程的模式例如，Hackl和Ledolter [4]提出的单个观测值的指数加权移动平均(EWMA)图是由观测值的“标准化秩”构建的，这是由受控(In-Control, IC)分布决定的。若分布未知，他们建议使用收集到的参考数据中的排名代替。该方法侧重于检测位置参数的变化，但是在监测整个多模态数据过程分布是否发生偏移时，其余的分布特征，如尺度参数，局部最大值的个数也是重要的质量指标。Qiu和Li [5]提出的对观测数据进行分类，然后利用一定的分类数据分析统计程序构造非参数SPC图。这种基于数据分类的方法会因为分类数量较多时容易导致列联表稀疏，监控方法在一定程度上会造成信息丢失。鉴于上述问题，针对单变量多模态过程开发既稳健又具有普适性的质量监控方法具有一定的研究意义。

本文提出了一种基于密度估计的数据驱动的监测方案，用于单变量多模态数据的在线监测。密度估计已广泛应用于各种应用中。通过估计给定数据集的概率密度函数(PDF)，可以了解数据的基本分布。基于GMM的方法对多模态过程进行建模，并且考虑了不同的模型估计方法[6] [7]。本文针对单变量多模态过程的质量监控问题，提出了一种基于高斯混合模型的负对数似然值的EWMA监控方案。这种数据驱动的监控方法，能够有效地检测模型的参数偏移和模型的组分变化。本文的其余部分组织如下。第2节介绍了数据驱动的方法，包括高斯混合模型，模型估计的方法和监控方案的构建。第3节通过蒙特卡洛模拟对所提出的监测方案在不同失控情景下进行了比较分析。在第4节中讨论了在实际案例中所提出的监控方案的有效性。最后对论文进行总结并给出了未来的可能的研究方向。

2. 方法

假设$X=\left\{x_1,x_2,\cdots ,x_t\right\}$ 为一组单变量多模态数据，$x_t$ 是在$t$ 时刻的观察值。在对多模态数据的第二阶段监测过程中，为了监测一个连续的一元多模态过程$\left\{x_t,t\ge 1\right\}$ ，假设存在一个IC数据集$\left\{x_1,x_2,\cdots ,x_n\right\}$ ，其样本量为$n$ ，预先从IC分布$f_0$ 中采集。我们的目的是确定过程是否已从IC状态变为失控(Out-of-Control, OC)状态，即可表示为假设检验：

$H_0:\left\{x_1,x_2,\cdots ,x_t\right\}$ 遵循受控分布$f_0$ ；

$H_1:\exists \tau \in \left[1,t\right],\left\{x_1,x_2,\cdots ,x_{\tau }\right\}$ 遵循受控分布$f_0$ ，$\left\{x_{\tau },\cdots ,x_t\right\}$ 遵循失控分布$f_{\text{1}}$ 。

其中τ为变化点，$f_{\text{1}}$ 为失控分布，满足$f_0\ne f_{\text{1}}$ 。对于含有多个峰的多模态过程[7]，指出其备择假设$H_1$ 可能存在两种类型的过程OC情景，第一类是关于分布中已有峰的参数发生偏移，第二类是关于分布的峰数发生变化。因此，在多模态过程质量监控方法应该具有可以及时准确地监测出这两种OC类型的能力。

本文提出的监测单变量多模态数据的方案包括以下三个步骤：

步骤1：使用高斯混合模型方法对多模态过程进行建模。

步骤2：讨论不同的参数估计方法，并进行比较分析。

步骤3：根据步骤2中得到的最优估计值，构建EWMA的统计量，并进行在线监测。

2.1. 高斯混合模型

高斯混合模型(Gaussian Mixture Model, GMM)是一种概率密度估计方法，它假设数据由若干个高斯分布混合生成。考虑数据集$X=\left\{x_1,x_2,\cdots ,x_n\right\}$ ，每个数据点都假定是由$K$ 个高斯分布的混合而成。高斯混合模型定义为：

$f\left(x\right)={\displaystyle \sum }_{j=1}^K{w}_j\mathcal{N}\left(x|{\theta }_j\right)$ (1)

其中，$w_j$ 为第$j$ 个高斯分量的权重，满足${\displaystyle \sum }_{j=1}^K{w}_j=1$ 且$0\le w_j\le 1$ ；${\theta }_j=\left\{{\mu }_j,{\sigma }_j^2\right\}$ 为第$j$ 个高斯分量的分布参数；$\mathcal{N}\left(x|{\theta }_j\right)$ 是高斯分布，表示对于给定的$x$ ，其均值为${\mu }_j$ ，方差为${\sigma }_j^2$ 。

2.2. 模型估计的方法

GMM中的关键问题是如何选择适当的分量数$K$ 。常用的模型评估准则包括赤池信息量准则(Akaike Information Criterion, AIC) [8]和贝叶斯信息量准则(Bayesian Information Criterion, BIC) [9]：

$AIC=2k-2\ln \left(L\right)$ (2)

$BIC=\ln \left(N\right)k-2\ln \left(L\right)$ (3)

其中，$L$ 为模型的最大似然估计；$N$ 为样本的数量；$k$ 为模型中自由参数的数量(GMM中包括了每个高斯分布的均值、方差和权重)，对于一维GMM，每个高斯成分有两个参数${\mu }_j,{\sigma }_j^2$ ，且每个成分的权重$w_j$ 有$K-1$ 个自由参数，因此自由参数的总数量为：

$k=K\left(2+1\right)-1=3K-1$ (4)

AIC和BIC通过惩罚复杂模型避免过拟合，选择对应最小值的$K$ 作为GMM最优解，分别记为AIC-GMM和BIG-GMM。在AIG-GMM和BIG-GMM模型中，采用期望最大化(Expectation-Maximization, EM)算法[10]估计GMM的参数：

E步：计算分模型$j$ 对观测数据$x_i$ 的响应度${\gamma }_{ij}$ ：

${\gamma }_{ij}=\frac{w_j\mathcal{N}\left(x_i|{\theta }_j\right)}{{\displaystyle \sum _{j=1}^K{w}_j\mathcal{N}\left(x_i|{\theta }_j\right)}},i=1,2,\cdots ,n,j=1,2,\cdots ,K$ (5)

M步：计算新一轮迭代的模型参数：

${\widehat{\mu }}_j=\frac{{\displaystyle \sum _{i=1}^n{\widehat{\gamma }}_{ij}{x}_j}}{{\displaystyle \sum _{i=1}^n{\widehat{\gamma }}_{ij}}},j=1,2,\cdots ,K$ (6)

${\widehat{\sigma }}_j^2=\frac{{\displaystyle \sum _{i=1}^n{\widehat{\gamma }}_{ij}}{\left(x_i-{\mu }_j\right)}^2}{{\displaystyle \sum _{i=1}^n{\widehat{\gamma }}_{ij}}},j=1,2,\cdots ,K$ (7)

${\widehat{\alpha }}_j=\frac{n_j}{n}=\frac{{\displaystyle \sum _{i=1}^n{\widehat{\gamma }}_{ij}}}{n},j=1,2,\cdots ,K$ (8)

最终迭代计算到参数收敛。

AIC和BIC需预设分量数K。在数据量大且符合混合高斯假设时，二者可以有效估计真实分量数。但对于非高斯数据，可能导致组分接近、产生冗余簇。为此，我们考虑了变分贝叶斯狄利克雷过程高斯混合模型(Variational Bayesian Dirichlet Process Gaussian Mixture Model, VBDP-GMM) [6] [11]。VBDP-GMM引入狄利克雷过程(Dirichlet Process, DP)作为先验来自动选择模型的分量数K。DP定义为$G~DP\left(\alpha ,G_0\right)$ ，其中，$\alpha$ 为浓度参数，控制新类别的产生频率，$G_0$ 为基准分布(通常为高斯分布)。在VBDP-GMM模型中，通过变分推断来估计模型的参数$\theta =\left\{{\mu }_j,{\sigma }_j^2,w_j,z_i\right\}$ (各高斯分量的均值${\mu }_j$ ，方差${\sigma }_j^2$ 、权重$w_j$ 以及隐藏变量(数据点$x_i$ 的分配分量$z_i$ )。在变分推断中，找到一个变分分布$q\left(\theta \right)$ 来近似真实的后验分布$p\left(\theta |X\right)$ ，变分目标是最小化变分分布$q\left(\theta \right)$ 与真实后验分布$p\left(\theta |X\right)$ 的KL散度：

$KL\left(q\Vert p\right)={\displaystyle \int q\left(\theta \right)\log \frac{q\left(\theta \right)}{p\left(\theta |X\right)}\text{d}\theta }$ (9)

由于$p\left(\theta |X\right)$ 包含棘手的边际似然$p\left(X\right)$ ，因此直接KL最小化是不可行的。经过式(9)的数学变形，我们得到以下结果：

$\begin{array}{c}L\left(q\right)=E_{q\left(\theta \right)}\left[\log \frac{p\left(X,\theta \right)}{q\left(\theta \right)}\right]\\ =E_{q\left(\theta \right)}\left[\log \frac{p\left(\theta \right)p\left(X|\theta \right)}{q\left(\theta \right)}\right]\\ =E_{q\left(\theta \right)}\left[\log p\left(X|\theta \right)\right]-KL\left(q\left(\theta \right)\Vert p\left(\theta \right)\right)\end{array}$ (10)

因此，最小化KL散度等价于最大化变分下界(Evidence Lower Bound, ELBO)。通过优化ELBO，我们可以得到模型的近似后验分布，从而估计模型参数。

2.3. 监控方案

在本小节中，提出了基于GMM的阶段II的监测方案，以实现对单变量多模态过程的OC状态的有效检测。在阶段II的过程监控中，我们目标是检测过程是否从IC状态偏移到OC状态。由于多模态数据呈现复杂的分布形态(如多峰分布)，直接使用传统的参数化控制图(如Shewhart图、CUSUM或EWMA)可能无法有效检测其变化。根据假设检验中的备择假设$H_1$ 表明该过程可能存在两种类型的过程OC情景，关于模型的参数发生偏移以及模型的峰数发生变化。因此，在多模态过程质量监控方法应该具有可以及时准确地监测出这两种OC类型的能力。

GMM通过对历史正常数据集进行训练，为每个新的观测值$x_i$ 提供密度估计$\widehat{f}\left(x_i\right)$ ，并通过计算对应的对数似然值(Log-Likelihood, LL)来监控过程中的异常变化。当观测值偏离IC状态时，其在当前密度估计模型中的出现概率显著下降，导致对应的LL值减小。为了在多模态过程中实现有效的在线监测，我们将计算新观测值的负对数似然值(NLL)统计量[12]，并与控制限比较，判断是否出现OC。对于第$i$ 个观测值的NLL值统计量定义如下：

$y_i=-{ℒ}_i=-\text{log}\left({\displaystyle \sum }_{j=1}^K{\widehat{w}}_j\mathcal{N}\left(x|{\widehat{\theta }}_j\right)\right)$ (11)

其中，${ℒ}_i=\text{log}\left(\widehat{f}\left(x_i\right)\right)$是第$i$ 个观测值的LL值。

为了适应多模态数据过程建模并进行公平比较，本文构建了基于公式(11)中负对数似然值(NLL)的非参数指数加权移动平均(EWMA)统计量

$Q_i=\lambda y_i+\left(1-\lambda \right)Q_{i-1},i=1,2,\cdots ,n$ (12)

其中，$E_0$ 是初始值，本文设置初始值$E_0=0$ 。$\lambda$ 是平滑系数，一般来说对于EWMA类型的控制图，选择一个相对较小的$\lambda$ 可以更好地检测到小的偏移，反之较大的$\lambda$ 则对较大的偏移更为敏感[13]。这一原理同样适用于本文控制图，为了进行比较研究，本文选择了$\lambda \in \left\{0.05,0.1,0.2\right\}$ 三个值。控制线$H$ 是一个关键阈值，用于判断多模态数据是否偏离了正常状态。我们可以采用二分搜索算法来计算IC数据的控制线$H$ 值。

3. 数值性能比较

3.1. 主要竞争方法

第二个对比方法我们考虑的是将无分布控制图或非参数控制图作为比较方法。Hackl和Ledolter提出了基于“标准秩”的EWMA图表(简称HLE) [4]，这些标准秩是使用受控分布作为历史参考样本计算的。也就是说，在受控分布未知的情况下，可以从受控过程中选取一个大小为$n$ 的历史样本$\left\{x_1,x_2,\cdots ,x_n\right\}$ 作为参考，并定义新观测值$x_t$ 的标准化秩为

$R_t=\frac{2}{n+1}\left(R_t^{\text{*}}-\frac{n}{2}\right)$ (13)

其中，$R_t^{\text{*}}$ 是$x_t$ 关于$n$ 个历史样本$\left\{x_1,x_2,\cdots ,x_n\right\}$ 的秩，即

$R_t^{\text{*}}=1+{\displaystyle \sum }_{i=1}^{n}I\left(x_i\le x_t\right)$ (14)

其中，$I\left(x_i\le x_t\right)$ 是指示性函数，表示当$x_i\le x_t$ 时，$I\left(x_i\le x_t\right)=1$ ，否则为0。于是HLE控制图的统计量为

$T_t=\lambda R_t+\left(1-\lambda \right)T_{t-1},t=1,2,\cdots ,n$ (15)

其中，$T_0=0$ 且$\lambda$ 是一个常数。当$T_t\ge H$ 时，则发出失控警报，$H$ 是由${\text{ARL}}_{\text{0}}$ 得到的控制线。

3.2. 控制图参数设置(模型估计的性能)

(1) 多模态模型的设置

在仿真实验中，首先需要生成服从单变量多模态过程的IC数据集。我们基于不同的分布类型的混合，构建了多模态过程分布模型。本文考虑了两种不同的模型，定义如下：

Model A $f=0.5p_1\left(x\right)+0.5p_2\left(x\right)$

Model B $f=0.2p_1\left(x\right)+0.2p_2\left(x\right)+0.2p_3\left(x\right)+0.2p_4\left(x\right)+0.2p_5\left(x\right)$

其中，$p_j\left(x\right)$ 是样本$x$ 的在第$j$ 簇的密度函数。对于多模态过程分布的每个峰，考虑了非高斯分布(均匀分布)以及高斯分布进行模拟分析，模型参数设置如表1所示。

Table 1. Model parameter settings

表1. 模型参数设置

分布类型	模型	$p_1\left(x\right)$	$p_2\left(x\right)$	$p_3\left(x\right)$	$p_4\left(x\right)$	$p_5\left(x\right)$
均匀分布 $U\left(a,b\right)$	A-1	$U\left(0,2\right)$	$U\left(5,7\right)$
	B-1	$U\left(0,1\right)$	$U\left(3,5\right)$	$U\left(6,7\right)$	$U\left(9,11\right)$	$U\left(12,15\right)$
高斯分布 $N\left(\mu ,{\sigma }^2\right)$	A-2	$N\left(0,1\right)$	$N\left(0,3\right)$

(2) 控制图参数设置和评价指标

在所有图表中，我们设置历史数据集个数$n=300$ ，$AR{L}_0=200$ 。平滑系数$\lambda \in \left\{0.05,0.1,0.2\right\}$ ，初始值E₀设为0。通过10,000次重复仿真，采用二分搜索算法，逼近不同控制图的控制限。控制图的监控性能评价指标，本文主要分析了控制图的OC性能指标——OC平均运行长度(ARL₁)和总体性能指标——相对平均指数RMI [14]。其中RMI定义如下：

$\text{RMI}=\frac{1}{m}{\displaystyle \sum _{i=1}^m\frac{{\text{ARL}}_{{\delta }_i}-\min {\text{ARL}}_{{\delta }_i}}{\min {\text{ARL}}_{{\delta }_i}}}.$ (16)

其中，$m$ 是所考虑的偏移数量，${\text{ARL}}_{{\delta }_i}$ 是给定控制图在偏移为${\delta }_i$ 时的${\text{ARL}}_1$ 值，而$\min {\text{ARL}}_{{\delta }_i}$ 是所有控制图在偏移为${\delta }_i$ 时的最小的${\text{ARL}}_1$ 值。在相同的偏移情景和相等的${\text{ARL}}_{\text{0}}$ 值下，具有较小${\text{ARL}}_1$ 的控制图表现出更佳的性能。从整体性能指标的角度来看，在给定的偏移范围内，具有最小RMI的控制图性能最优。

3.3. 数值模拟比较分析

3.3.1. GMM的模型估计效果的比较结果

在对模型的监控性能比较之前，我们首先对GMM的估计性能进行了对比分析。具体而言，当模型的真实分量数未知时，对于AIC-GMM和BIC-GMM方案，我们预先设定初始分量数K的范围为2到40，而VBDP-GMM利用DP自适应地选择模型中的分量数无需预设分量数。然后对表1中给出的方案进行了估计，最终算法收敛时的分量数如表2所示。

由表2可知，在估计由高斯分布生成的多模态数据时，三种方案都收敛到了真实分量数。然而，上述三种方案都假设每个组分均由高斯分布生成的，因此它们在捕捉高斯分布混合的多模态数据时表现较好。然而，在拟合非高斯数据时，均出现了不同程度的偏差。其中AIC-GMM与BIC-GMM出现的过拟合尤为严重。由于无法准确估计分量数可能会导致GMM不能有效地描述这些复杂数据，于是在后续讨论中，将不再考虑。基于VBDP-GMM的NLL值的EWMA监控方案，记为VGNE。

对于模拟，为了验证本文所提出的控制图的有效性，我们考虑了基于非高斯和高斯多模态过程的两类具有代表性的OC情景：多模态过程中的参数偏移，包括位置参数和尺度参数；多模态过程中的分量数变化。

Table 2. Number of components at convergence for different algorithms

表2. 不同算法收敛时的分量数

序号	真实分量数	AIC-GMM	BIC-GMM	VBDP-GMM
A-1	K = 2	K = 8	K = 4	K = 2
B-1	K = 5	K = 18	K = 9	K = 5
A-2	K = 2	K = 2	K = 2	K = 2

3.3.2. 非高斯多模态过程偏移时的性能

对于模型A-1，OC情景我们分别考虑了：(a) 位置参数均偏移$U\left(a,b\right)+\delta$ ，$\delta \in \left\{0.1,0.2,0.3,0.4,0.5,0.75,1,2\right\}$ ，以模拟过程发生不同幅度的位置偏移；(b) 为了保证在尺度参数发生偏移时，位置参数保持不变，我们设定尺度参数均偏移的失控情况为${\delta }_1\cdot \left(U\left(a,b\right)-\frac{a+b}{2}\right)+\frac{a+b}{2}$ ，${\delta }_1\in \left\{1.1,1.2,1.3,1.4,1.5,1.75,2,3\right\}$ ，表示尺度参数的缩放系数，对应过程标准差放大至原来的1.1倍、1.2倍等；(c) 添加新的组分。表3给出了多个分量的位置参数偏移的性能比较结果；表4给出了多个分量的尺度参数偏移下控制图的性能比较结果；表5给出了对于模型分量数变化的OC示例，表6给出了模型分量数变化OC示例下控制图的性能比较结果。从上述表中可得，无论OC情景是模型参数偏移还是模型组分数变化，VGNE控制图的ARL1和RMI值都达到最低，都可以针对二阶的非高斯多模态过程进行有效且稳健的在线监控。

Table 3. Performance of control charts under location parameter shifts in multiple components

表3. 对于多个分量的位置参数偏移下控制图的性能比较

δ	VGNE			HLE
	λ = 0.05	λ = 0.1	λ = 0.2	λ = 0.05	λ = 0.1	λ = 0.2
0	200.99	199.60	200.48	200.06	200.66	200.63
0.1	149.12	124.36	105.85	95.00	103.22	107.30
0.2	96.49	67.55	47.75	55.46	60.07	65.64
0.3	68.83	42.86	27.60	37.81	40.92	43.85
0.4	52.06	31.08	18.86	28.23	29.94	32.89
0.5	41.65	23.99	14.07	22.34	23.91	25.66
0.75	26.47	14.70	8.43	15.20	15.34	15.91
1	18.28	10.09	5.75	11.92	11.77	11.65
2	6.90	3.86	2.28	8.45	8.49	7.77
RMI	0.521	0.195	0.0086	0.214	0.239	0.257
H	1.670	1.757	1.882	0.165	0.275	0.432

Table 4. Performance of control charts under scale parameter shifts in multiple components

表4. 对于多个分量的尺度参数偏移下控制图的性能比较

δ1	VGNE			HLE
	λ = 0.05	λ = 0.1	λ = 0.2	λ = 0.05	λ = 0.1	λ = 0.2
1	200.99	199.60	200.48	200.06	200.66	200.63
1.1	82.57	58.57	44.40	186.38	183.07	170.92
1.2	56.82	35.40	23.73	178.64	169.30	152.84
1.3	44.86	26.30	16.41	170.28	156.17	137.24
1.4	37.44	21.40	12.84	161.15	148.13	130.32
1.5	32.22	18.05	10.65	155.73	142.56	122.41
1.75	23.70	13.15	7.52	146.95	131.54	107.86
2	18.52	10.25	5.87	141.71	124.68	102.79
3	9.70	5.49	3.25	130.34	108.65	85.26
RMI	0.123	0.043	0.001	0.808	0.729	0.618
H	1.670	1.757	1.882	0.165	0.275	0.432

Table 5. OC case for varying number of model components

表5. 模型分量数变化的OC示例

OC case	wj	模型分量数变化
1-1	0.3, 0.5, 0.2	U(0,2), U(5,7), U(0.5,2.5)
1-2	0.3, 0.5, 0.2	U(0,2), U(5,7), U(1,3)
1-3	0.3, 0.5, 0.2	U(0,2), U(5,7), U(2,4)
2-1	0.5, 0.3, 0.2	U(0,2), U(5,7), U(5.5,7.5)
2-2	0.5, 0.3, 0.2	U(0,2), U(5,7), U(6,8)
2-3	0.5, 0.3, 0.2	U(0,2), U(5,7), U(7,9)

Table 6. Performance of control charts under the OC case with varying numbers of model components

表6. 对于模型分量数变化OC示例下控制图的性能比较

OC case	VGNE			HLE
	λ = 0.05	λ = 0.1	λ = 0.2	λ = 0.05	λ = 0.1	λ = 0.2
1-1	97.80	67.25	47.13	117.25	134.42	159.13
1-2	50.55	30.38	19.45	81.82	102.75	129.81
1-3	25.70	14.60	8.88	62.02	79.62	112.18
2-1	84.92	56.23	37.85	88.37	87.72	86.06
2-2	46.70	27.63	17.63	57.36	56.11	51.54
2-3	22.08	12.79	7.93	43.92	41.75	37.42
RMI	0.345	0.126	0.000	0.580	0.638	0.714
H	1.670	1.757	1.882	0.165	0.275	0.432

除了考虑二阶的情况，本文也考虑了基于均匀分布五阶的多模态过程模型B-1，进一步说明本文所提出的VGNE控制图在多模态过程的监控性能。OC情景仍然考虑模型参数偏移和模型分量数变化，失控情况分别设置为：(a) 位置参数均偏移$U\left(a,b\right)+\delta$ ，(b) 尺度参数均偏移$\delta \cdot \left(U\left(a,b\right)-\frac{a+b}{2}\right)+\frac{a+b}{2}$ ，(c)添加新的组分$U\left(12,15\right)+\delta$ 。由表7给出了基于非高斯多模态过程模型B-1下控制图的监控性能结果，在位置参数偏移较小的情况下，HLE图监控性能优于VGNE图，这是因为HLE图的设计针对位置参数变化更为敏感。但当偏移量增大时，VGNE图的监控性能更优，并且在尺度参数发生偏移以及注入新的组分时，VGNE控制图都可以针对五阶的非高斯多模态过程进行有效且稳健的在线监控。

Table 7. Monitoring performance of control charts under the non-Gaussian multimodal process model B-1

表7. 基于非高斯多模态过程模型B-1下控制图的监控性能

OC case	δ		VGNE			HLE
			λ = 0.05	λ = 0.1	λ = 0.2	λ = 0.05	λ = 0.1	λ = 0.2
位置参数	a	0.2	184.01	167.54	138.08	96.33	105.21	118.22
		0.5	121.73	91.23	67.41	50.54	58.76	69.11
		1	70.80	42.97	28.15	32.45	35.76	41.31
尺度参数	b	1.2	142.86	118.63	96.39	189.74	184.46	187.62
		1.5	107.17	75.40	53.30	186.81	183.81	175.91
		2	82.27	51.32	33.60	183.76	183.70	173.82
组分数	c	0.5	136.67	110.15	89.12	170.37	162.31	159.45
		1	104.83	72.65	50.85	147.97	143.74	136.76
		2	73.50	45.07	27.98	127.27	123.92	112.36
H			2.840	2.923	3.041	0.165	0.275	0.432

3.3.3. 高斯多模态过程偏移时的性能

在本节中我们比较了VGNE和HLE控制图在高斯多模态过程中的监控性能。以表1的模型A-2为例，分析了模型参数偏移和模型分量数变化，OC情况分别设置为：(a) 位置参数整体偏移$N\left(\mu +\delta ,{\sigma }^2\right)$ ，(b) 尺度参数整体偏移$N\left(\mu ,{\delta }^2\right)$ ，(c) 新增组分$N\left(3,{\delta }^2\right)$ 。表8给出具有代表性的结果。

Table 8. Monitoring performance of control charts under the Gaussian multimodal process model A-2

表8. 基于高斯多模态过程模型A-2下控制图的监控性能

OC case			VGNE			HLE
			λ = 0.05	λ = 0.1	λ = 0.2	λ = 0.05	λ = 0.1	λ = 0.2
位置参数	a	0.2	175.24	169.39	165.52	81.19	89.26	95.74
		0.5	108.59	93.49	90.41	33.59	36.31	41.06
		1	54.50	37.97	31.20	15.76	15.63	16.83
尺度参数	b	1.2	120.33	106.75	101.35	199.51	188.86	182.48
		1.5	80.83	62.94	54.98	184.53	171.85	163.73
		2	56.44	38.82	31.09	167.22	153.37	137.51
组分数	c	1.2	177.58	174.46	165.99	209.12	199.89	199.10
		1.5	149.85	139.58	131.20	203.20	201.34	194.02
		2	120.88	103.32	94.07	213.08	199.75	186.10
H			2.102	2.283	2.595	0.165	0.275	0.432

由表8可知，由于HLE图被设计成只对位置变化敏感，HLE图在检测位置参数发生偏移时更有效。但在尺度参数发生偏移时，由于VGNE图善于捕捉高斯混合分布的特征，在检测尺度参数中小偏移方面表现优秀的检测性能，而HLE图的ARL₁却是VGNE控制图的几倍。在组分发生变化时，VGNE图有更低的ARL₁，而HLE图的ARL₁在预设的ARL₀ = 200附近变化，难以检测到组分的变化。因此，VGNE图能够对高斯多模态过程进行有效且稳健的在线监控。

综合上述分析，本文系统评估了所提出的VGNE控制图与基准方法HLE控制图在多个候选平滑参数$\lambda \in \left\{0.05,0.1,0.2\right\}$ 下的性能表现。结果表明，在所有失控场景下，当$\lambda =0.2$ 时，VGNE控制图均取得了最优的RMI值。因此，本文选取$\lambda =0.2$ 作为最优平滑参数，该取值在保持对小且持续偏移的敏感性与快速检测较大幅度偏移的能力之间取得了最佳平衡。后续的实例研究均基于该参数设置展开。

4. 实例分析

半导体是一种具有介于导体和绝缘体之间的电导特性的材料。主要应用在集成电路、通信系统、光伏发电、大功率电源转换等领域。在统计过程控制中，监测半导体制造过程主要通过收集和分析生产过程中的数据来实现。本节给出了一个实际数据示例，以说明VGNE图的有效性。本文选用加州大学欧文分校的机器学习库中维护了一段从2008年7月至10月期间由自动管理半导体制造过程的计算机系统收集的半导体制造过程数据，其原始数据可以在UC Irvine Machine Learning Repository (http://archive.ics.uci.edu/ml/datasets/SECOM)网站上找到。该组半导体制造过程数据共包含1567个观测数据，每个观测数据有590个质量特性，其中前1463个观测数据为合格品，剩余104个为不合格品。从各种属性中选择一个变量，记作“X528”。考虑到数据量大，去除异常值后，我们将前600个数据作为训练样本，后104个不合格数据作为测试样本。Shapiro-Wilk正态性检验结果显示，X528的IC样本的p值为$7.20\times {10}^{-11}$ ，拒绝原假设，表明数据显著偏离高斯分布。

由图1可知，由测试样本得到的经验概率密度图存在两个峰值，这表示由这600个数据作为IC数据集并非都来自同一个分布，而是由两种分布不同的过程混合而成。因此这个半导体制造过程可能是多种操作模式的组合，并且过程数据可能遵循多模态分布。当${\text{ARL}}_0=200$ ，$\lambda =0.2$ 时，绘制VGNE以及HLE控制图。

由图2可知，对于后104个测试样本，VGNE图偏移发生的第一个信号出现在第7次观测，HLE图未发出报警信号，结果显示我们的方案略优于HLE。在实际应用中，当控制图发出第一个信号时，我们认为过程已经失控了，因为在此信号之后所有的统计数据几乎都超过了控制线。此时，操作员应立即检查工艺设置和设备校准情况，努力将工艺恢复到IC状态。如果初步调整未能解决问题，则需要及时联系工艺工程师，以进一步分析和解决潜在的工艺问题。

Figure 1. Empirical PDF of the IC dataset

图1. IC数据集的经验PDF

Figure 2. Monitoring of Phase II observations using VGNE and HLE control charts

图2. VGNE和HLE控制图用于阶段II观测数据的监控

5. 总结

多模态数据在许多复杂应用中都有出现，这类数据不满足单模态和高斯分布的基本假设，因此难以使用某一具体的参数分布对其进行建模。传统的统计过程控制方法在处理多模态数据时可能表现不佳。为解决这一限制，拟议的VGNE图表采用了数据驱动的方法，利用DP自适应地选择模型中的分量数无需预设分量数，该方法有效捕捉了数据固有的多模态性，展示了在高斯和非高斯多模态过程中稳健的监控性能。

最后，尽管我们提出的方法为多模态过程设计得复杂，VGNE图表展现了显著的鲁棒性和广泛的适应性。其对未知分布的过程进行有效监控，而无需事先了解流程知识。目前的主要限制在于当IC观测不足时，在线监测的实施。因此，未来的改进将聚焦于多元扩展和自启动能力开发，以促进更灵活的实时应用。