国家标准 GB/T40859一2021 流式数据监测控制图 Controlchartforstreaming datamonmitoring 2021-10-11发布 2022-05-01实施国家市场监督管理总局发布国家标涯花管理委员会国家标准
GB/T40859一2021 目次前言引言范围 2 规范性引用文件术语和定义符号应用条件基于离散小波变换的数据预处理单变量流式数据控制图的构建 7.1概述 7.2单变量流式数据EWMA控制图的构建 7.3单变量流式数据EWMS控制图的构建 7.!单变量流式数据控制图应用步骤多变量流式数据控制图的构建 8.1概述 8.2小波系数重组 8.3多变量流式数据MEwMA控制图的构建 8.4多变量流式数据MEwMC控制图的构建 8.5多变量流式数据控制图应用步骤附录A资料性)离散小波变换及多分辨率分析原理 10 附录B(资料性)单变量流式数据应用实例 12 附录c(资料性多变量流式数据应用实例 15 附录D(资料性)代码实现 18 参考文献 29
GB/T40859一2021 前言本文件按照GB/T1.1一2020<标准化工作导则第1部分;标准化文件的结构和起草规则》的规定起草请注意本文件的某些内容可能涉及专利本文件的发布机构不承担识别专利的责任本文件由全国统计方法应用标准化技术委员会(SAC/TC21)提出并归口本文件起草单位;清华大学、标准化研究院本文件主要起草人:王凯波、杨洛、丁文兴、张帆、赵静
GB/T40859一2021 引言在现代复杂生产系统和装备运行过程中,由于自动化传感器的广泛应用,所获取的数据往往是流式数据的形式所谓流式数据,是一组数量庞大、顺次有序、快速高频到达的数据序列现有的控制图往往基于较大间隔进行离散采样获得的数据进行过程控制由于流式数据具有快速高频等特点,现有控制图通常无法适用于控制产生流式数据的过程传统的数据采集与分析往往是低频的,如每小时对一条生产水杯的生产线进行产品采样和关键质量指标测量,以判断生产过程的状态是否正常而在复杂系统的生产制造和运行过程中,通常会布局很多设备状态传感器和自动化检测仪器,以秒甚至更小的时间单位为间隔,对设备运行状态及产品质量指标进行实时、自动和持续的数据收集,由此产生流式数据这些流式数据包含比低频数据更加复杂,全面的信息,但其信息密度也相对较低,因此,需要合适的特征提取手段,挖掘流式数据中包含的信息 20世纪,有学者提出了多分辨率分析方法这是一种基于离散小波变换的数据分析方法,可以将数据分解到不同的频率尺度下进行分析此后小波分析的理论方法不断完酱多分辨率小波分析十分适用于处理高频的流式数据将高颊的流式数据分解到不同的频率尺度后,可在不同的尺度下分别分析,提取流式数据中复杂的信息本文件可对高速铁路、装备制造,复杂产品制造等系统在运行过程中所产生的实时状态数据和产品质量数据进行监测,对关键信号的异常变化进行预警,有可能在早期发现危险故障,提升系统安全可靠性,提高产品质量 IN
GB/T40859一2021 流式数据监测控制图范围本文件给出了对单变量及多变量流式计量型数据进行监测的控制图方法本文件适用于对产生流式数据的一个或多个有一定关联关系的过程或质量变量,在均值和基本波动特征已知的条件下,对偏移和波动进行控制规范性引用文件下列文件中的内容通过文中的规范性引用而构成本文件必不可少的条款其中,注日期的引用文件,仅该日期对应的版本适用于本文件;不注日期的引用文件,其最新版本(包括所有的修改单)适用于本文件 GB/T3358.1统计学词汇及符号第1部分;一般统计术语与用于概率的术语 GB/T3358.2统计学词汇及符号第2部分:应用统计术语和定义 r358.1和GB/T358.2界定的以及下列术语和定义适用于本文件 GB/T 3.1 流式数据streamingdata -组数量庞大、顺次有序、快速高频到达的数据序列注1传感器监测数据,交通监测数据等往往都是以流式数据的形式被收集注2:流式数据通常是由一个或多个数据源持续产生,随时间延续而增长的动态数据集合,具有采样间隔小、,采样频率高,到达速度快,数据量大等特点 3.2 EwMMA控制图EwMIAchart averagecontrolchart 指数加权移动平均控制图exponentiallyweightedmoving 用指数平滑加权平均评估和监测过程水平的计量控制图 [来源:GB/T3358.2一2009,2.3.16,有修改] 3.3 EWMS控制图EWMSchart 指数加权均方误差控制图exponentiallyweightedmeansquarederrorcontrolchart 用指数平滑加权平均评估和监测过程偏移与变异的计量控制图 3.4 MEWMA控制图MEWMAchart 多变量指数加权移动平均控制图multivariateeponentialyweightedmovingaverageeontrolehart 将两个或两个以上相关的变量合成一个样本统计量,并用指数平滑加权平均评估和监测过程水平的计量控制图
GB/T40859一2021 3.5 MEWC控制图MIEwMcchart 多变量指数加权移动协方差矩阵控制图 multivariateexponentiallyweightedmm0vingcovariance matrixcontrolchart 将两个或两个以上相关的变量合成一个样本统计量,并用指数平滑加权平均评估和监测过程变异的计量控制图 3.6 小波变换waveletransform 用样本数据与小波函数的内积计算得到小波系数的过程 3.7 离散小波变换disecretewavelettransform" 用样本数据与正交小波基的内积计算得到小波系数的过程 3.8 小波基函数waveletbasisfunetion 母小波motherwavelet 均方可积且均值为0的函数注1:常用小波基函数有Haar小波,Meyer小波等注2:小波基函数通常表示为() 3.9 小波函数wavelets 用小波基函数伸缩与平移得到的一组丽数注:小被丽数可表示为的a.(o)=团-((a-1一b).其中“为伸缩参数.南为平移参数，为尺度参数.为位置参数 3.10 orthonoralwaveletbases 正交小波基用某些小波基函数进行特定的离散化伸缩与平移,得到的可作为正交基的小波函数注1，Haar小波是可用于产生正交小波基的一种小波基丽数注2:通常选取伸缩参数a=2,平移参数6=1 3.11 小波系数waveleteoefficient 用样本数据与小波函数内积的结果注内积过程可表示为d,=()9,.()e"d,其中)为样本数据 3.12 父小波fatherwavelet 尺度函数scalingfunction 可生成正交基的范数为1的规范化函数,与母小波对应注1:并不是所有的母小波都有对应的父小波注2:可用于离散小波变换的母小波都有对应的父小波注3，母小波用于表现样本波动信息,父小波用于表现样本近似信息 3.13 多分辨率分析mutiresolutionanalysis 用离散小波变换得到的小波系数进行数据分析的方法注离散小波变换得到的小波系数包含不同频率尺度下的信息
GB/T40859一2021 3.14 近似层小波系数approximatewaveletcefieient 用样本数据与父小波内积的结果注1:近似层小波系数包含样本趋势信息注2:样本经离散小波变换通常仅保留一层近似层小波系数 3.15 细节层小波系数detailwavelecoeffieient 用样本数据与母小波内积的结果注1:细节层小波系数包含样本波动信息注2样本经离散小波变换会产生若干层近似层小波系数符号下列符号适用于本文件 ARL 平均链长 ARL 过程受控状态下的平均链长第」层第k个近似层小波系数第」层第人个近似层小波系数向量 c1.A d 第层第个细节层小波系数 d 第层第人个细节层小波系数向量 MEwMA控制图控制限 he MEwMC控制图控制限力维的单位矩阵离散小波变换分解层数 EwMA控制图控制限系数多变量数据的维度第个EwMS控制图统计量 La 下控制限上控制限多变量流式数据单变量流式数据第i维变量 .r MEwMA控制图统计量第j个MEwMC控制图统计量 y EwMA控制图统计量 EwMA控制图平滑系数 EwMs控制图平滑系数 MEwMA控制图平滑系数 MEwMC控制图平滑系数近似层小波系数均值的目标值近似层小波系数向量均值的目标值展第层细节层小波系数均值的目标值 Al 第层细节层小被系数向量均值的目标值从a
GB/T40859一202 近似层小波系数标准差的目标值第层细节层小波系数标准差的目标值近似层小波系数向量均值协方差矩阵的目标值第层细节层小波系数向量均值协方差矩阵的目标值母小波应用条件 S 现实的生产过程中,由于各种传感器的大量使用,过程样本数据往往以流式数据的形式被收集流式数据具有采样频率高、采样间隔短、数据量大等特点由于数据获取频率高,短时间内便能获得大量数据,提高了数据计算和存储的难度,故流式数据需采用新的监测方法此外,生产过程中的过程信号往往十分复杂,可能同时存在着多种不同的特征信号这些复杂的信号组合在一起会导致难以对故障的原因进行诊断,更严重的是有时多种特征信号的叠加可能会导致一些特征信号难以被监测多种特征信号叠加的数据如图1所示 wwM W人 wM V \ 说明: 白噪声; .ZN -包含均值突变的数据工B 包含波动突变的数据" .zC -rA,r、,r叠加后的数据 rD 复杂的特征信号从图1中可看出.发生了均值的偏移,r发生了波动幅度的变化虽然这二者变化特征明显但是对于由rA、rg、re叠加得到的数据rp,其均值偏移与波动变化两种特征信号互相影响若使用单 -控制图监测rp,则对均值偏移与波动变化都难以进行有效识别若可以将图1中叠加在一起的复杂信号分开进行监测,那么便会大大提高监测的效果离散小波变换是一种对数据进行分解的方法,原始数据经过离散小波变换后,可在不同的频率尺度上监测不同的特征信号此外,即使原始数据存在一定自相关性,经离散小波变换得到的小波系数的自相关性会显著
GB/T40859一2021 降低或消失,亦即离散小波变换所得到的系数对原始数据中存在的自相关性不敏感,因此,离散小波变换的方法可适用于流式数据的监测中在使用本控制图对过程进行监测时,应有一定量能够代表过程正常运行状态的历史数据,用于计算控制图参数的目标值基于离散小波变换的数据预处理为了实现对流式数据样本的监测,需要首先对所获取数据应用基于Haar小波的离散小波变换算法进行预处理,生成小波系数,然后将小波系数用于控制图监测 Haar小波是一种常用的正交小波,其定义为: (01) 结合附录A中介绍的小波函数的公式(A.1),可知其在第层的小波函数的定义域长度仅为2',因此仅使用有限的数据点即可计算一个小波系数假定分解层数为了,流式数据序列为x={r,r , 近似层小波系数和细节层小波系数的计算公式分别为 =1x当 cJ. .r(-1)2十 =2+xC d. r(-1)2/,E1,, -1)2+ 在监测过程中,需对新获取的数据点进行计数假定计数为i,当i可以整除2时,计算前层小波系数表1给出了应用小波变换得到细节层小波系数的示例表中每一列表示一个新获取的数据每一行表示小波分解的层,表内的d;表示第层的第人个小波系数表中共展示了8个数据点,进行 3层小波分解的系数表中可看出.当计数为2时,可通过1、2两个数据点计算获得一个第一层的小波系数d.;计数为8时,可通过18八个数据点计算第三层小波系数d.,同时通过5,6、7,8四个数据点计算第二层小波系数da ,通过7,8两个数据点计算第一层小波系数d.a 当获取新的数据时,可重复该过程,不断获取新的小波系数表1小波变换样例计数层数 d. dl1，？ d. d1. dl2， d2， da. 通过该算法,可将实时获得的数据及时分解为各层小波系数,共可得到1层近似层系数{c.1cya 以及」层绸节层小波系数(d,.,da)jEl),用于对过程的监测在构建控制图的过程中,分解层数的选择取决于流式数据特征信号的频率范围,若要监测较低频的波动信号则选取较大的分解层数J 但分解层数」过大会导致计算复杂度增加,影响监测效率通常分解层数建议选为 2或3即可在应用控制图时,往往假设被监测特征服从正态分布在流式数据服从正态分布的假设下,基于 Haar小波的数据预处理得到的各层小波系数仍然满足正态分布的假设,可使用控制图进行监测在过
GB/T40859一202 程受控状态下获取的数据经预处理得到的各层小波系数,可用于计算控制图参数的目标值单变量流式数据控制图的构建 7.1概述流式数据经小波变换分解为小波系数后,使用EwMA控制图监测近似层小波系数,使用EwMs 控制图监测各层细节层小波系数附录B中给出了构建控制图的案例 7.2单变量流式数据EwMIA控制图的构建 EwMA控制图适用于监测数据的均值变化应使用经离散小波变换得到的近似层小波系数 c.,ce.3,}构建该控制图 EwMA统计量可表示为义，=Ac,+(1一A))- 其中,入E(0,l)为平滑系数,文0=A 将公式(4)展开,可以得到 (5 =;,=A,c,十a,QI一A c,-I十A，(1-Ai)'c-;十可以看出,EwMA控制图将所有以往的信息都包括进来计算统计量距离当前时间越远,该数据点的信息在此时刻的统计量中的占比越小 EwMA控制图的控制限为十lo 6 Ua=从一La Lc=A c 其中,为控制限系数应用ARL分析的方法可以确定控制图的两个参数;控制限系数!和平滑系数入I 通常,A的取值取决于偏移量的大小,一般推荐范围为0.05<入1<0.5 ARL=370时的参数选择如表2所示表2ARL =370时EwMIA控制图参数选择标准 0.05 0.1 0.2 0.3 0.4 0.5 2.488 2.715 2.864 2.928 2.961 2.979 7.3单变量流式数据EwMs控制图的构建 EwMS控制图适用于监测数据的方差变化应使用经离散小波变换得到的各层细节层小波系数 ).jE(1),共构建个EwMs控制图 d.1,d,2, 其中第个控制图的EwMS统计量可表示为 si，=入d,一M小)'(1一入2)s- 8 其中AE(0.1)为平滑系数二i 通常,绘制关于s,的控制图,而非关于s.的控制图,即将迭代序列的平方根绘制人控制图 EwMS控制图的控制限为: Ue,= =dd
GB/T40859一2021 Lc!=o 10 山其中,表示发生第一类错误的概率,X-=()和X()分别表示自由度为的卡方分布的(I a/2)×100%和(a/2)×100%分位数,=(2-入,/入入，的取值取决于偏移量的大小,一般推荐范围为0.01<入,<0.4 ARL=370时(a=0.0027),平滑系数入，与控制限的对应关系如表3所示表3ARL=370时EwMs控制限选择标准 0.01 0.02 0.05 0.10 0.20 0.33 Ua 1.15a 1.22o 1.35d小 1.50a小 1.73aa 1.99a小 Lct 0.85d 0.79a 0,67o 0.,54a 0.37da 0,22a 7.4单变量流式数据控制图应用步骤单变量流式数据控制图构建的具体步骤如下步骤1确定控制图参数;选取离散小波变换分解层数1、,EwMA控制图的控制限L和平滑系数入1,EwMS控制图的第一类错误的概率a和平滑系数入g,并计算各控制图的控制限步骤2计算小波系数;对新获取的数据r，进行计数,代人公式(2)、公式(3)获取新的近似层或细节层小波系数步骤3计算统计量;使用新获得的小波系数,近似层小波系数代人公式(4)计算EwMA统计量，各层细节层小波系数分别代人公式(8)计算各层的EwMS统计量步骤4绘制控制图:将各统计量及相应的控制限分别绘制在控制图上,共得到一张EwMA控制图及J张EwMS控制图根据控制图判断受控状态当任何一张控制图出现控制点超过控制限时，表明过程可能“失控”,需要对过程采取必要的措施,予以纠正多变量流式数据控制图的构建 8.1概述假设某一时刻开始获得多变量流式数据序列x=(x1x),变量的维度为声，即x1- cP,,c"')' 为达到监测目的,首先使用离散小波变换将每一个变量的数据分解为近似层小波系数和细节层小波系数,再将不同变量同一层的小波系数重组，构成多维的小波系数向量使用 MEwMA控制图监测近似层小波系数向量,使用MEwMC控制图监测各细节层小波系数向量附录 C中给出了构建控制图的案例 8.2小波系数重组使用离散小波变换将每一维数据都分解至层以2维数据为例;第1维数据{rP,rg,}可分解得到1层近似层系数c',c,}以及层细节层小波系数d,,d,),je{1,,); 第2维数据{.zP,r'P,)可分解得到1层近似层系数(c',c'3,}以及」层细节层小波系数 d,d3,,jE1,, 将2维数据的近似层系数组成向量,可以得到1层近似层系数向量{ce.1,c2,},其中c.1= e,c);将2维数据对应的细节层系数组成向量可得到1层细节层小波系数d1,dy2,)， E(1,,J},其中d.1=(d,d. 图2展示了变量维度为力时小波系数的重组过程
GB/T40859一202 r0 小波变换小波变换小波变换近似层小波细节层小波细节层小波系数向量e 系数向量D 系数向量D MEwMA控制图 MEwNMC控制图 MEWN ne教制间图2多变量流式数据控制图小波系数重组示意图 8.3多变量流式数据MEwMA控制图的构建 MEwMA控制图适用于监测多变量流式数据的均值变化应使用8.2中经离散小波变换并重组得到的近似层小波系数向量eea)构建控制图1 定义向量序列Z,Z=u. 每获取一个新的小波系数向量e.,便可通过迭代计算得到2 (11 Z=入sc.+(1一入)Z- 其中,常数入E(0,1)是平滑系数判断数据失控的条件为 Y=(z一M);'(z一M.)>h 12 其中,Y为该控制图的统计量,h为MEwMA控制图的控制限;习，为Z，的协方差矩阵,其可表示为: 13 》,=一工应用ARL分析的方法可以确定控制图的两个参数;控制限h，和平滑系数入入，的取值取决于偏移量的大小,一般推荐范围为0.05GB/T40859一2021 8.4多变量流式数据MIEwC控制图的构建 MEwMC控制图适用于监测多变量流式数据的方差及协方差变化应使用8.2中经离散小波变换并重组得到的各层细节层小波系数向量'd.1,d2,},jE1,,J}构建控制图构建第个MEwMC控制图时,定义向量序列s4,s.0=1, 每获取一个新的小波系数向量 d4,先进行一步转换: 14 U)1=A(d.i一M山其中,矩阵A使得A2A'=1, 使用U.通过迭代计算得到sa 15 ;=入,U.,U'.4+(1一入,)s sj i-1 其中,常数入E(0,l)是平滑系数判断数据失控的条件为: 16 =tr(s))一lnls一p>he yj 其中,y,,为该控制图的统计量，力为变量的维度,he为MEwWMC控制图的控制限，/r()为求矩阵迹的诵数构建该控制图同样需要确定两个参数;控制限A 和权重系数A. A，的取值取决于偏移量一般推荐范围为0.01<入的大小, <0.4 应选择合适的h 使得控制图可达到一定的受控ARL,ARL =370时hc选择的原则如表5所示表5ARL=370时MIEwMC控制图控制限选择标准 0.05 0,15 0,30 0.10 0.,20 0,.25 0,.27 0.63 1.01 1.43 1.85 2.32 0.59 1.29 2.04 2.85 3.71 4.61 0.99 2.16 3,41 4.74 6,14 7.63 l.49 3.25 5.ll 7.12 9,.23 ll1.51 2.12 4.56 7.17 l0 l0,01 13.04 16.34 8.5多变量流式数据控制图应用步骤多变量流式数据控制图构建的具体步骤如下步骤1确定控制图参数;选取离散小波变换分解层数J、MEwMA控制图的平滑系数入，、 MEwMC控制图的平滑系数入,并选择各控制图的控制限步骤2计算小波系数向量;对新获取的数据x，进行计数,并分别代人公式(2)、公式(3)获取新的近似层或细节层小波系数,再将对应层小波系数重组为小波系数向量步骤3计算统计量;使用新获得的小波系数向量,近似层小波系数向量代人公式(l1)、公式(12 和公式(l3)计算MEWMA统计量,各层细节层小波系数向量分别代人公式(14)、公式(15)和公式(16) 计算各层的MEwWMC统计量步骤4绘制控制图;将各统计量及相应的控制限分别绘制在控制图上,共得到一张MEwMA控制图及张MEWMC控制图根据控制图判断受控状态当任何一张控制图出现控制点超过控制限时,表明过程可能“失控”,需要对过程采取必要的措施,予以纠正以上步骤往往需要使用计算机辅助实现,实现控制图的Python代码在附录D中列出
GB/T40859一2021 附录 A 资料性离散小波变换及多分辨率分析原理多分辨率数据分析的原理是利用离散小波少换将数据分解到不同的尺度上进行分析,这对于流式数据的监测十分有效离散小波变换是一种时频分析算法,与传统的傅里叶变换不同它能保留原始数据在时域上的信息,即可将原始数据分解到不同的频率尺度上用于过程监测的流式数据通常可通过小波函数分解成近似信号和细节信号的形式 (A.1 f(t= 之'f(+之.妈,.G) 其中,和大分别为尺度参数和位置参数,中是低频的小波基父函数,必是高频的母小波;c.，和d,分别为近似小波系数和细节小波系数其表达式为: 中.()=2(2'一k),e，=(,内.,》 r).()ed (A.2 ,.()=2t(e一),d小=(,的.)-[)9,.(o)e"山d A.3) 当k一m2',mZ时,这一小波变换的方法被称作离散小波变换在分解时,首先确定分解层数 J,之后依次计算各层的小波系数,在每层中以k=i2',iE乙逐步平移小波基函数共可得到1层的近似小波系数和」层细节小波系数,如图A.1所示第一层细节小波系数第二友绷节小波系数原始数豁第三层第一层细节小波系数近似小波系数第二层近似小谈系数第三层近似小波系数图A.1离散小波变换小波系数每一个小波系数都是由一段原始数据与该层小波基函数的内积得到,即该小波系数可以描述原始数据在一段时间内与该层小波基函数的相似程度因此不同层的小波系数描述了原始信号在不同频域上的信息,如图A.2所示 10
GB/T40859一2021 rM儿wwMm wv从小wm w 八M" 说明: 原始数据; 第一细节层小波系数第二细节层小波系数 d 近似层小波系数图A.2小波变换效果从图A.2中可看出,原始信号包含有均值偏移和方差变化两种不同的特征信号利用离散小波变换将其分解到2层后,可看出均值偏移的信息表现在近似小波系数上,而波动变化的信息表现在细节系数上因此,可对分解后的各层小波系数分别进行监测,这便是多分辨率分析的原理 11
GB/T40859一2021 录附 B 资料性) 单变量流式数据应用实例在某直拉单晶硅生产过程中,需要对单晶的关键质量参数进行监测记该质量参数为变量x 单晶拉制过程中,数据的获取间隔约1s 生产系统已经平稳运行了一段时间,获取了足够长的平稳数据本案例中相关数据已进行脱敏处理在建立控制图前首先要基于历史平稳数据对控制图参数的目标值进行估算当所获得历史数据存在异常值时,基于工程经验判断和统计分析,对异常值进行剔除,仅使用平稳数据进行参数的目标值估计将过程平稳状态下获取的平稳数据x应用离散小波变换分解到2层,得到一层近似层小波系数 }和两层细节层小波系数d.,d.a,}和'da.,d.a,} 根据这些小波系数可估计出各 c2.1,c2,2，层小波系数的均值和方差的目标值,计算结果如表B.1所示表B.1各层小波系数均值与方差小波系数均值方差 0.l49" 0,00 0.00 d 0.0038 0.007 5 0.00 现从该过程获取新的1000组数据x=.r1,,r1mm),如图B.1所示下面对这组新的数据实施控制图监测 0.4 0.2 0.0 -0.2 -0.4 一0.6 200 400 600 800 1000 时间/s 图B.1流式数据x 12
GB/T40859一2021 对数据x=(ri,,r1m}应用公式(2),公式(3),并代人j=2,即可使用离散小波变换将数据分解到2层,同样得到一层近似层小波系数(ce,1,cs.,,c?,2)和两层细节层小波系数'd.,d.2,,d,m" 和(ds.1,de.2,,d3.2m} 对于近似层小波系数,使用EwMA控制图监测变量的均值,选择平滑系数入=0.1,控制限系数 L=3 将(cesa,e,es)和从,=0.00代人公式(4)即可得到EwMA控制图的统计量= 将从 0.00和:=0.l497代人公式(6)公式(7)可计算EwMA控制图的控制限将统计量-，和控制限绘制成图,如图B.2所示 o---- 0.2 0.0 -0.2 -0.4 -0.6 -0.8 20o 40o 60o 800 10o0 时间/s 图B.2近似层小波系数EWMA控制图从控制图中可看出,EwMA统计量在第580904秒超出了控制限,可判断在此段时间内流式数据的均值发生了变化对于两层细节层小波系数,分别使用EwMs控制图监测变量的方差,选择平滑系数A;=0.1,控制限系数a=0.0027 将{d.1,d.a,d.o、从小=0.00,、d，=0.0038和d ，d 2,l,d2,2, 2,25o、'd 0.00,.; =0.0075分别代人公式(8),即可得到EwMS控制图的统计量和s 将入,=0.1代人公式(9)、公式(10)可计算EWNMS控制图的控制限将统计量si、sa.和控制限绘制成图,如图B.3和图B.4所示 13
GB/T40859一202 0.12 0.10 0.08 . 0.06 0.04 0.02 200 400 600 800 000 时间/s 图B.3第1层细节层小波系数EWMIS控制图 0.11 0.10 0.09 气10.08 0.07 0.06 0.05 0.04 600 200 400 800 1000 时间/s 图B.4第2层细节层小波系数EwMS控制图从图B.3中可看出,第1层细节层小波系数的EwMS统计量在第158一240秒,第278一298秒第490678秒,第738一786秒和第892962秒内出现了超出控制限的情况;从图B.4中可看出,第2 层细节层小波系数的EwMS统计量在第468500秒,第560600秒和第676一708秒内出现了超出控制限的情况由于两个控制图分别包含原始数据在不同频域上的信息,应将两个控制图的结果结合在一起考虑,最终可判断流式数据在上述时间段内变量的方差发生了变化 14
GB/T40859一2021 录附 C 资料性多变量流式数据应用实例在某直拉单晶硅生产过程中,涉及两个质量相关的变量A和B(该数据已进行数据脱敏),过程中数据的获取间隔是1s 生产系统已经平稳运行了一段时间,获取了足够长的平稳数据X 在建立控制图前首先要对控制图参数的目标值进行估算当所获得历史数据存在异常值时,基于工程经验判断和统计分析,对异常值进行剔除,仅使用平稳数据进行参数的目标值估计将过程平稳状态下获取的平稳数据x的每一个变量都用离散小波变换分解到2层,并重组得到一层近似层小波系数向量(ea1,ca2,}和两层细节层小波系数向量d,.1,d12,}和da1,da.a,} 根据这些小波系数向量可估计出各层小波系数的均值和方差的目标值,近似层小波系数向量均值为M=[0.00,0.00],协方差矩阵为习,=[0.1497,一0.0102;一0.0102,0.0579] 两层细节层小波系数向量均值为M=M [O0.00,0.00],协方差矩阵为=[o.0038,0.0000;0.0000,0.0405],习，=[0.0075，一0.0017; -0.0017,0.0520 现获取了1000组新的数据x={xi,x2,,x1o),其中xi=[zI,r'"] 两个变量如图C.1所示 0.5 0.0 -0.5 200 400 600 800 000 200 400 600 800 000 时间/s 图C.1 变量示意图对于获取的数据X={x1,x2,x10o}进行控制图监测应用公式(2)、公式(3),并代人=2,即可使用小波变换分别将两个变量的数据分解到2层,每个变量都可得到一层近似层小波系数cs,1,cs.a }和'd.a,d..d.m) 将两个变量对应层的 ,es.3o}和两层细节层小波系数'dI,d.a,,d.o 小波系数重新组合,则每一层都包含两个变量的小波系数,得到小波系数向量ea.1ea.a,,ca.25、 d.1,d12,d,.s0)和(da1,da.2 ,d2.,250} 对于近似层小波系数向量,使用MEWMA控制图监测变量的均值,选择平滑系数入，=0.l 将 lea1,e.2,,e;.】和关代人公式(11),即可得到迭代向量2 再将2，和工代人公式(12)公式(13),即可得到MEWMA控制图的统计量Yi 根据表4可知,该控制图的控制限hA=10.08 将统计量Y和控制限hA绘制成图,如图C.2所示 15
GB/T40859一202 35 30 25 20 15 10 200 400 600 800 1000 时间/s 图c.2近似层小波系数MEwMIA控制图从控制图中可看出,MEwMA统计量在第428一656秒,第684806秒和第888904秒内超出了控制限,可判断在此段时间内数据中的某一个或某几个变量发生了均值的变化对于两层细节层小波系数向量,使用MEwMC控制图监测变量的方差和协方差,选择平滑系数入，=0.1 将'd.1,da,,d.s0、Mu,》，和(da.1,da.a,,d.20),uu>分别代人公式(14),即可得到转换后的向量序列U.，和Ua1 再将入A，=0.1、U,和U.代人公式(15、公式16)，即可计算 MEwMC控制图的统计量y,和y 根据表5可知,该控制图的控制限h =0.63 将统计量y y..和控制限h 绘制成图,如图C.3和图C.4所示 2.0 1.5" 1.0 y Ue 0.5 0.0 200 00 600 800 1000 时间/s 图c.3第一层细节层小波系数\IEwMc控制图 16
GB/T40859一2021 1.2 1.0 0.8 y o.6 0.4 0.2 0.0 600 800 1000 200 400 时间/s 图c.4第二层细节层小波系数MIEwMC控制图从图C.3,图C.4中可看出,第一层细节层小波系数的MEwMC统计量在第246250秒,第340 358秒、第478一588秒,第610~654秒,第774~812秒和第956988秒内出现了超出控制限的情况; 第二层细节层小波系数的MEwMC统计量在第560一592秒,第624一680秒和第964一976秒内出现了超出控制限的情况由于两个控制图分别包含原始数据在不同频域上的信息,因此将两个控制图的结果结合在一起考虑,最终可判断原始数据在上述时间段内某一个或某几个变量的方差或几个变量间的协方差发生了变化 17
GB/T40859一2021 附录 D 资料性代码实现 D.1概述为了实现对生产过程的监测,需要使用计算机的辅助且在本文件介绍的算法中,小波变换、矩阵计算等复杂的数学计算往往也可通过编程软件的内置丽数进行 Python是目前应用最广泛的数据处理语言之一,以下是使用Python实现绘制两种控制图的主要代码 D.2单变量流式数据控制图构建算法建立EwMA类 importmath importnumpyasnp matplolib,pyplatas叫t import classEWMA 1I 这是用于建立EwMA控制图的类 definitself，lambda1=0.l，muc=0， sigma-c=1，level=1，L=2.715): 111 初始化EwMA类参数: anmbda_1;平滑系数均值目标值 mu_C 方差目标值 sigmLc: 小波系数的分解层数 level: 控制限系数 ttt self.lambda_1 lambda_1 self.mu_c mu_e self.sigmma_c sigmma_c self.level level C]井记录所有的统计量的值 sel.value self.oldnumber=mu_c井用于表示z(i一1),其中z(0)=mu_c self.new_number=0#用于表示z(i) self.uel=self.muLc十1self.sigma_c”sel.lambda_1/(2-self.lambda_1))0.5井上控制限 (sel.、 self.lcl self. L .sigma_e”self.lambda_1/2-self.lambda1)"0.5井下控 .mu_c 制限 defupdate(sel,data): 18
GB/T40859一2021 获取新的小波系数后计算对应的控制图统计量参数: 新获取的小波系数 data; 1 len(data井新获取的小波系数的总长度 n foriin e(0,n): rangel sel. fnew-number=sdlI.lambdaL_1'dataC门十(sedl.lambda_1)'sdlf.oldLnumber 对应公式(4 self.value. eappend(sel.new- _number self.old_number self.new_number defplot(self、length): I111I 将最新获取的部分数据绘制成控制图参数: ength;绘制控制图的长度原始数据的长度) otaL_length=len(self.value井目前已经计算出的统计量的长度 adjusted_length=int(math.floor(length'2(-self.level))井原始数据的长度对应到该层小波系数的长度 plot_length=min(totaL_length，adjusted_length#井若目前计算出的统计量不足则绘制已计算的统计量 L_label一np.arange(2”sel.level,plot_length'2”self.level十1,2”self.level井该层小波系数对应到原始数据的时间戳 plt.plot(t_label,np.asarray(self.value[plot_length:]).reshape(plot_length,),'k',linewidth =0.6#绘制统计量 plt.plot(t_label,self.ucl”np.ones(plot_length),'k-',linewidth=0.6井绘制上控制限 plt.text(0,sef.ucl,'$\mathregular(U_(CL}$') plt.plot(t_label,self.lcl”np.ones(plot_length),'k',inewidth=0.6#绘制下控制限 plt.text(0,self.lel,'$\mathregular(1_CL$') plt.xlabel'时间(秒)'，fontproperties="SimHe") plt.ylabel('z$_i$',rotation=360,fontsize=16) plt.show(O 建立EwMS类 importmath importnumpyasnp importmatplotlib.pyplotasplt EwMs() class 这是用于建立EwMS控制图的类 19
GB/T40859一2021 def _init_(self,lambda2 =0.1，muLd=0，sgmaLd=1,level=1，uel=1.32,ll=0.67): 初始化EwMS类参数 lamm mbda_2:平滑系数均值目标值 tmmu_C; sigma_c:方差目标值 level: 小波系数的分解层数 ucl: 上控制限参数(标准差的倍数) lcl 下控制限参数(标准差的倍数 self.lambda_2 lambda_2 self.mu d mud self.sigmma_d d SIgma self.level level 井记录所有的统计量的值 self.value [ sigma_d井用于表示s^2(i-1),其中s2(o)=sigma_d self.old_number 井用于表示、2(i) self.new_number 井上控制限 self.ucl=ucl”self.sigma_d"0.5 #下控制限 self.lel lcl'self.sigma_d"0.5 defupdate(self，data): 获取新的小波系数后计算对应的控制图统计量参数新获取的小波系数 data: m 井新获取的小波系数的总长度 n len(data foriinrange(0,n): elf.lambda_2'(data[门-self.mu_d)2十1-self.lambda_2)”sel. sel.new_umber oldLnumber井对应公式(8) self.value.append(math.sqrt(self.new_number))井计算统计量s() self.old_number=self.new_number defplot(self，length): " 将最新获取的部分数据绘制成控制图参数 ength;绘制控制图的长度(原始数据的长度 " totalength=len(sell.value井目前已经计算出的统计量的长度 dustedLIength= int(math.floor(length'2(-self.level))井原始数据的长度对应到该层小波系数的长度 20
GB/T40859一2021 lot_length=min(totaL_length，adusted_ength#若目前计算出的统计量不足则绘制已计算的统计量 Llabel一np.arange(2”selL.level,plotL_length2selL.level+1,2”sel.evel井该层小波系数对应到原始数据的时间戳 gth,),'k',linewidth plt.plot(t_label,np.asarray(self.value[-plot_length]).reshape(plot_lengt 0.6 井绘制统计量 plt.plot(t_label,self.ucl'np.ones(plot_length),'k-'",linewidth=0.6井绘制上控制限 plt.text(O,self.ucl,'$\mathregular(U_(CL)$' plt.plot(t_label,self.Iel np.ones(plot_length)，k--,linewidth=0.6井绘制下控制限 plt.text(0,self.lel,$\mathregular(I_ICL)}$' plt.xlabel'时间(秒'，fontproperties="SimHei") 360,fontsize=16 plt.ylabel('s$_i$',rotation plt.show( 构建控制图假设构建控制图前获取了一组受控数据'data' import numpyasnp Imp0rt matplotlib.pyplotasplt romEwMAimportEwMIA井引用前面定义的EwMA类 romEwMsimportEwMS井引用前面定义的EwMs类 inmportpywt井用于小波变换的程序包井井这里应给定用于计算控制图参数目标值的平稳历史数据,赋值给data变量井 data np.array([datal,data2,] 井也可以应用下行代码生成1000个仿真样本数据进行测试 np.random.normal(0,1,1000) data 井通过受控数据计算控制图参数的目标值 sample_size=len(data 井对应7.4步骤1,将受控数据分解为各层小波系数,设置控制图参数并计算各参数的目标值井分解层数为2 wavelet_coefficient=pywt.wavedec(data,haar',level=J井使用离散小波变换将受控数据分解为各层小波系数井设置EwMA控制图的平滑系数,并计算近似层小波系数的均值,方差的目标值 A2 amda 0.l np.mean(waveletLcoefficient[o],axis=0 m =0 np.cov(wavelet_coeffieent[0],rowvar= sigma #设置第二层EwMS控制图的平滑系数,并计算第二层细节层小波系数的均值、方差的目标值 lamdaD2=0.1 muD2 oefeient[1],axis- 0 np.mean(wavelet_coe D2 np.cov(wavelet_coeffieient[1],rowvar=0) sigmaL 井设置第一层EwMS控制图的平滑系数,并计算第一层细节层小波系数的均值、方差的目标值 lamda_D1=0.1 21
GB/T40859一202 D1 0) rudle-oefarcenm[2]axsi mu np.mean(wave 0 igma_D1 np.cov(wavelet_coeffieient[2],rowvar 井查表2设置EwMA控制限系数L=2.715,创建EwMA类 ontrol_chart_A2=EwMA(lamda_A2,mu_A2,sigma_A2,2,2.715 controlchartA2. .coeffieient[0] .update(wavelet 井查表3设置上控制限为1.50倍标准差、下控制限为0.54倍标准差,创建EwMS类 controlchartD2 EWMSlamdaD2,muD2， D2,2,1.50,0.54 ,sigma control_chartD2. te(wavelet_coefficient[1] .updat 井查表3设置上控制限为1.50倍标准差、下控制限为0,54倍标准差,创建EwMS类 controlchartD1 EWMS(lamdaD1,muD1， ,,sigma_D1,l,l.50,0.54) ontrol_chart_D1.update(wavelet-coefficient[2] controlchart fcontrol_chart_A2control_chartD2,control_chart_D1] 井井数据监测 #这里应给定用于监测的样本数据,赋值给test_data变量井 test_data np.array([datal,data2,] 井也可以用下行代码生成1000个仿真样本数据进行测试 test_data" np.random.normal0,l,l000 len(testdata井用于监测的样本数据的长度 oldLdata=np,.array(]井保存用于计算小波系数的数据井计数为0,该计数用于判断是否更新小波系数 whilenGB/T40859一2021 =0 old_data=np.array([口井清空用于计算小波系数的数据 n=n十1 foriin range(0,J+1): pll.rcParams['igure.dpi]=300 ontrol_chart[].plot(N井分别绘制各个控制图 D.3多变量流式数据控制图构建算法建立MEwMA类 importmath mmportnumpyasnp importmatplotlib. b,pyplot“叫 classMEWMA(): 11 这是用于建立MEwMA控制图的类 " def_init_self,lanmbda_3=0.1，muLc=0，sigma_c=1，level=1，uel=8.64) 初始化MEwMA类参数: ambda_3;平滑系数均值向量目标值 mu_c: 协方差矩阵目标值 sigma_c: 小波系数的分解层数 level: 控制限 ucl: sel.lambda_3 lambda_3 sel.mu_c mu_e self.sigma_c np.mat(sigma_c) self.leve level ]井记录所有的统计量的值 self.value muLc井用于表示Z(i一1),其中Z(o) self.old_number mu_c self.new_number=np.array([#Z(i) self.uel 井控制限 ucl defupdate(sel,data): 获取新的小波系数后计算对应的控制图统计量参数: data;新获取的小波系数向量 1 [n,p data.shape井n为小波系数的总长度,p为变量维数 23
GB/T40859一202 foriin range(O0. ),n: el.ambdha-3》"data[i]十(ld.ambda_3)sell.old走Lnunber= sdlf.new_number 对应公式(11) seL.ambda_3/(2.sedl.ambad_3)'eell.sgmaL- #对应公式q3) sigma_z trans self.new_numberself.mu-e z)),trans.T)) 书 dot( sedlI.sadlue.append(G .do(ras,p.lhmlk.in(sigm np.lloatnp.d np. 对应公式(12)计算统计量并保存 self.oldnumber self. faew-numbes defplot(self，length): 111 将最新获取的部分数据绘制成控制图参数 ength;绘制控制图的长度(原始数据的长度 111 #目前已经计算出的统计量的长度 total_length=len(self.value) adjusted_length=int(math.floor(length*2(-self.level))井原始数据的长度对应到该层小波系数的长度 plot_length=min(totaL_length，adjusted_length#若目前计算出的统计量不足则绘制已计算的统计量 L-label=np.arange(2”self.level,plot_length*2”self.level+1,2self.level)井该层小波系数对应到原始数据的时间截 plt.plot(_label,np,.asarray(self.value[plotL_length;]).reshape(plot_length,)'k',linewidth 0.6井绘制统计量 plt,plot(_label,self.ucl'np.ones(plot_length),'k',inewidth=0.6)井绘制上控制限 plt.text(0,self.ucl,'$\mathregular'U_CL$') plt.xlabel('时间(秒)'，fontproperties="SimHei") plt.,ylabel('$\mathregular(Y_iM2}$',rotation=360,fontsize 16 plt.show( 建立MEwMC类 importmath importnumpyasnp importmatplotlib.pyplotasplt dlassMEwMC(. 转" 这是用于建立MEwMC控制图的类转" def self =0.1， lambd-4 =1，udl-l.169) _init ，mu_d=0，sigma_d=l，leve m 初始化MEwMC类参数 ambda-_4:平滑系数均值向量目标值 tmmu_d: 24
GB/T40859一2021 sigma_d:协方差矩阵目标值 leel 小波系数的分解层数控制限 ucl 一 self.lambda lambda4 _4 self.mud np.matmu_d self. d .sigma_ sigma_d self.level level self.value [ #记录所有的统计量的值 self.old_number=np.eye(Ien(mu_d)井用于表示s(i-1),其中s(0=Ip aray(》井用于表示s(D) self.new_number=np.ar self.ucl ucl 井控制限 nat(sigma_d) eigVals,eigVeets=np.linalg,.eig(np.m foriinrange(0,len(eigVals) egVeets[,]=eieVects[,](a/mp.sr(eigVals[]) self.Amatrix=eigVects.T井计算公式(14)中的A矩阵 defupdate(self，data): 获取新的小波系数后计算对应的控制图统计量参数: data;新获取的小波系数向量 data.shape井n为小波系数的总长度,p为变量维数 [n,p= foriinrange(0,n): trans=data[i,]-self.mu_d np.dot(self.Amatrix,trans.T UUT=np.dot(U,U.T井对应公式(15)中的U(i)"U(i)AT sel.new_number=self.lambda_4UUT十(1-self.lambda_4)”self.oldLnumber#对应公式(15 self.value.append(np.trace(self.new_number)-np,.log(np.linalg,det(self.new_number)) p井对应公式(16)计算统计量并保存 self.old_number=self.new_number defplot(sel，length): 将最新获取的部分数据绘制成控制图参数 ength;绘制控制图的长度原始数据的长度 otaL_length=len(self.value井目前已经计算出的统计量的长度 int(math.floor(le ength2( adljusted_length= -self.level)井原始数据的长度对应到该层小波系数的长度 25
GB/T40859一202 plot_length=min(totaL _ength，adjusted_length井若目前计算出的统计量不足则绘制已计算的统计量 _labe e(2”self.level,plot_length'2”self.level十1,2”self.level井该层小波 =np.arange 系数对应到原始数据的时间戳 plt.plot(t_label,np.asarray(self.value[-plot_length]).reshape(plot_length,),'k',linewidth =0.6井绘制统计量 np.ones(plo_length),'k-',imewidth plt.plot(_label,self.ucl'" =0.6井绘制上控制限 ext(0,self.ucl,$\mathregular(U_(CL$" plt.tex plt.xlabel('时间(秒' ontproperties="SimHei") 360,fontsize 16 plt.ylabel'$\mathr {y_i)$',rotation regular plt.show( 构建控制图: 假设构建控制图前获取了一组p=2维的受控数据,分别记为a和b 1mportnumpyasnp importmatplotlib. pyplatas plt fromMEwMAimmportMEwMA井引用前面定义的MEwMA类 romMEwMcimportMEwMC井引用前面定义的MEwMC类 importpywt井用于小波变换的程序包井这里应给定用于计算控制图参数目标值的平稳历史数据,赋值给data变量 np.array(L[al,a2, np.array([Cb1,b2,] np.vstack((a,b)).T data 也可以应用下行代码生成1000个仿真样本数据进行测试 np.random.normal(0,l,1000) np.random.normalo,1,1000 np.vstack((a,b)).T data 井通过受控数据计算控制图参数的目标值 data.shape[o]井受控数据长度 sample_size data.shape[1]井变量维度井对应8.5步骤1,设置控制图参数并计算各参数的目标值分解层数为2 井 pywt.wavedec(a,'haar',level=J eficdiemE a wavelet_coef wavelet_coeffieent_b=pywt.wavedec(b,"haar',level=J 井对应8.2,分别对各个变量进行小波变换,并将小波系数重组为系数向量 wavelet_coefficient oriinrange(0,len(wavelet_coefficient_a)):井将对应层系数组成矩阵 np.vtack((waveldeL.coefiaien_a[门,wavdletL_oedfeien_b[门) temp p.T temp temp. wavelet-coefficient.apend(temp) 井设置MEwMA控制图的平滑系数,并计算近似层小波系数向量的均值向量,协方差矩阵的目标值 lamda_A2=0.1 26

流式数据监测控制图GB/T40859-2021详解

随着数据时代的到来，数据的收集、处理和分析已经成为企业决策的重要依据。而在这其中，流式数据监测控制图作为质量控制和监测的一个重要工具，受到了越来越多企业的重视。

流式数据监测控制图是指通过实时获取生产过程中所涉及的各项参数数据，并结合统计方法和控制技术，对数据进行分析和判断，从而实现对生产过程的监控和控制。它可以帮助企业及时发现生产过程中的异常情况，并采取相应措施加以纠正，从而保证产品质量。

GB/T40859-2021是我国发布的最新标准，旨在规范流式数据监测控制图的应用，规定了流式数据监测控制图的基本概念、原理、方法和应用等方面的内容。下面我们将从以下几个方面对该标准进行详细解读：

一、基本概念

GB/T40859-2021标准中提到了流式数据监测控制图的基本概念和术语，包括：特性参数、不良、样本、子组、控制线、中心线等。对于这些概念和术语，我们需要深入了解，并在实际应用中准确使用，以保证控制图的有效性。

二、控制图的建立

GB/T40859-2021标准中规定了控制图的建立方法，包括：选择特性参数、确定子组大小、计算中心线和控制限、作图等。在使用控制图之前，我们需要按照标准中的要求进行建立，以确保控制图的可靠性和有效性。

三、控制图的应用

GB/T40859-2021标准中明确规定了控制图的应用范围和方法，包括：判断过程是否稳定、识别过程异常、优化生产过程等。在实际应用中，我们需要根据生产过程的具体情况，选择合适的控制图，并结合标准中的方法进行分析和判断。

四、控制图的维护

GB/T40859-2021标准中还对控制图的维护和管理提出了要求，包括：监测数据质量、更新控制线和中心线等。在使用控制图的过程中，我们需要定期检查和维护控制图，以确保其长期有效性。

总之，流式数据监测控制图是企业质量控制和监测的重要工具，而GB/T40859-2021标准则为其应用提供了明确的指导。我们需要深入理解和掌握该标准的要求和应用方法，并在实际应用中加以运用，才能更好地保证生产过程的质量和稳定性。希望本文对您有所帮助。