-
摘要: 多模态是复杂工业生产过程的普遍特性.不同模态具有不同的过程特性,需要建立不同的模型,因此离线建模数据的模态划分与识别是整个多模态过程建模的关键问题之一.目前,常用的聚类算法需要对其结果进行人工分析和后续处理,无法真正实现多模态过程的全自动模态识别.因此,本文提出一种全自动的多模态过程离线模态识别方法.首先通过宽度为H的大切割窗口对数据进行切割,利用改进的K-means聚类算法对窗口单元进行聚类;根据聚类结果,对稳定模态淹没现象进行处理,得到模态的初步划分结果;最终,利用小滑动窗口L,对稳定模态及过渡模态交接区域进行细划分,准确定位稳定模态与过渡模态的分割点.算法实现了多模态过程的全自动离线识别,并给出合理有效的识别结果.仿真分析表明此方法能够实现模态的自动识别,且识别结果准确.Abstract: Multimode is a general characteristic of complex industrial processes. Different modes have different process characteristics, and different models should be established. Therefore, offline mode identification is one of the critical problems for multimode processes modelling. Presently, the commonly used clustering methods cannot realize offline mode identification of multimode processes automatically because human analysis and further processing are needed to gain the final identification result. A fully automatic offline mode identification method is proposed in the paper. First, the data is divided into a series of data segments by a cutting window with the designated width H. The improved K-means clustering method is used to assign the segments into different clusters. According to the clustering result, the missing stable modes are dealt with to obtain the preliminary mode identification results. Finally, the regions between the stable modes and transitional modes are further analyzed by a small moving window L to determine the accurate boundaries between different modes. The mode identification of multimode processes can be realized automatically by the method for a reasonable and effective identification result. Feasibility and practical value of the method are evaluated by case study.
-
Key words:
- Mode identification /
- multi-mode process /
- transition mode /
- stable mode /
- fully automatic
-
随着社会的进步,工业生产更加倚重于可以生产多种产品的高效过程,过程生产方案的变动或者产品类型的改变会导致生产过程出现具有不同过程特性的多种模态,多模态也成为复杂工业生产过程中普遍存在的一种特性.因此,多模态过程成为研究热点之一.
多模态生产过程包含多个稳定生产模态,不同的稳定模态之间存在不同的过渡模态[1].不同模态的过程特性差别很大,变量相关性也有着显著的差异,故无法用一个模型对多个不同的模态特性同时进行描述.为了解决多模态特性对过程监测、故障诊断等问题的影响,Yu等[2-3] 提出了基于高斯混合模型(Gaussian mixture model,GMM)的多模态处理方法.在此基础上,Xie等[4]提出了滑动GMM的方法,解决多模态过程中单个模态的动态问题. Ge等[5-6]利用混合概率PCR的思想,巧妙解决了多模态过程的质量回归问题.这些基于贝叶斯思想的方法将所有模态放在一起建立一个大的混合模型,无需进行模态识别.但此类方法只考虑了稳定模态,并没有考虑过渡模态的影响.如果将这些方法直接应用到含有过渡模态的过程中,则会导致过渡模态描述信息不完全,在监测的过程中可能给出错误的结论.但是过渡模态是多模态过程的一大特征,在许多过程中是无法忽略的.因此,对于多模态过程,尤其含有过渡模态的多模态过程,无论是过程优化、过程监测、故障诊断或者过程状态评价等等,都需要对差异明显的多个模态采用不同的方法分别进行分析.
为了解决该问题,Ng等[7]提出了一种基于重叠PCA的联合监测方法,但是该方法无法分辨稳定模态与过渡模态.Lu等[8]和Zhao等[9]提出了多模型的监测思想并将基于该思想的一系列监测方法[10-13]成功应用到了多阶段批次过程和多模态连续过程中.Zhang等[14-15]有效提取所有模态之间的公共信息,将每个模态分成公共子空间和独立子空间,提高了监测效果.这些多模型方法能够有效解决过渡模态所带来的影响,但是在应用时首先需要对离线数据进行模态识别和划分,然后才能针对不同模态建立不同模型.也就是说,离线建模数据的模态划分与识别是实现多模态过程建模、优化、监测以及状态评价的基础和首要前提,是研究多模态过程的关键问题之一.
离线建模数据的模态划分与识别就是将一组正常操作数据按照其生产模态的不同进行类别划分,同时辨别出生产数据对应的模态信息(某一稳定模态或某过渡模态)[1],对于简单的有模态指示变量的生产过程,可以根据已知的模态指示变量准确地判断出建模数据所对应的运行模态.但大多工业生产过程相对复杂,通常没有模态指示变量,因此需要提出一种合理有效的方法来实现过程模态的自动区分与识别.
目前,已经有许多聚类算法应用到模态识别中,例如K-means聚类算法[16]、减聚类算法[17]、模糊C均值聚类算法[18-19]等.但这些算法在处理多模态过程时存在以下问题:1) 一个较完整的多模态过程的操作数据,尤其一些连续生产的大工业过程,包含了多个稳定模态以及稳定模态间的过渡模态,样本容量大、复杂度高.如果数据不加处理,直接应用聚类算法,则计算量大,且容易给出错误的聚类结果.2) 这些聚类算法只能将数据分成不同的类,无法给出具体的模态信息.尤其是当多模态过程含有过渡模态时,聚成的子类属于哪个模态,这些子类之间有什么关系,这些信息都无法根据聚类结果直接获取.3) 过程数据经常会出现一些奇异点,这些奇异点对模态的划分有着很大的影响,因此需要对奇异点现象进行处理.4) 由于过渡模态具有一定的动特性,多模态过程数据复杂,聚类算法的参数选取比较困难.不准确的参数有可能导致某一稳定模态被划分到不同子类中,从而出现稳定模态被淹没的现象.因此,这些聚类算法只适用于多模态过程的初始分类,要得到具体的模态信息,需要对其聚类结果进行进一步分析.
目前,大多数实际过程都是对聚类结果进行简单人为处理,而且没有考虑到初始聚类结果不理想的情况,因此当前的一些处理算法存在两大主要问题:1) 不能实现完全自动的模态识别; 2) 给出错误识别结果的可能性很大.
为了解决上述问题,本文提出一种系统的全自动离线模态识别方法.首先,选取长度为 $H$ 的切割窗口对离线建模数据进行切割,将各窗口作为分析的基本单元,提取每个窗口内数据的均值,通过改进的K-means聚类算法[16],对窗口均值向量进行聚类;针对由于聚类算法参数不准确等原因造成的稳定模态被埋没的现象进行处理,实现稳定模态和过渡模态的初步划分;在初步确定的模态基础上,选定一个小滑动窗口 $L$ ,对稳定模态及过渡模态交接区域进行进一步细划分,准确定位稳定模态与过渡模态的分割点,实现多模态建模数据的模态识别.
1. 多模态过程的离线模态识别
多模态过程中,稳定模态是生产过程运行状态比较稳定的过程,变量基本在某一种状态下稳定运行,因此时间窗口内的变量均值在稳定模态期间变化不大.相比于稳定模态,过渡模态是从一种稳定运行模式向另一种稳定运行模式渐近式量变的动态过程,出现了较多的动态变量,变量的均值随时间方向也呈现出一种动态的渐变趋势.根据多模态过程的这一特性,本文对多模态过程进行离线模态识别.
假设一个多模态生产过程,其建模数据表示为最常用的二维矩阵形式 $\hat{X}\in{\bf R}^{N× J}$ ,其中 $N$ 和 $J$ 分别表示过程的采样个数以及过程变量个数.首先对数据进行无量纲化处理[20],即将不同变量的方差归一,从而使每一个变量在数据模型中都具有同等的权重,消除了过程变量测量值的量程差异对聚类结果的影响,无量纲处理后的数据记为 $X$ .
1.1 模态初步确定
首先,寻找过程数据中的稳定模态和过渡模态,实现过程模态的初步确定.
步骤 1. 窗口切割
选取长度为 $H$ 的切割窗口对二维数据矩阵 $X$ 沿着采样方向进行分割,即将矩阵 $X$ 转置后沿横轴方向分割.
设 $N$ 可表示为 $N=K× H+d$ $(0\leq d＜H)$ .这样 $N$ 个采样数据被分割为 $K$ 个窗口,每个切割窗口按着采样时间方向顺次排列.第 $k$ , $k=1,2$ , $\cdots$ , $K$ 个窗口的数据表示为 $X^{\rm T}_k\in {\bf R}^{J× H}$ .求取这些二维矩阵的均值向量,记为 $\overline{x}_k$ ,这些均值向量表征了窗口内变量的相关信息,成为模态划分的基本单元.
步骤 2. 聚类分析
下面利用改进的K-means聚类算法对均值向量进行聚类.该聚类算法的具体步骤简述如下:
算法的输入是均值向量集合 ${\overline{x}_1,\overline{x}_2,\cdots,\overline{x}_K}$ ,以及两个子类中心的最小距离阈值 $\theta$ .算法的输出是子类数量 $C_st$ .子类中心 ${w_1,w_2,\cdots,w_{C_{st}}}$ ,以及每个窗口属于不同子类的隶属关系 $u(k):\overline{x}_k\to$ $1,2$ , $\cdots$ , $C_{st}$ .变量 $i$ , $k$ 及 $c$ 分别是算法中迭代次数、分类模式以及聚类中心的索引.
1) ~从 $K$ 个被分类单元中,任意选择 $C_0$ 个单元向量作为初始聚类中心 $W_{i,c}$ $(c=1,2,\cdots,C_0) $ .如果 $C_0$ 个初始聚类中心选取的过于集中或者个数较少,有可能造成局部最优,导致最终形成的子类比真正的子类要少,给出错误结果.因此,一般常用方法是从被分类单元向量中均匀抽取 $C_0$ 个模式,建议 $C_0$ $=$ $K/3$ $\sim$ $K/2$ .
2) ~若两个子类中心的距离 $dist(w_{i,c_1},w_{i,c_2})$ 小于预定的阈值 $\theta$ ,剔除其中一个聚类中心,更新子类数量,记为 $C_{i+1}$ .
3) ~计算每个单元向量 $\overline{x}_k$ $(k=1,2,\cdots,K)$ 到所有聚类中心的距离,若 $\overline{x}_k$ 和第 $c^*$ 类的中心的距离最小,则将 $\overline{x}_k$ 的隶属关系定义为 $u(k)=c^*$ .
4) ~根据单元向量的隶属关系重新计算新的聚类中心 $dist(w_{i,c}$ , $c=1,2,\cdots,C_{st})$ .
5) ~如果算法满足收敛条件则结束,否则返回2) ,进行下一次迭代计算.收敛条件有:两次迭代中的聚类中心距离的变化小于一个很小的阈值 $\epsilon$ ，或者每个子类中各单元向量 $\overline{x}_k$ 到子类中心的距离平方和以及子类之间的距离平方和达到最小.
上述分类算法最终可以将整个多模态过程分成 $C_{st}$ 子类.
步骤 3. 时段划分
将窗口数据沿着时间轴展开,时间上连续且属于同一个子类的窗口单元划分到同一个时段.
假设整个过程被初步划分为 $M_0$ 个时段,各个时段按照时间方向有序排列.第 $m$ 个子时段记为 $\hat{X}_m$ , $m=1,2,\cdots,M_0$ . $\hat{X}_m$ 属于第 $c$ 个子类,即 $u(\hat{X}_m)$ $=$ $c$ ,其中 $c=1,2,\cdots,C_{st}$ .第 $m$ 个子时段时间长度记为 $\hat{S}_m$ .
由于聚类算法的输入是按照时间顺序排列的窗口均值向量,因此按照时间顺序,将时间上相邻且属于同一子类的窗口单元划分到同一个时段中,这样就使得那些具有相似过程状态的连续时间窗口被包含在同一个时段中,时段内部也呈现出时间的连续性.也就是说,如果过程变量在某段时间内保持一致,这段时间上的均值向量将会被聚类算法划分在同一个时段中.当过程变量的状态发生改变时,对应的均值窗口被分到不同的时段.
步骤 4. 稳定模态确定
根据经验知识定义稳定模态最短运行时间 $S_{\min}$ ,将时段长度大于稳定模态最短运行时间 $S_{\min}$ 的时段定义为稳定模态,且隶属于同一子类的稳定模态定义为同一种稳定模态.
如果两个相邻的稳定模态之间所有时段的运行长度之和大于 $S_{\min}$ ,说明在两个稳定模态之间仍然含有稳定模态.这是由于多模态过程既含有稳定模态,又含有动态特性明显的过渡模态,数据特性复杂.聚类算法在聚类时由于参数设定不恰当或者陷入局部最优等问题,可能使得某稳定模态被聚成比较多的子时段,导致该模态在初步识别时无法满足稳定模态的识别条件,造成稳定模态识别上的缺失,这里将这种现象定义为稳定模态的淹没现象.我们需要对该现象进行分析处理.
假设现在已确定有 $A$ 和 $B$ 两种稳定模态,其中稳定模态 $A$ 和 $B$ 相邻,且 $A$ 和 $B$ 之间所有时段的长度之和大于 $S_{\min}$ ,则做以下处理:
1) ~在稳定模态 $A$ 和模态 $B$ 之间的所有子时段中寻找运行时间最长的一段,作为中间新稳定模态的基准时段.如果运行时间最长的时段不只一个,则选择最靠近中间位置的时段.这是由于两个稳定模态之间一定含有过渡模态,且稳定模态运行时间远远大于过渡模态运行时间,所以若稳定模态 $A$ 和稳定模态 $B$ 之间含有稳定模态,则该模态应该位于靠近模态 $A$ 和 $B$ 之间的中间位置.假设第 $k$ 个时段是该稳定模态的基准时段,记为 $\hat{X}_{\rm base}=\hat{X}_k$ ,其长度 $\hat{S}_{\rm base} =\hat{S}_k$ .与 $\hat{X}_{\rm base}$ 相邻的两个时段分别定义为 $\hat{X}_{k_1}$ 和 $\hat{X}_{k_2}$ ,其中 $k_1=k-1$ , $k_2=k+1$ .
2) ~以 $\hat{X}_{\rm base}$ 为基准,比较其左右两个相邻的时段与 $\hat{X}_{\rm base}$ 的相似度,即计算 $\gamma(\hat{X}_{\rm base},\hat{X}_{k_1})$ 以及 $\gamma(\hat{X}_{\rm base},\hat{X}_{k_2})$ .时段 $\hat{X}_1$ 和 $\hat{X}_2$ 的相似度定义为
\begin{align}\gamma_{1,2}=\gamma_{\hat{X}_1,\hat{X}_2}=\exp\left(-\frac{\Arrowvert{\overline{X}_1-\overline{X}_2}\Arrowvert}{2J}\right)\end{align}
(1) 其中, $\overline{X}_1$ 和 $\overline{X}_2$ 分别为时段 $\hat{X}_1$ 和 $\hat{X}_2$ 的均值向量.
3) ~如果两个时段相似度相等,则将两个时段同时加入到 $\hat{X}_{\rm base}$ 中,否则选择相似度高的时间段加入到 $\hat{X}_{\rm base}$ 中,组成新的 $\hat{X}_{\rm base}$ .
\begin{align}&\hat{X}_{\rm base}=\begin{cases} [\hat{X}_{\rm base}\ \ \hat{X}_{k_2}],& \gamma_{{\rm base},k_1}＜\gamma_{{\rm base},k_2}\\ [\hat{X}_{k_1}\ \ \hat{X}_{\rm base}\ \ \hat{X}_{k_2}],& \gamma_{{\rm base},k_1}=\gamma_{{\rm base},k_2}\\ [\hat{X}_{k_1}\ \ \hat{X}_{\rm base}],& \gamma_{{\rm base},k_1}>\gamma_{{\rm base},k_2}\end{cases}\end{align}
(2) \begin{align} &\hat{S}_{\rm base}=\begin{cases} \hat{S}_{\rm base}+\hat{S}_{k_2},& \gamma_{{\rm base},k_1}＜\gamma_{{\rm base},k_2}\\ \hat{S}_{k_1}+\hat{S}_{\rm base}+\hat{S}_{k_2},& \gamma_{{\rm base},k_1}=\gamma_{{\rm base},k_2}\\ \hat{S}_{k_1}+\hat{S}_{\rm base},& \gamma_{{\rm base},k_1}>\gamma_{{\rm base},k_2}\end{cases}\end{align}
(3) 4) ~如果时段满足终止条件则认为稳定模态已经确定,否则返回2) ,进行下一次迭代计算.终止条件为:迭代中最高相似度小于阈值 $\sigma$ ,而且 $\hat{X}_{\rm base}$ 的时间长度 $\hat{S}_{\rm base}$ 大于稳定模态最短运行时间 $S_{\min}$ .
5) ~计算新确定稳定模态与所有已知的稳定模态的中心距离 $d^i_{\rm new}=\Arrowvert{\overline{X}_{\rm new}-\overline{X}_i}\Arrowvert$ ,其中 $i=A,B$ , $\overline{X}_{\rm new}$ 和 $\overline{X}_{i}$ 分别为新确定稳定模态和第 $i$ 个已知稳定模态的样本均值.如果最小的距离小于预定的阈值 $\theta$ ,则将该稳定模态可定义为相应的稳定模态 (例如, $d^A_{\rm new}=\min(i^A_{\rm new})$ 且小于阈值 $\theta$ ,则将新确定稳定模态定义为 $A$ 稳定模态),否则说明该稳定模态是没有出现过的新模态,定义为 $C$ 稳定模态.
6) ~计算新确定稳定模态与相邻的前一稳定模态 $A$ 之间所有时段的运行长度 $S_{A\rm new}$ .判断 $S_{A\rm new}$ 是否大于 $S_{\min}$ .如果 $S_{A\rm new}$ $＜$ $S_{\min}$ ,认为新确定模态与模态 $A$ 之间不存在其他稳定模态.否则,新稳定模态与模态 $A$ 之间仍然存在稳定模态,按照新稳定模态的确定方法,继续寻找稳定模态,即重复1) $\sim$ 5) 所描述的提取过程,直到找出新确定模态与模态 $A$ 之间所有的稳定模态.
7) ~对新确定稳定模态与相邻稳定模态 $B$ 之间的时段进行处理,处理方法与6) 类似,直到找出新确定模态与模态 $B$ 之间所有的稳定模态.
除了考虑两个相邻稳定模态之间是否存在新稳定模态,仍有两种情况需要考虑:1) 过程起始状态为非稳定模态时,分析过程起点即第1个时刻到第一个稳定模态起点之间的数据是否含有新稳定模态; 2) 过程终止状态为非稳定模态时,分析最后一个稳定模态的终止时刻到过程终点之间是否含有新稳定模态.分析过程与两个稳定模态( $A$ 和 $B$ 为例)之间的分析过程类似.
经过上述处理,多模态过程的所有稳定模态都已确定,而且稳定模态之间的长度均小于稳定模态最短运行时间 $S_{\min}$ .
步骤 5. 处理两个相邻稳定模态之间的时段
对于两个相邻的稳定模态,有以下两种情况: 1) 相邻的两个稳定模态分别为两种不同的模态 $A$ 和 $B$ ,则认为两个稳定模态之间的所有子时段属于稳定模态 $A$ 到模态 $B$ 的过渡模态. 2) 相邻两个稳定模态为同一种模态 $A$ ,则将两个稳定模态之间的所有子时段定义为稳定模态 $A$ 的扰动过程.这种现象的发生是由于实际过程在正常运行的时候,出现较长时间的扰动,使得一些变量发生变化,经过控制系统的调节,过程再次达到稳定时的新状态与扰动前的状态相同,因此扰动前后的两种模态可以看作同一种模态.
通过步骤1 $\sim$ 5,算法实现了多模态过程中所有的稳定模态以及过渡模态的识别与初步划分.模态的初步划分过程示意图如图 1所示.
1.2 模态准确划分
模态划分的关键是确定当前模态的结束时间,即确定下一种模态的开始时间.上一节实现了多模态过程各模态的初步划分,但是聚类单元为时间长度一定的切割窗口,因此模态的初步划分只能确定模态开始转变的起始窗口,却无法得知窗口内模态开始发生变化的具体时刻,因此需要对稳定模态与过渡模态相邻的窗口进行更深层次的分析.
如图 2所示,假设从第1个窗口到第 $k_1-1$ 个窗口为稳定模态 $A$ ,从 $k_1$ 个窗口开始进入过渡模态.此时可能有两种情况:1) 过渡模态的开始时间发生在第 $k_1$ 个窗口的前半段.那么第 $k_1-1$ 个窗口完全属于稳定模态,第 $k_1$ 个窗口数据以过渡模态为主,其均值向量到稳定模态 $A$ 中心距离偏大,因此没有聚到稳定模态 $A$ 所在的类中. 2) 过渡模态的开始时间发生在第 $k_1-1$ 个窗口的后半段.稳定模态的过程数据占据第 $k_1-1$ 个窗口的大部分,因此其均值向量离稳定模态 $A$ 所在的类中心距离最近,故将此窗口归到稳定模态中.第 $k_1$ 个窗口内数据完全属于过渡模态,到稳定模态类中心的距离很大,因此将其分到了过渡模态起始子时段中.
图 2 多模态准确识别过程(以 $AC$ 过渡模态起始时刻为例)Fig. 2 Exact mode identification for multimode processes (i.e.beginning of the transitional mode $AC$ )基于以上分析,为了准确判断稳定模态的结束时间,需要从第 $k_1-1$ 个窗口开始用较短的滑动窗口 $L$ 重新分析判断,滑动步长为 $h$ .如图 2所示,通过比较每个滑动窗口与稳定模态 $A$ 的相似度来确定过渡模态的起始时间,定义滑动窗口与稳定模态的相似度,如式(4) 所示:
\begin{align}\gamma_{\widetilde{x}_t,\overline{X}_A}=\exp\left(-\frac{\Arrowvert{\widetilde{x}_t-\overline{X}_A}\Arrowvert}{2J}\right)\end{align}
(4) 其中, $\widetilde{x}_t$ 为第 $t$ 个小滑动窗口的均值, $\overline{X}_A$ 为稳定模态 $A$ 的变量均值.
相似度定量指示了各窗口与稳定模态的相似度.当多模态过程从稳定模态逐渐进入过渡状态,其变量越来越偏离前一稳定模态所在的稳定状态,相应地,与稳定模态的相似度也逐渐减少.本文引入可以调节的相似度阈值 $\alpha$ 作为边界参数,按照如下的准则分析模态识别问题:如果 $t_1$ 是第一个满足以下条件的窗口:从第 $t_1$ 个小窗口开始,连续有 $r$ 个窗口都满足 $\gamma_t＜\alpha$ ,即满足时,认为过程从 $t_1$ 个小窗口进入到过渡模态.将第 $t_1$ 个小窗口的起始位置作为过渡模态的开始,则此过渡模态起始时刻为 $(k_1$ $-$ $2) × H+{(t_1-1) × h}+1$ .
\begin{align}\begin{cases} \gamma_{t_1-1}\ge\alpha,&\\ \gamma_t＜\alpha,&t=t_1,t_1+1,\cdots,t_1+r-1\end{cases}\end{align}
(5) 过渡模态终止时间的确定与起始时刻的确定类似,不同的是考虑过渡模态与下一个稳定模态 $C$ 的相似度,而且当多模态过程从过渡模态逐渐进入稳定模态时,滑动窗口与稳定模态 $C$ 的相似度逐渐增大. $AC$ 过渡模态到第 $k_2-1$ 个窗口结束,从第 $k_2$ 个窗口进入下一稳定模态 $C$ .此时需要对从第 $k_2-1$ 个窗口开始用较短的滑动窗口 $L$ 进行分析.通过比较每个滑动窗口与稳定模态 $C$ 的相似度来确定过渡模态的终止时间.如果 $t_2$ 是第一个满足以下条件的窗口:从第 $t_2$ 个小窗口开始,连续有 $r$ 个窗口都满足 $\gamma_t\ge\alpha$ ,认为过程从 $t_2$ 个小窗口进入到稳定模态.将第 $t_2$ 个小窗口的起始位置记为下一稳定模态的开始,则过渡模态终止时刻为 $(k_2-2) × H +(t_2-1) $ $×$ $h+1$ .
\begin{align}\begin{cases} \gamma_{t_2-1}＜\alpha,\\ \gamma_{t}\ge\alpha,&t=t_2,t_2+1,\cdots,t_2+r-1\end{cases}\end{align}
(6) 由于规定窗口长度 $H$ 不得大于最短过渡时间,因此,一般不会出现过渡模态被淹没的现象.即使出现这一情况,也不会影响最终的识别结果.例如,由于窗口长度 $H$ 取值过大,导致初步识别时某一过渡模态被淹没,即两个稳定模态直接相连,中间没有过渡模态.这种情况下,说明被淹没的过渡模态一定位于两个稳定模态相连接的两个 $H$ 窗口内.根据本节提出的准确划分方法,对这两个相邻的 $H$ 窗口进行分析处理,便可确定过渡模态的起始时间和终止时间.
1.3 重要参数分析
在模态识别过程中,涉及到几个重要参数的选取,影响模态的辨识结果.
1) 大切割窗口 $H$
一个比较完整的多模态过程操作数据,包含了不同的稳定模态以及稳定模态间具有动特性的过渡模态,样本容量大,复杂度高.如果不经过任何处理,直接对其进行聚类,奇异点等异常数据的存在以及过渡模态的动特性会大大影响聚类效果.因此,大切割窗口的设定一方面可以快速处理大量的建模数据,降低复杂度,另一方面可以有效减少异常数据对聚类的影响.
$H$ 的选取会直接影响模态的初步确定.首先, $H$ 不得大于最短过渡时间长度,否则可能造成过渡模态的淹没现象. $H$ 越大,对噪声、奇异点的冗余能力越强.但同时,各模态起始终止时间的初步确定范围越宽. $H$ 越小,各模态起始终止时间的初步确定更加准确,但噪声奇异点等随机扰动的影响越大,聚类复杂度越高.在实际应用中,可以根据实际过程选取恰当的切割窗口.例如,在实际过程中,离线数据的样本数很大,数据质量较差,可选取较大的 $H$ .如果离线数据样本数较少,或者过渡时间较短,便选取较小的 $H$ .
2) 小滑动窗口 $L$
为了更加精确地确定各模态起始终止时间,需要对稳定模态与过渡模态相交接的两个窗口进行分析,此段数据长度为 $2H$ .因此需要一个小的滑动窗口 $L$ 来进行分析. $L$ 越大,去噪能力越强,但过渡模态起始、终止时刻的确定误差有可能越大. $L$ 越小,噪声或奇异点等随机扰动导致错误结果的可能性越大.
滑动步长 $h$ 一般不能取得太大,否则会造成漏点现象. $h$ 越大,滑动速度越快,结果误差可能越大. $h$ 越小,滑动速度越慢,但结果越精确.在小滑动窗口 $L$ 和滑动步长 $h$ 的实际选取过程中,需要根据实际情况,选取合适的窗口长度,保证识别的有效准确.例如, $H$ 较小时, $h$ 完全可以取一个较小的值,甚至取为1. $H$ 很大时,考虑到算法的运算速度, $h$ 可以取一个较大的值,但是如果过程对识别精度要求非常高, $h$ 取值也不能太大.
3) 相似度阈值 $\alpha$
在模态识别的过程中,阈值 $\alpha$ 的选取会直接影响到模态识别精度和监测效果. $\alpha$ 越大,要求窗口数据与稳定模态相似度越高,稳定模态的数据离稳定模态中心越近,但过高的阈值会导致原本属于稳定模态的数据也划分到过渡模态中,在增大过渡模态范围的同时丢失了部分稳定模态信息. $\alpha$ 太小,会将过渡模态数据划分到稳定模态中,从而无法准确定位模态的起始终止位置.
因此,针对不同的模态,结合过程数据提供的信息,定义一个弹性相似度阈值 $\alpha$ .以图 2为例,在初步划分时,从第一个窗口到第 $k_1-1$ 个窗口为稳定模态 $A$ ,从 $k_1$ 个窗口开始进入过渡模态. $AB$ 过渡模态到第 $k_2-1$ 个窗口结束,从第 $k_2$ 个窗口进入下一稳定模态 $B$ .
经过分析知道,过渡模态确切起始位置应该在第 $k_1-1$ 和第 $k_1$ 这两个窗口的某个时刻,而第 $k_1-2$ 以及往前几个窗口一定属于稳定模态 $A$ .如图 2所示,第 $k_1-2$ 往前的 $q$ 个大的分割窗口数据记为 $\widehat{X}$ ,即
\begin{align}\widehat{X}=[X_{k_1-1-q} \cdots X_{k_1-2}]\end{align}
(7) 用小滑动窗口 $L$ 对 $\widehat{X}$ 进行分析,计算各个小滑动窗口与稳定模态 $A$ 的相似度 $\gamma$ .利用核密度估计[21],计算上述相似度的密度分布,这个分布给出了当前数据属于稳定模态 $A$ 时相似度的概率分布情况.
在一个事先给定的置信水平 $1-\beta$ 下,计算相似度 $\gamma$ 的控制限,定义为阈值 $\alpha_1$ .此阈值可以作为分析第 $k_1-1$ 和第 $k_1$ 这两个窗口的边界参数.
确定过渡模态结束时间时,同上述方法类似.分析 $k_2+1$ 往后 $q$ 个大的分割窗口数据与稳定模态 $B$ 的相似度,求取控制限,定义为阈值 $\alpha_2$ ,作为分析第 $k_2-1$ 和第 $k_2$ 两个窗口的边界参数.
2. 田纳西-伊斯曼过程的仿真研究
2.1 TE 过程描述及实验设定
田纳西-伊斯曼(Tennessee Eastman,TE)过程是由伊斯曼化学品公司创建的,该过程是基于一个真实工业过程的仿真,作为比较各种方法的数据源,已在故障监测与诊断方面得到了广泛的应用[22].
TE过程包括五个主要单元[23]: 反应器、冷凝器、压缩机、分离器和汽提塔,共有41个测量变量和12个控制变量. 研究表明,TE过程的反应器压力和反应器液位设定值大小对生产成本有较大影响.TE过程的操作点可以根据生产要求进行调整,这种情况下TE过程就是一个典型的多模式过程.本文以反应器压力和液位设定值变化这两种操作模式改变的过程为对象,设计三种稳定模态和一种扰动,如表 1所示.设定值的改变会对很多过程变量产生较大影响,共选择15个过程连续变量作为被监控变量,并将其列于表 2.
表 1 TE过程操作模态列表Table 1 Operating modes of Tennessee Eastman process模态类型 反应器压力设定值(kPa) 反应器液位(%) 稳定模态 $A_0$ 2800 65 稳定模态 $B_0$ 2600 65 稳定模态 $C_0$ 2400 65 扰动1 2800 75 表 2 TE过程变量表Table 2 Variables of Tennessee Eastman process序号 变量名称 1 $A$ 进料量(流1) 2 $D$ 进料量(流2) 3 $E$ 进料量(流3) 4 $A,C$ 混合物料流量 5 再循环流量(流8) 6 反应器进料速度(流6) 7 反应器温度 8 排放速度(流9) 9 产品分离器温度 10 产品分离器压力 11 产品分离器塔底流量(流10) 12 汽提塔压力 13 汽提塔温度 14 反应器冷却水出口温度 15 分离器冷却水出口温度 模拟180个小时的仿真过程,每小时采样100次.稳定模态最短运行时间为15小时,即最小稳定模态长度 $S_{\min}$ 为15h.改变两次设定值,得到三种稳定模态的运行数据.反应器压力初始设定值为2800kPa,过程为 $A_0$ 稳定模态.当运行至 $T=30$ h时,将反应器压力设定值降低为2600kPa,过程进入过渡模态 $A_0B_0$ .在 $T=35.20$ h时,过渡模态结束,过程开始进入下一稳定模态 $B_0$ .过程运行至 $T$ $=$ $70$ h时,改变反应器压力设定值2400kPa,过程状态再次发生变化,此时进入过渡模态 $B_0C_0$ .在 $T$ $=$ $80.28$ h时,过渡模态结束,过程进入下一稳定模态 $C_0$ . $T=110$ h时,将反应器压力重新设定为2800kPa,此时过程从稳定模态 $B_0$ 向稳定模态 $A_0$ 过渡.在 $T=122.36$ h时,过程再次进入稳定模态 $A_0$ .在 $T=140$ h时,维持反应器压力不变,对过程加了一个扰动,即将反应器液位设定值提高 $10%$ ,经过调节,在 $T=147.21$ h时重新回到稳定模态 $A_0$ 运行.设计的实验如表 3所示,仿真共得到18000个采样数据,包括三种稳定模态和三种过渡模态,以及一次扰动.1 $\sim$ 3000个样本点属于模态 $A_0-1$ ,从3001个样本点进入到 $A_0B_0$ 过渡模态,其中 $A_0-1$ 表示过程第一次进入稳定模态 $A_0$ ,其他数据的含义一样.图 3显示了15个过程变量的变化曲线,其中横轴为采样点,纵轴为各个变量的值.
表 3 TE过程实验设计Table 3 Experimental design of TE process模态 聚类单元 $A_0-1$ 1 $\sim$ 3000 $A_0B_0$ 3001 $\sim$ 3520 $B_0$ 3521 $\sim$ 7000 $B_0C_0$ 7001 $\sim$ 8028 $C_0$ 8029 $\sim$ 11000 $C_0A_0$ 11001 $\sim$ 12235 $A_0-2$ 12236 $\sim$ 14000 扰动 14001 $\sim$ 14720 $A_0-3$ 14721 $\sim$ 18000 
图 3 TE过程15个变量变化曲线 ( $X$ 轴表示采样点)Fig. 3 15 variables of the TE process ( $X$ axis indicates sample points.)2.2 TE 过程初步识别
在离线识别模态前,首先对数据进行去量纲处理,然后利用分割窗口对数据进行切割.本文在第2节中提到,分割窗口 $H$ 的不同取值会影响最终聚类结果.这里针对三种情况进行讨论分析:1) ~ $H=1$ ,即不设置窗口,直接对数据进行聚类; 2) ~ $H=50$ ,即窗口宽度选择合适的情况;3) ~ $H=500$ ,即窗口宽度过大的情况.两个子类中心的最小距离阈值设为 $\theta$ $=$ $2.5$ ,阈值 $\varepsilon=2.5$ .
当 $H=1$ 时,即利用改进的K-means聚类算法对数据直接进行聚类,聚类单元共18000个.初始聚类中心 $C_0$ 为6000个.聚类结果如图 4所示.将时间相邻隶属类别相同的窗口划分到一个时段中,计算各时段的运行时间.
从图 4可以看出,只有三个时段运行时间大于稳定模态最短运行时间,而且属于同一稳定模态,定义为稳定模态 $A$ ,与仿真设定对比可知,稳定模态 $A$ 对应预设模态中的稳定模态 $A_0$ .而设定中其余两个稳定模态由于出现稳定模态淹没现象而没被识别出来.从聚类结果可以看出,如果不对聚类结果进行处理分析,会给出错误结论.而目前并没有发现相关算法对不理想的聚类结果进行讨论分析,也再次证明了本文所提算法的重要性.
图 5表示前后两次迭代聚类中心的距离的变化( $H=1$ ),在第52次迭代时,此值为0.05674,小于阈值 $\varepsilon$ ,算法结束.迭代过程中聚类中心个数的变化如图 6所示,最终确定的聚类中心个数为20个.
图 5 前后两次迭代聚类中心的距离( $H=1$ )Fig. 5 Distance between clustering centers in two iterations ( $H=1$ )当 $H=50$ 时,采样点被分割为360个窗口,每个窗口求取均值,并利用改进K-means算法对窗口均值向量进行聚类.聚类单元共360个,初始聚类中心为120个.聚类结果如图 7所示.
算法初步划分过程将TE过程分成了8个类,20个子时段,每个时段的运行时间以及所属类别如表 4所示.
表 4 各时段运行时间( $H=50$ )Table 4 Runtime of each period ( $H=50$ )时段 窗口个数 运行时间 (h) 所属类别 1 60 30 1 2 3 1.5 3 3 7 3.5 6 4 70 35 4 5 3 1.5 2 6 1 0.5 5 7 7 3.5 6 8 9 4.5 8 9 11 5.5 5 10 19 9.5 7 11 10 5 5 12 6 3 7 13 9 4.5 5 14 11 5.5 7 15 7 3.5 3 16 10 5 4 17 37 18.5 1 18 3 1.5 3 19 12 6 2 20 65 32.5 1 从表 4可以发现,共有4个时段运行时间超过了稳定模态最短运行时间,分别为第1个、第4个、第17个和第20个时段.根据定义,将这4个时段定义为稳定模态.其中第1个、第17个和第20个时段都属于第1个子类,因此将其定义为稳定模态 $A$ .第4个时段隶属于第4个子类,因此将其定义为稳定模态 $B$ .
稳定模态 $A$ 与稳定模态 $B$ 之间所有时段的运行时间之和为 $5$ h,小于稳定模态最短运行时间,因此认为稳定模态 $A$ 与稳定模态 $B$ 之间不存在新的稳定模态.稳定模态 $B$ 与其相邻的后一稳定模态 $A$ 之间的所有时段长度之和为38h,远远大于稳定模态最短运行长度 $S_{\min}$ ,因此判断这两个稳定模态之间仍含有稳定模态,即存在稳定模态淹没现象.
图 8表示前后两次迭代聚类中心的距离的变化( $H$ $=$ $50$ ),在第4次迭代时,此值为0.0697,小于阈值 $\varepsilon$ ,算法结束.迭代过程中聚类中心个数的变化如图 9所示,最终确定的聚类中心个数为8个.
图 8 前后两次迭代聚类中心的距离( $H=50$ )Fig. 8 Distance between clustering centers in two iterations ( $H=50$ )当 $H=500$ 时,采样点被分割为36个窗口,每个窗口求取均值利用对窗口均值向量进行聚类.聚类单元只有36个,初始聚类中心 $C_0$ 为12个.聚类结果如图 10所示.从图 10中可以看到,算法将过程分成了3个类,4个子时段,其中3个时段运行时间大于稳定模态最短运行时间,因此分别定义为稳定模态 $A$ 、 $B$ 、 $C$ .而稳定模态 $C$ 与稳定模态 $A$ 之间有一个单元属于过渡模态.而 $AB$ 过渡模态, $BC$ 过渡模态及扰动被淹没,这是由于 $H$ 取值过大的缘故.
图 10 K-means算法聚类结果 ( $H=500$ )Fig. 10 Clustering result of K-means algorithm ( $H=500$ )图 11表示前后两次迭代聚类中心的距离的变化( $H=500$ ),在第3次迭代时,此值为0,小于阈值 $\varepsilon$ ,算法结束.迭代过程中聚类中心个数的变化如图 12所示,最终确定的聚类中心个数为3个.
图 11 前后两次迭代聚类中心的距离 ( $H=500$ )Fig. 11 Distance between clustering centers in two iterations( $H=500$ )将 $H$ 取不同值时的聚类结果进行比较,如表 5所示.当 $H=1$ 时,即直接对数据进行聚类时,迭代次数为52次,运行时间为519.2秒.而用切割窗口的均值作为聚类单元时,运行时间大大较少,基本都在1秒左右,而且在很少迭代次数内就能达到最优.
表 5 聚类结果比较(不同 $H$ 值)Table 5 Comparison of clustering results (different $H$ )$H$ 值 聚类单元 初始类别数量 最终类别数量 迭代次数 运行时间(s) 1 18000 6000 20 52 519.2 50 360 120 8 4 0.906 500 36 12 3 3 0.745 这是由于窗口的设定能够成倍减少聚类单元的个数,大大降低了算法的复杂度,另外可以有效去除噪声等对数据的影响.图 13～15展示了三种情况下,初始聚类中心与第一个聚类中心的距离变化趋势.
图 13 初始聚类中心到第一个聚类中心的距离( $H=1$ )Fig. 13 Distance between original clustering center and the first clustering center ( $H=1$ )
图 14 初始聚类中心到第一个聚类中心的距离( $H=50$ )Fig. 14 Distance between original clustering center and the first clustering center ( $H=50$ )
图 15 初始聚类中心到第一个聚类中心的距离( $H=500$ )Fig. 15 Distance between original clustering center and the first clustering center ( $H=500$ )在 $H=1$ 时,前1000个初始聚类中心实际上是从第一个稳定模态中抽取出来的.我们希望通过聚类中心的剔除,相同的模态可以由同一个聚类中心表示.但从图 13中可以看出,由于噪声的影响,即使初始聚类中心是从同一个稳定模态下抽取出来的,它们的距离有可能大于阈值 $\theta=2.5$ ,所以第一次迭代中同一模态下仍然保留了很多的聚类中心.另外,由于过渡模态变化比较大,而数据比较多,因此过渡模态会保留比较多的聚类中心.而在每次迭代过程中,计算聚类单元隶属度的时候,算法复杂度为 ${\rm O}(KN)$ ,在 $H=1$ 时,从图 6中可以看到,在第一次迭代中保留的聚类中心个数为 $N=18000$ ,迭代过程共有52次,且聚类中心个数均在20以上.这样,其复杂度远超于其他两种情况,这也是算法运行时间远远大于其他情况的原因.如果为了大大降低聚类中心的个数,将阈值 $\theta$ 增大,如图 13中所示,将 $\theta$ 取为3.8,这会导致稳定模态 $B$ 的某些数据被划分到代表稳定模态 $A$ 的聚类中.随着迭代过程的进行,聚类中心不断更新,可能导致将两种稳定模态划分到同一类的情况.通过分析,说明了如果对数据不加处理,直接进行聚类,由于奇异点等噪声的影响会导致算法复杂度高,参数选取困难的问题.
由图 14可知,当 $H=50$ 时,由于窗口的设定,窗口均值代替原数据作为聚类单元,这样会大大减少噪声的影响,因此属于同一个稳定模态的聚类中心大量被剔除.例如图 14中,从稳定模态 $A$ 抽取的聚类中心到第一个聚类中心的距离均在阈值 $\theta$ 之下,因此通过删减,稳定模态 $A$ 只用一个聚类中心表示.从图 9中可以看出,在第一次迭代时,聚类中心个数就骤减到18个,而聚类单元 $N$ 也只有360个,这是算法在4次迭代就可以达到最优的原因.
图 15表示了当 $H=500$ 时,初始聚类中心与第一个聚类中心的距离变化趋势.对比实验设定,发现其中第 $i$ 个聚类中心 $(i=1,2,9,11,12) $ 都是从稳定模态 $A$ 中抽取的原始聚类中心, 因此到第一个聚类中心的距离都在阈值之下,因此这些聚类中心在第一次迭代中就会被删减掉,只剩下第一个聚类中心.而第10个聚类中心应该表示扰动,但由于窗口设置过大,大大降低了扰动的影响,使其也归到了稳定模态 $A$ 中.这也再次证明了窗口设置过大,会使得聚类中心个数太少,且由于过渡时间远小于稳定模态时间,从而导致过渡模态的淹没.
另外,在模态细划分时,需要在相邻的两个窗口内寻找稳定模态与过渡模态的切割点,而大 $H$ 的设定会增大精确定位的范围.例如 $H=500$ 时,相较于 $H=50$ 的情况,算法运行时间提高的并不明显,但是精确定位的范围却扩大了10倍.
这里,我们选择 $H=50$ 作为初步识别的窗口设定.对图 7中稳定模态的淹没现象进行处理.根据第2节中步骤4的处理方法,确定新稳定模态,直到所有稳定模态之间的运行长度小于稳定模态最短运行时间 $S_{\min}$ .模态初步识别结果如图 16所示.
从第1 $\sim$ 60个窗口定义为稳定模态 $A$ ,第61 $\sim$ 70个窗口为稳定模态 $A$ 到模态 $B$ 的过渡模态,记为 $AB$ 过渡模态.第71 $\sim$ 140个窗口定义为稳定模态 $B$ ,从第141 $\sim$ 160个窗口为稳定模态 $B$ 到模态 $C$ 的过渡模态,记为 $BC$ 过渡模态.从第161 $\sim$ 226个窗口定义为稳定模态 $C$ ,从第227 $\sim$ 243个窗口定义为稳定模态 $C$ 到模态 $A$ 之间的过渡模态,记为 $CA$ 过渡模态.从第244 $\sim$ 280个窗口以及第296 $\sim$ 360 个窗口定义为稳定模态 $A$ .两个相邻稳定模态 $A$ 之间的所有时段,即从281 $\sim$ 295个窗口定义为稳定模态 $A$ 的扰动过程.与实验设定比较发现,稳定模态 $A$ , $B$ , $C$ 分别对应设定中的稳定模态 $A_0$ , $B_0$ , $C_0$ .
2.3 TE 过程准确划分
接下来需要准确定位各模态的起始终止位置.利用滑动窗口, 对稳定模态与过渡模态相邻的窗口进行分析.为了比较各种参数对识别结果的影响,小滑动窗口 $L$ 设置了两种情况: 1) $L=1$ ; 2) $L=5$ .滑动步长 $h$ 设置两种情况: 1) $h=1$ ; 2) $h=6$ .为了与根据数据特性计算的弹性阈值 $\alpha$ 进行比较,设置两个 $\alpha$ : 1) $\alpha=0.95$ ;2) $\alpha=0.85$ .
图 17和图 18表示在滑动步长 $h=1$ , $L$ 分别取不同值时,过渡模态 $AB$ 起始时刻的识别过程,分析了稳定模态 $A$ 与过渡模态 $AB$ 相连的5个切割窗口的数据与稳定模态 $A$ 的相似度变化趋势.横坐标表示滑动窗口初始位置所对应的采样点,纵坐标表示小滑动窗口均值与稳定模态的相似度.在模态初步识别时,过程从第61个切割窗口进入过渡过程.因此,过渡模态的精确起始时间应该位于第60和第61个切割窗口,对应图中第2951采样点到第3050个采样点.而第1个窗口到第29个切割窗口一定属于稳定模态 $A$ .为了得到相似度阈值,分析了第27、28和29共3个窗口数据与稳定模态 $A$ 的相似度.利用小滑动窗口 $L$ 对3个大窗口数据进行分析,对应图中2800采样点到2950采样点,根据核密度估计求取相似度阈值.
图 17 $AB$ 过渡过程起始位置识别 ( $L=1$ )Fig. 17 Identification at the beginning of $AB$ mode( $L=1$ )
图 18 $AB$ 过渡过程起始位置识别 ( $L=5$ )Fig. 18 Identification at the beginning of $AB$ mode( $L=5$ )图 17表示小窗口 $L=1$ 的情况,即在模态准确定位时,直接对数据进行分析,而不加以处理.过程实际是从第3001个采样点进入 $AB$ 过渡模态.从图中可以看出,当 $\alpha=0.95$ 时,识别起始位置为第2983个样本点,误差达到18个样本点.这是由于相似度阈值 $\alpha$ 设置过大,即使是稳定模态 $A$ 中的数据到稳定模态 $A$ 的相似度也有很多在阈值之下.而当 $\alpha=0.85$ 时,识别起始位置为第3007个样本点,误差为6个样本点.利用前3个大 $H$ 窗口计算的弹性阈值,识别起始位置为第3006个样本点.造成误差的原因是数据含有大量噪声,使得数据到稳定模态 $A$ 的相似度波动比较大.
图 18表示小窗口 $L=5$ 的情况,即在模态准确定位时,利用小滑动窗口对数据进行处理,求取窗口均值到前一个稳定模态 $A$ 的相似度.从图中可以看出,当 $\alpha=0.95$ 时,识别起始位置为第3001个样本点,误差为0.而当 $\alpha=0.85$ 时,识别起始位置为第3006个样本点,误差为5个样本点.而利用前3个大 $H$ 窗口计算的弹性阈值 $\alpha$ ,识别起始位置为第3001个样本点,误差为0.
图 19和图 20表示 $AB$ 过渡模态与下一个稳定模态 $B$ 相交接的5个大切割窗口数据与稳定模态 $B$ 的相似度变化趋势.根据初步识别, $AB$ 过渡过程终止时刻应该位于第69和第70个大切割窗口内.而71 $\sim$ 73这3个大窗口一定属于下一个稳定模态 $B$ ,因此分析这5个窗口,进行离线识别.
图 19 $AB$ 过渡过程终止位置识别 ( $L=1$ )Fig. 19 Identification at the end of $AB$ mode ( $L=1$ )
图 20 $AB$ 过渡过程终止位置识别 ( $L=5$ )Fig. 20 Identification at the end of $AB$ mode ( $L=5$ )从图 19中可以看出,由于没用小窗口对数据进行处理,噪声导致数据的相似度波动大,而且变化趋势非常不明显.当相似度阈值 $\alpha=0.95$ 时,前2个大 $H$ 窗口内数据到下一个稳定模态 $B$ 的相似度均在阈值之下.因此,将 $2H$ 窗口最后一个采样点作为其识别的过渡模态终止位置,即 $t=3550$ .而 $\alpha$ $=$ $0.85$ 时,识别终止时刻为 $t=3458$ .利用后3个大窗口数据计算的弹性相似度阈值 $\alpha$ ,识别终止时刻为第3482个采样点.图 20表示当 $L=5$ 时, $AB$ 过渡模态终止时刻的识别情况,分析过程类似.表 6表示了在其他参数相同,小滑动窗口 $L$ 取不同值时的识别结果的对比.通过比较发现,当 $L=5$ 时,识别效果远好于不设置小窗口( $L=1$ )的情况.这说明小窗口 $L$ 的设置,会大大减少噪声对数据的影响.另外,在其余参数一样的情况下,根据数据特性计算的弹性阈值 $\alpha$ , 效果明显好于另外两个硬性定义的相似度阈值参数.
表 6 识别结果比较 (不同 $L$ 值)Table 6 Comparison of clustering results (different $L$ )模态 小窗口 $L$ 相似度阈值 识别位置(采样点) 识别误差(采样点) $AB$ 起始时刻 $L=1$ $\alpha=0.95$ 2983 18 $\alpha=0.85$ 3007 6 弹性阈值 $\alpha$ 3006 5 $L=5$ $\alpha=0.95$ 3001 0 $\alpha=0.85$ 3006 5 弹性阈值 $\alpha$ 3001 0 $AB$ 终止时刻 $L=1$ $\alpha=0.95$ 3550 30 $\alpha=0.85$ 3458 62 弹性阈值 $\alpha$ 3482 38 $L=5$ $\alpha=0.95$ 3550 30 $\alpha=0.85$ 3451 31 弹性阈值 $\alpha$ 3518 2 图 21和图 22表示在小窗口长度 $L=5$ ,滑动步长 $h$ 分别取不同值时,过渡模态BC起始时刻的识别过程.图 23和图 24表示过渡模态BC终止时刻的识别过程.表 7表示了在其他参数相同,滑动步长 $h$ 取不同值时的识别结果的对比.通过比较可以发现,当 $h=6$ 时,需要分析的窗口个数明显小于 $h$ $=$ $1$ 时,提高了算法的运行速度.但由于滑动较大,可能真正的分割点被错过,导致误差增大.例如,在图 22中,实际起始点为7001,但由于步长较大,当第7003个采样点仍然在阈值之上时,下一个窗口对应的采样点直接滑动到第7009个采样点,从而增大了误差.实际上,由于大切割窗口 $H=50$ ,需要准确定位的范围并不大,因此滑动步长 $h$ 完全可以取1.
图 21 $BC$ 过渡过程起始位置识别 ( $L=1$ )Fig. 21 Identification at the beginning of $BC$ mode ( $L=1$ )
图 22 $BC$ 过渡过程起始位置识别 ( $L=5$ )Fig. 22 Identification at the beginning of $BC$ mode ( $L=5$ )
图 23 $BC$ 过渡过程终止位置识别 ( $h=1$ )Fig. 23 Identification at the end of $BC$ mode ( $h=1$ )
图 24 $BC$ 过渡过程终止位置识别 ( $h=6$ )Fig. 24 Identification at the end of $BC$ mode ( $h=6$ )表 7 识别结果比较(不同 $h$ 值)Table 7 Comparison of clustering results (different $h$ )模态 小窗口 $L$ 相似度阈值 识别位置(采样点) 识别误差(采样点) $BC$ 起始时刻 $h=1$ $\alpha=0.95$ 7001 0 $\alpha=0.85$ 7006 5 弹性阈值 $\alpha$ 7001 0 $h=5$ $\alpha=0.95$ 7009 8 $\alpha=0.85$ 7015 14 弹性阈值 $\alpha$ 7009 8 $BC$ 终止时刻 $h=1$ $\alpha=0.95$ 8050 22 $\alpha=0.85$ 7951 77 弹性阈值 $\alpha$ 8027 1 $h=5$ $\alpha=0.95$ 8045 17 $\alpha=0.85$ 7951 77 弹性阈值 $\alpha$ 8033 5 通过比较可知,在其余参数一样的情况下,根据数据特性计算的弹性阈值 $\alpha$ ,识别效果最好.通过图 17～24,可以发现由于各段数据的特性不尽相同,稳定模态上的数据到稳定模态中心的相似度也有较大的差别.因此,硬性的规定一个相似度阈值是不合理且难以确定的.而弹性相似度的定义很好地解决了该问题,且达到了非常好的效果.表 8给出了算法的最终识别结果,其中大切割窗口 $H=50$ ,小滑动窗口 $L=5$ ,滑动步长 $h=1$ ,相似度阈值 $\alpha$ 为计算的弹性阈值.
表 8 TE过程识别结果Table 8 Identification result of TE process模态 聚类单元 模态 聚类单元 模态 聚类单元 $A-1$ 1～3 000 $BC$ 7 001～8 027 $A-2$ 12 238 ～14 000 $AB$ 3 001～3 518 $C$ 8 028～11 250 扰动 14 001 ～14 725 $B$ 3 519～7 000 $CA$ 11 251～12 237 $A-3$ 14 726 ～18 000 比较表 3和表 8,可以发现算法识别误差基本都在10个采样点以下.而 $CA$ 过渡模态起始时刻与原设定相差较大.原实验设定在 $T=110$ h的时候,增大反应器压力,使稳定模态 $C$ 进入过渡模态,也就是说进入 $CA$ 模态时对应的采样点应为11001,但是算法得到的是11251,滞后一段较长的时间.图 25表示了 $CA$ 过渡模态的起始时刻的识别过程.分析这段滞后数据与稳定模态 $C$ 的相似度可以发现,虽然理论上从11001采样时刻就进入了过渡模态,但从11001 $\sim$ 11250的采样数据非常接近稳定模态 $C$ 的稳定工作状态.从图 25可以看出,此段数据与模态 $C$ 的相似度最高达到0.96,平均值在0.94左右,因此将这段数据归到模态 $C$ 中是非常合理的.
图 25 $CA$ 过渡过程起始位置识别( $L=5,h=1$ )Fig. 25 Identification at the beginning of $CA$ mode ( $L=5,h=1$ )通过以上仿真结果和过程分析不难发现,本文提出的算法可以有效的识别出各模态的起始终止时刻,且达到了比较理想的效果.
3. 结论
本文提出了一种多模态过程模态识别方法.对于离线模态识别的问题,先利用大分割窗口 $H$ 对过程数据进行粗分割,提取均值向量进行聚类,并针对聚类出现的稳定模态淹没现象作了相关处理,实现了离线模态的初步划分.然后用一个较短长度的滑动窗口对稳定模态与过渡模态相邻的几个窗口数据进行分析,通过分析相似度的变化趋势,定义了相似度阈值,从而精确定位各模态的起始终止时间,最终实现了多模态过程的自动模态识别."先粗后细"的模态识别方法提高了运算速度,对噪声、奇异点等随机扰动也有较强的冗余能力,同时提高了算法的精确性.该算法无需人为分析以及后续处理,能够真正实现多模态过程的全自动识别,并给出合理有效的识别结果.最后通过TE仿真,定量分析了不同参数取值对识别结果的影响,验证了该方法的有效性与准确性.
-
图 2 多模态准确识别过程(以 $AC$ 过渡模态起始时刻为例)
Fig. 2 Exact mode identification for multimode processes (i.e.beginning of the transitional mode $AC$ )
图 3 TE过程15个变量变化曲线 ( $X$ 轴表示采样点)
Fig. 3 15 variables of the TE process ( $X$ axis indicates sample points.)
图 5 前后两次迭代聚类中心的距离( $H=1$ )
Fig. 5 Distance between clustering centers in two iterations ( $H=1$ )
图 8 前后两次迭代聚类中心的距离( $H=50$ )
Fig. 8 Distance between clustering centers in two iterations ( $H=50$ )
图 10 K-means算法聚类结果 ( $H=500$ )
Fig. 10 Clustering result of K-means algorithm ( $H=500$ )
图 11 前后两次迭代聚类中心的距离 ( $H=500$ )
Fig. 11 Distance between clustering centers in two iterations( $H=500$ )
图 13 初始聚类中心到第一个聚类中心的距离( $H=1$ )
Fig. 13 Distance between original clustering center and the first clustering center ( $H=1$ )
图 14 初始聚类中心到第一个聚类中心的距离( $H=50$ )
Fig. 14 Distance between original clustering center and the first clustering center ( $H=50$ )
图 15 初始聚类中心到第一个聚类中心的距离( $H=500$ )
Fig. 15 Distance between original clustering center and the first clustering center ( $H=500$ )
图 17 $AB$ 过渡过程起始位置识别 ( $L=1$ )
Fig. 17 Identification at the beginning of $AB$ mode( $L=1$ )
图 18 $AB$ 过渡过程起始位置识别 ( $L=5$ )
Fig. 18 Identification at the beginning of $AB$ mode( $L=5$ )
图 19 $AB$ 过渡过程终止位置识别 ( $L=1$ )
Fig. 19 Identification at the end of $AB$ mode ( $L=1$ )
图 20 $AB$ 过渡过程终止位置识别 ( $L=5$ )
Fig. 20 Identification at the end of $AB$ mode ( $L=5$ )
图 21 $BC$ 过渡过程起始位置识别 ( $L=1$ )
Fig. 21 Identification at the beginning of $BC$ mode ( $L=1$ )
图 22 $BC$ 过渡过程起始位置识别 ( $L=5$ )
Fig. 22 Identification at the beginning of $BC$ mode ( $L=5$ )
图 23 $BC$ 过渡过程终止位置识别 ( $h=1$ )
Fig. 23 Identification at the end of $BC$ mode ( $h=1$ )
图 24 $BC$ 过渡过程终止位置识别 ( $h=6$ )
Fig. 24 Identification at the end of $BC$ mode ( $h=6$ )
图 25 $CA$ 过渡过程起始位置识别( $L=5,h=1$ )
Fig. 25 Identification at the beginning of $CA$ mode ( $L=5,h=1$ )
表 1 TE过程操作模态列表
Table 1 Operating modes of Tennessee Eastman process
模态类型 反应器压力设定值(kPa) 反应器液位(%) 稳定模态 $A_0$ 2800 65 稳定模态 $B_0$ 2600 65 稳定模态 $C_0$ 2400 65 扰动1 2800 75 
下载: 导出CSV
表 2 TE过程变量表
Table 2 Variables of Tennessee Eastman process
序号 变量名称 1 $A$ 进料量(流1) 2 $D$ 进料量(流2) 3 $E$ 进料量(流3) 4 $A,C$ 混合物料流量 5 再循环流量(流8) 6 反应器进料速度(流6) 7 反应器温度 8 排放速度(流9) 9 产品分离器温度 10 产品分离器压力 11 产品分离器塔底流量(流10) 12 汽提塔压力 13 汽提塔温度 14 反应器冷却水出口温度 15 分离器冷却水出口温度
下载: 导出CSV
表 3 TE过程实验设计
Table 3 Experimental design of TE process
模态 聚类单元 $A_0-1$ 1 $\sim$ 3000 $A_0B_0$ 3001 $\sim$ 3520 $B_0$ 3521 $\sim$ 7000 $B_0C_0$ 7001 $\sim$ 8028 $C_0$ 8029 $\sim$ 11000 $C_0A_0$ 11001 $\sim$ 12235 $A_0-2$ 12236 $\sim$ 14000 扰动 14001 $\sim$ 14720 $A_0-3$ 14721 $\sim$ 18000
下载: 导出CSV
表 4 各时段运行时间( $H=50$ )
Table 4 Runtime of each period ( $H=50$ )
时段 窗口个数 运行时间 (h) 所属类别 1 60 30 1 2 3 1.5 3 3 7 3.5 6 4 70 35 4 5 3 1.5 2 6 1 0.5 5 7 7 3.5 6 8 9 4.5 8 9 11 5.5 5 10 19 9.5 7 11 10 5 5 12 6 3 7 13 9 4.5 5 14 11 5.5 7 15 7 3.5 3 16 10 5 4 17 37 18.5 1 18 3 1.5 3 19 12 6 2 20 65 32.5 1
下载: 导出CSV
表 5 聚类结果比较(不同 $H$ 值)
Table 5 Comparison of clustering results (different $H$ )
$H$ 值 聚类单元 初始类别数量 最终类别数量 迭代次数 运行时间(s) 1 18000 6000 20 52 519.2 50 360 120 8 4 0.906 500 36 12 3 3 0.745
下载: 导出CSV
表 6 识别结果比较 (不同 $L$ 值)
Table 6 Comparison of clustering results (different $L$ )
模态 小窗口 $L$ 相似度阈值 识别位置(采样点) 识别误差(采样点) $AB$ 起始时刻 $L=1$ $\alpha=0.95$ 2983 18 $\alpha=0.85$ 3007 6 弹性阈值 $\alpha$ 3006 5 $L=5$ $\alpha=0.95$ 3001 0 $\alpha=0.85$ 3006 5 弹性阈值 $\alpha$ 3001 0 $AB$ 终止时刻 $L=1$ $\alpha=0.95$ 3550 30 $\alpha=0.85$ 3458 62 弹性阈值 $\alpha$ 3482 38 $L=5$ $\alpha=0.95$ 3550 30 $\alpha=0.85$ 3451 31 弹性阈值 $\alpha$ 3518 2
下载: 导出CSV
表 7 识别结果比较(不同 $h$ 值)
Table 7 Comparison of clustering results (different $h$ )
模态 小窗口 $L$ 相似度阈值 识别位置(采样点) 识别误差(采样点) $BC$ 起始时刻 $h=1$ $\alpha=0.95$ 7001 0 $\alpha=0.85$ 7006 5 弹性阈值 $\alpha$ 7001 0 $h=5$ $\alpha=0.95$ 7009 8 $\alpha=0.85$ 7015 14 弹性阈值 $\alpha$ 7009 8 $BC$ 终止时刻 $h=1$ $\alpha=0.95$ 8050 22 $\alpha=0.85$ 7951 77 弹性阈值 $\alpha$ 8027 1 $h=5$ $\alpha=0.95$ 8045 17 $\alpha=0.85$ 7951 77 弹性阈值 $\alpha$ 8033 5
下载: 导出CSV
表 8 TE过程识别结果
Table 8 Identification result of TE process
模态 聚类单元 模态 聚类单元 模态 聚类单元 $A-1$ 1～3 000 $BC$ 7 001～8 027 $A-2$ 12 238 ～14 000 $AB$ 3 001～3 518 $C$ 8 028～11 250 扰动 14 001 ～14 725 $B$ 3 519～7 000 $CA$ 11 251～12 237 $A-3$ 14 726 ～18 000
下载: 导出CSV
-
[1] 谭帅. 多模态过程统计建模及在线监测方法研究[博士学位论文], 东北大学, 中国, 2012.Tan Shuai. Statistical Modeling and Online Monitoring for Multiple Mode Processes[Ph.D. dissertation], Northeastern University, China, 2012. [2] Yu J, Qin S J. Multimode process monitoring with Bayesian inference-based finite Gaussian mixture models. AIChE Journal, 2008, 54(7):1811-1829 doi: 10.1002/(ISSN)1547-5905 [3] 王静, 胡益, 侍洪波. 基于GMM的间歇过程故障检测. 自动化学报, 2015, 41(5):899-905 http://www.aas.net.cn/CN/abstract/abstract18664.shtmlWang Jing, Hu Yi, Shi Hong-Bo. Fault detection for batch processes based on Gaussian mixture model. Acta Automatica Sinica, 2015, 41(5):899-905 http://www.aas.net.cn/CN/abstract/abstract18664.shtml [4] Xie X, Shi H B. Dynamic multimode process modeling and monitoring using adaptive Gaussian mixture models. Industrial and Engineering Chemistry Research, 2012, 51(15):5497-5505 doi: 10.1021/ie202720y [5] Ge Z Q, Gao F R, Song Z H. Mixture probabilistic PCR model for soft sensing of multimode processes. Chemometrics and Intelligent Laboratory Systems, 2011, 105(1):91-105 doi: 10.1016/j.chemolab.2010.11.004 [6] Ge Z Q. Mixture Bayesian regularization of PCR model and soft sensing application. IEEE Transactions on Industrial Electronics, 2015, 62(7):4336-4343 doi: 10.1109/TIE.2014.2385042 [7] Ng Y S, Srinivasan R. An adjoined multi-model approach for monitoring batch and transient operations. Computers and Chemical Engineering, 2009, 33(4):887-902 doi: 10.1016/j.compchemeng.2008.11.014 [8] Lu N Y, Gao F R, Wang F L. Sub-PCA modeling and on-line monitoring strategy for batch processes. AIChE Journal, 2004, 50(1):255-259 doi: 10.1002/(ISSN)1547-5905 [9] Zhao C H, Wang F L, Lu N Y, Jia M X. Stage-based soft-transition multiple PCA modeling and on-line monitoring strategy for batch processes. Journal of Process Control, 2007, 17(9):728-741 doi: 10.1016/j.jprocont.2007.02.005 [10] Zhao C H, Zhang W D. Reconstruction based fault diagnosis using concurrent phase partition and analysis of relative changes for multiphase batch processes with limited fault batches. Chemometrics and Intelligent Laboratory Systems, 2014, 130:135-150 doi: 10.1016/j.chemolab.2013.10.014 [11] Tang X C, Li Y, Xie Z. Phase division and process monitoring for multiphase batch processes with transitions. Chemometrics and Intelligent Laboratory Systems, 2015, 145:72-83 doi: 10.1016/j.chemolab.2015.04.007 [12] Wang F L, Tan S, Peng J, Chang Y Q. Process monitoring based on mode identification for multi-mode process with transitions. Chemometrics and Intelligent Laboratory Systems, 2012, 110(1):144-155 doi: 10.1016/j.chemolab.2011.10.013 [13] Tan S, Wang F L, Peng J, Chang Y Q, Wang S. Multimode process monitoring based on mode identification. Industrial and Engineering Chemistry Research, 2012, 51(1):374-388 doi: 10.1021/ie102048f [14] Zhang Y W, Zhang H L. Fault detection for time-varying processes. IEEE Transactions on Control Systems Technology, 2014, 22(4):1527-1535 doi: 10.1109/TCST.2013.2273498 [15] Zhang Y W, Li S. Modeling and monitoring between-mode transition of multimodes processes. IEEE Transactions on Industrial Informatics, 2013, 9(4):2248-2255 doi: 10.1109/TII.2012.2220977 [16] Alguwaizani A. Degeneracy on K-means clustering. Electronic Notes in Discrete Mathematics, 2012, 39:13-20 doi: 10.1016/j.endm.2012.10.003 [17] 潘天红, 薛振框, 李少远. 基于减法聚类的多模型在线辨识算法. 自动化学报, 2009, 35(2):220-224 doi: 10.3724/SP.J.1004.2009.00220Pan Tian-Hong, Xue Zhen-Kuang, Li Shao-Yuan. An online multi-model identification algorithm based on subtractive clustering. Acta Automatica Sinica, 2009, 35(2):220-224 doi: 10.3724/SP.J.1004.2009.00220 [18] Chiu S L. Fuzzy model identification based on cluster estimation. Journal of Intelligent and Fuzzy Systems:Applications in Engineering and Technology, 1994, 2(3):267-278 http://www.oalib.com/references/16296749 [19] Jain A K, Murty M N, Flynn P J. Data clustering:a review. ACM Computing Surveys, 1999, 31(3):264-323 doi: 10.1145/331499.331504 [20] 王惠文. 偏最小二乘回归方法及其应用. 北京:国防工业出版社, 1999.Wang Hui-Wen. Partial Least-Squares Regression-Method and Applications. Beijing, China:National Defence Industry Press, 1999. [21] 钱鹏江, 王士同, 邓赵红. 快速核密度估计定理和大规模图论松弛聚类方法. 自动化学报, 2011, 37(12):1422-1434Qian Peng-Jiang, Wang Shi-Tong, Deng Zhao-Hong. Fast kernel density estimate theorem and scaling up graph-based relaxed clustering method. Acta Automatica Sinica, 2011, 37(12):1422-1434 [22] 蒋浩天, 拉塞尔 E L, 布拉茨 R D. 工业系统的故障检测与诊断. 北京:机械工业出版社, 2003.Chiang L H, Russell E L, Braatz R D. Fault Detection and Diagnosis in Industrial Systems. Beijing, China:China Machine Press, 2003. [23] Larsson T, Hestetun K, Hovland E, Skogestad S. Self-optimizing control of a large-scale plant:the Tennessee Eastman process. Industrial and Engineering Chemistry Research, 2001, 40(22):4889-4901 doi: 10.1021/ie000586y 期刊类型引用(18)
1. 季新芳，张勇，巩敦卫，郭一楠，孙晓燕. 异构集成代理辅助的区间多模态粒子群优化算法. 自动化学报. 2024(09): 1831-1853 . 
本站查看2. 任超，叶泽甫，程兰，乔铁柱，阎高伟. 基于测地线流式核的隐空间多工况软测量建模. 传感器与微系统. 2023(10): 28-31 . 百度学术3. 韩鹏东，阎高伟，任密蜂，程兰，叶泽甫. 基于迁移子空间学习的偏最小二乘回归软测量方法. 控制与决策. 2023(11): 3147-3155 . 百度学术4. 曹跃，余冲，李智，杨明磊. 工业数据驱动的加氢裂化装置多工况切换过渡状态检测. 化工学报. 2023(09): 3841-3854 . 百度学术5. 宋冰，郭涛，侍洪波，谭帅，陶阳，马浴阳. 基于双子空间并行回归的化工过程质量相关故障检测方法. 化工学报. 2023(11): 4600-4610 . 百度学术6. 董峰，李昭，李凌涵，张淑美. 多模态动态核主成分分析的气液两相流状态监测. 自动化学报. 2022(03): 762-773 . 本站查看7. 徐志强，王芳，程兰，陈泽华，阎高伟. 基于工况差异度量的多源即时学习软测量. 控制工程. 2022(10): 1768-1773 . 百度学术8. 朱明睿，纪杨建，甘红宇，张念. 基于近邻转移约束规则的非确定工业过渡过程的模态识别方法. 计算机集成制造系统. 2022(11): 3576-3587 . 百度学术9. 任超，阎高伟，程兰，王芳. 基于最大均值差异的多模态过程过渡模态识别方法. 浙江大学学报(工学版). 2021(03): 563-570 . 百度学术10. 来颜博，阎高伟，程兰，陈泽华. 基于动态独立成分分析和动态主成分分析的测地线流式核无监督回归模型. 上海交通大学学报. 2020(12): 1269-1277 . 百度学术11. 王通，段泽文，张文喜. 基于窗口样本相似因子分析的油井工况识别方法. 沈阳工业大学学报. 2019(06): 681-686 . 百度学术12. 孔祥玉，曹泽豪，杜柏阳，罗家宇. 基于偏最小二乘的质量相关多模态故障检测技术. 控制与决策. 2019(12): 2547-2557 . 百度学术13. 王通，段泽文. 基于模糊评估自适应更新的油井动液面软测量建模. 化工学报. 2019(12): 4760-4769 . 百度学术14. 熊伟丽，郭校根. 一种基于多工况识别的过程在线监测方法. 控制与决策. 2018(03): 403-412 . 百度学术15. 郑皓，熊伟丽. 基于多阶段ICA-SVDD的间歇过程故障监测. 南京理工大学学报. 2018(02): 195-203 . 百度学术16. 陈晓方，钱荧灿，王雅琳，阳春华. 基于主元导数特征聚类的加氢裂化动态调整区间识别. 清华大学学报(自然科学版). 2018(01): 81-86 . 百度学术17. 王署东，尹柏强，何怡刚，李兵，佐磊，侯金波. 高压局部放电定位模型转换求解与改进K-means聚类优化方法. 电子测量与仪器学报. 2018(11): 178-186 . 百度学术18. 南男，杨健，赵晶晶，侍洪波. 基于谱聚类特征向量分析的模态划分方法. 华东理工大学学报(自然科学版). 2017(05): 669-676 . 百度学术其他类型引用(18)
-
下载:
计量
- 文章访问数: 2370
- HTML全文浏览量: 352
- PDF下载量: 1184
- 被引次数: 36
