基于变量选择的工业过程故障诊断方法研究

《基于变量选择的工业过程故障诊断方法研究》由会员上传分享，免费在线阅读，更多相关内容在学术论文-天天文库。

分类号：TP277UDC：621.3学号：15451186195密级：公开温州大学硕士学位论文基于变量选择的工业过程故障诊断方法研究作者姓名：张申波学科、专业：计算机应用技术研究方向：电力电子装备信息化指导教师：舒亮副教授指导教师：闫正兵讲师完成日期：2018年5月温州大学学位委员会 分类号：TP277UDC：621.3学号：15451186195密级：公开温州大学硕士学位论文基于变量选择的工业过程故障诊断方法研究作者姓名：学科、专业：计算机应用技术研究方向：电力电子装备信息化指导教师：指导教师：完成日期：2018年5月温州大学学位委员会 温州大学学位论文独创性声明本人郑重声明：所呈交的学位论文是我个人在导师指导下进行的研究工作及取得的研究成果。尽我所知，除了文中特别加以标注和致谢的地方外，论文中不包含其他人已经发表或撰写过的研究成果，也不包含为获得温州大学或其它教育机构的学位或证书而使用过的材料。与我一同工作的同志对本研究所做的任何贡献均已在论文中作了明确的说明并表示了谢意。论文作者签名：日期：年月日温州大学学位论文使用授权声明本人完全了解温州大学关于收集、保存、使用学位论文的规定，即：学校有权保留并向国家有关部门或机构送交论文的复印件和电子版，允许论文被查阅和借阅。本人授权温州大学可以将本学位论文的全部或部分内容编入有关数据库进行检索，可以采用影印、缩印或扫描等复制手段保存和汇编本学位论文。本人在导师指导下完成的论文成果，知识产权归属温州大学。保密论文在解密后遵守此规定。论文作者签名：导师签名：日期：年月日日期：年月日 基于变量选择的工业过程故障诊断方法研究摘要故障诊断技术可以及时识别和定位出工业过程中的故障，提高工业生产的安全性与可靠性。现代工业过程自动化与智能化程度的日渐提高，产生并存储了大量的过程数据，而这些数据包含了许多有用的过程信息，因此数据驱动的故障诊断方法受到了广泛的关注。本文以变量选择技术为基础，针对实际工业过程中故障方向未知、数据具有随机性和不确定性等问题，进行了相应的研究。具体工作和创新如下：（1）将难以求解的故障重构问题转化为易于求解的变量选择问题，可以在故障方向未知的情况下，同时计算出故障的大小和方向。（2）实际工业数据通常具有随机性和不确定性，为了提高故障诊断的鲁棒性，本文通过构造预测矩阵，将故障诊断问题转化为参数回归问题，并结合贝叶斯理论与变量选择算法，计算出所有变量发生故障的概率，然后利用故障概率找出过程中最有可能的故障变量。最后，实验研究表明了本文所提工业过程故障诊断方法的有效性。本文方法不仅可以准确地识别造成故障发生的变量，找出故障发生的根本原因，还可以给出变量发生故障的概率，为接下来的过程恢复提供更多的有用信息。关键词：故障诊断，故障重构，变量选择，贝叶斯理论I II ResearchonFaultDiagnosisMethodofIndustrialProcessBasedonVariableSelectionABSTRACTFaultdiagnosistechnologyisabletoidentifyandlocatethefaultsintheindustrialprocess,andimprovethesafetyandreliabilityofindustrialprocess.Withtheincreasingdegreeofautomationandintellectualizationinmodernindustrialprocess,alargeamountofdatathatcontainsusefulprocessinformationiscollectedandobtained.Therefore,thedata-drivenfaultdiagnosismethodshasreceivedextensiveattention.Theresearchofthisthesis,basedonthevariableselectiontechnique,aimstosolvetheproblems,whichincludeunknowndirectionsofthefault,therandomnessanduncertaintyofthedataintheactualindustrialprocess.Thepracticalworkandinnovationsareasfollows:(1).Theproblemoffaultreconstruction,whichisdifficulttobesolved,istransformedintoaneasiervariableselectionproblem.Withoutfaultydirections,thesizeanddirectionsofthefaultcanbecalculatedsimultaneously.(2).Thepracticalindustrialdatausuallyhasthecharacteristicofrandomnessanduncertainty.Inordertoimprovetherobustnessoffaultdiagnosis,thefaultdiagnosisproblemistransformedintoparametricregressionproblembyconstructingpredictionmatrix.CombiningwiththeBayesiantheoryandthevariableselectionalgorithm,thefaultyprobabilitiesofthevariablescanbecalculated.Withthefaultyprobability,thevariableswhicharemostlikelytobethefaultyvariablescanbeidentified.Theexperimentalstudiesshowtheeffectivenessoftheproposedfaultdiagnosismethod.TheproposedmethodscannotonlyaccuratelyidentifyIII thefaultyvariables,butalsogivethefaultyprobabilityofthevariables,providingmoreusefulinformationfortheprocessrecovery.KEYWORDS:faultdiagnosis,faultreconstruction,variableselection,BayesiantheoryIV 目录摘要....................................................................................................................................IABSTRACT.....................................................................................................................III第一章绪论......................................................................................................................11.1研究背景及意义..................................................................................................11.2过程监控的概述..................................................................................................11.3国内外工业过程故障诊断的研究概况...............................................................21.3.1诊断方法..................................................................................................21.3.2数据驱动的故障诊断方法研究现状........................................................41.4本文的主要研究内容和章节安排.......................................................................6第二章基本原理...............................................................................................................72.1变量选择.............................................................................................................72.1.1变量选择的定义......................................................................................72.1.2变量选择问题求解...................................................................................82.2PCA监测统计量.................................................................................................92.3蒙特卡洛方法的概述........................................................................................102.3.1蒙特卡洛数值积分.................................................................................102.3.2蒙特卡洛方法的基本思路.....................................................................102.3.3Gibbs抽样.............................................................................................112.4核密度估计.......................................................................................................112.5本章小结...........................................................................................................13第三章基于变量选择的故障重构方法..........................................................................153.1引言...................................................................................................................153.2故障重构...........................................................................................................153.3Lasso重构.........................................................................................................163.3.1求解故障方向Ξ...................................................................................173.3.2参数μ的选择........................................................................................183.3.3计算故障幅值e......................................................................................183.3.4故障诊断................................................................................................193.4实例研究...........................................................................................................203.4.1数值案例................................................................................................203.4.2TEBenchmark实例...............................................................................213.5本章小结...........................................................................................................29第四章基于参数回归的工业过程故障诊断方法...........................................................314.1引言...................................................................................................................314.2参数回归...........................................................................................................314.3BayesianLasso故障诊断..................................................................................324.3.1参数后验密度函数的推导.....................................................................344.3.2Gibbs抽样的构建..................................................................................354.3.3回归系数β控制限的计算......................................................................364.3.4变量故障概率的计算.............................................................................36V 4.3.5实例研究................................................................................................374.4BayesianFusedLasso故障诊断........................................................................424.4.1实例研究................................................................................................464.5本章小结...........................................................................................................53第五章总结与展望...........................................................................................................555.1总结...................................................................................................................555.2展望...................................................................................................................55参考文献............................................................................................................................57致谢....................................................................................................................................61攻读硕士学位期间发表的论文.........................................................................................63VI 第一章绪论第一章绪论1.1研究背景及意义随着自动化技术和计算机技术的不断发展与进步，现代工业过程的复杂度、自动化和智能化程度变得越来越高，过程中的耦合也变得越来越强。因此，工业过程中的一些微小故障也可能在过程中传播和扩大。如果不及时排除，可能导致整个生产过程不能正常的运行，甚至可能引发安全事故，造成难以估计的后果。例如：1979年，美国三里岛发生的核泄漏事故，不仅造成放射物质外流，严重影响居民的日常生活，而且还造成了高达数十亿美元的经济损失；1984年印度博帕尔农药厂的化工有毒气体泄漏，造成了3500多人死亡，近百亿美元的经济损，博帕尔也因此而成为一座恐怖之城；1984年，前苏联的切尔诺贝利核电站核事故造成了2000多人死亡，直接的经济损失高达30亿美元，并对自然环境产生了巨大的影响；2011年日本福岛核电站的事故，使二十公里外受到严重辐射、附近所有的居民撤离、周边地区的蔬菜与鱼类的核辐射指标严重超标。在国内，工业事故也数见不鲜。例如：1996年的天津津西大华化厂爆炸事件造成19人死亡、14人受伤，直接的经济损失达120多万元；2013年，中石油天然气公司大连石化分公司的三苯罐区爆炸事故，则造成了4人死亡，直接经济损失大约697万元；2016年，万华化工集团有限公司的烟台工业园区“9.20”爆炸事件，造成4人死亡、4人受伤，直接的经济损失约573万元；2017年，河北利兴特种橡胶股份有限公司的氯气泄漏事件，造成公司的现场员工以及附近人员中毒，紧急疏散周边群众近1000人，还造成2人死亡与25人住院。因此，提高工业过程的可靠性以及稳定性，增强工业过程中的实时监控，具有重要的实际意义。与此同时，现代工业过程因其广泛采用了集散控制系统和现场总线技术，产生并存储了大量的过程数据。在这些数据中，包含有能够反映工业过程运行状态的有用信息，而这信息对于过程监控非常重要。因此，如何有效地利用这些过程数据来监控工业过程已成为目前关键的研究方向[1]。1.2过程监控的概述过程监控主要包括四个步骤：故障检测、故障诊断、故障根源分析与过程恢复[2]。在进行监控之前一般先收集工业过程数据，并进行相应的预处理，所以一1 温州大学硕士学位论文个完整的过程监控系统大致如图1-1所示。由于大多数现代工业过程运行环境苛刻，控制策略复杂，传统单变量故障诊断方法难以进行有效地故障诊断。同时，在过程监控领域，对不同类型的故障检测方法进行了较为深入的研究。相比之下，因难以分析多变量对监测统计指标的影响，针对故障诊断的研究相对较少。所以，深入地进行故障诊断方法研究具有实际的应用价值。故障根源工业过程数据采集数据预处理故障检测故障诊断过程恢复分析图1-1过程监控示意图Figure1-1Theschematicdiagramofprocessmonitoring在现代工业过程中，故障一般指的是一个过程或者系统的特征（变量）出现了不同于正常状态的偏差[3]。工业过程中的故障类型主要有：过程的参数变化、干扰参数变化、执行器问题及传感器问题[4]。过程监控中故障检测用来找出发生故障的样本，确认系统是否有故障发生。故障诊断的任务是在故障样本中找出发生故障的观测变量。故障根源分析则用以确定过程中发生的故障类型、准确地找出导致故障发生的变量，为接下来的过程恢复提供支撑。过程恢复则利用诊断得到的信息使过程或者系统恢复到正常状态，消除故障对其的影响。1.3国内外工业过程故障诊断的研究概况1.3.1诊断方法目前，已经有了很多不同类型的故障诊断方法，学者和专家从不同的角度将这些故障诊断方法进行了不同的分类。这其中的Frank[5]分类方法得到学术界的普遍认同。他将故障诊断的方法大致分成三类：基于解析模型[3]、基于数据驱动[6]以及基于知识[7-8]的方法。图1-2给出了一些常用的故障诊断方法。其中，解析模型诊断方法属于一种定量的分析方法。基于构建的系统数学模型，该方法通过参数模型估计以及等价空间方程来提取残差。然后，通过设定的特定阈值对其进行评价和分析，以判断系统是否发生了故障，并找出故障发生的原因[9]。若残差信号值的偏离程度越大，则表示系统发生的故障越严重，反之，则认为故障较小，系统受到了扰动因素的影响。这类诊断方法主要包括参数估计法、状态估计法与等价空间法三种类型[3]。因该方法可以在那些被控对象是输入和输出较少的过程（系统）中建立准确的数学模型，实时深入地分析系统的动态特性，所以能够对这些过程进行有效的故障诊断。但现代复杂的工业过程大多是输入输出较多、参数不定的系统，因而基于解析模型的方法难以对其建立准确的2 第一章绪论数学模型。此外，该方法因过于依赖建立的数学模型，导致建模时产生的误差、噪声及参数的扰动对诊断结果的影响极大，得到的诊断结果准确性不高。因此，为了解决这些问题，近年来该方法的研究主要集中于抗干扰问题的解决。参数估计方法基于解析模型的方法状态估计方法等价空间方法模糊推理方法故神经网络方法障诊基于知识的方法断方专家系统方法法定型模型分析方法部分最小二乘方法费舍尔判别分析法主成分分析方法基于数据驱动的方法支持向量机方法基于重构的方法基于变量选择的方法图1-2故障诊断方法Figure1-2Classificationoffaultdiagnosisapproaches随着工业过程控制系统的复杂化和大型化，精准的数学模型也变得越来越难以有效建立，这使得基于解析模型的方法在实际工业过程中难以充分发挥作用。基于知识的方法具有易于模块化、修改的优点[10]。故障诊断时，基于知识的方法既不需建立系统模型，又能利用被控对象的信息建立故障检测规则特性，因而能够对一些过程(系统)进行及时的故障检测。该方法主要适用于非线性、复杂控制系统和过程的故障诊断，受到了越来越多研究者们的关注[11]。并且在诊断过程中，该方法主要通过启发式症状、过程历史知识和相应的统计知识对故障进行诊断。3 温州大学硕士学位论文这类方法可以分成为神经网络[12]、专家系统[13]、定性模型分析[14]、模糊推理[15]和故障树[16]等方法。同时，基于知识的诊断方法大多是基于人工智能原理，并借助计算机模拟人的思维过程，对系统进行相应的推理和诊断。而人类的智能技术则相对复杂，因此要使计算机具备和人相近的推理能力还需要一定的时间。与前两种方法不同，基于数据驱动的故障诊断方法则不需要完整的系统先验知识。该方法在获得过程产生的数据基础上，利用数据处理和统计分析方法挖掘数据中的隐含信息来进行故障诊断。其基本思想是先利用过程历史数据建立数据预测模型，然后通过比较预测模型和系统实时数据建立的模型来实现故障诊断。通常利用降维技术将高维数据投影到低维的特征空间，以得到新的、便于常规处理的特征集数据。因此，基于数据驱动的方法既可以最大限度地将冗余信息去除，保留有用的原始数据集信息，又能减少计算的复杂度。这类方法适用于那些解析模型难以建立的复杂工业过程。主要基于数据驱动的方法有主成分分析（PrincipalComponentAnalysis,PCA）[17,18]、支持向量机（SupportVectorMachine,SVM）[19]、部分最小二乘（PartialLeastSquares，PLS）[20]、费舍尔判别分析（Fisherdiscriminantanalysis,FDA）[21]、基于重构的方法[22]和基于变量选择的方法[23]。因工业过程的复杂性及难以建模的特性，前面两类方法难以对其进行有效的故障诊断。而传感器技术与计算机技术的快速发展为实时的数据获取与分析提供了可能。所以，数据驱动的方法能够有效地诊断出工业过程中的故障，为现代工业过程系统的稳定运行提供有力的保障。本课题选择此类方法进行故障诊断研究，接下来将阐述这类方法的大致发展概况。1.3.2数据驱动的故障诊断方法研究现状现代工业过程的复杂度和自动化程度越来越高，过程中的耦合也变得越来越强，产生和存储了大量的过程数据。然而，解析模型的故障诊断方法和基于过程知识的故障诊断方法均难以对复杂工业过程进行有效的故障诊断。因此，目前的现代工业过程故障诊断大多采用数据驱动的故障诊断方法。在这些数据驱动故障诊断方法中，主成分分析法（PCA）[24-26]在工业过程的故障诊断中被广泛采用，特别是在大规模复杂工业过程的故障诊断研究领域[27,28]。因其可以通过变换从高维数据中提取包含整个数据信息的特征集，所以PCA具有计算简单和易于执行的优点。在故障诊断过程中，PCA利用HotellingT2和预测误差（SPE）[29]两个统计监测指标来判别工业过程是否发生故障。但是，PCA也存在一些难以克服的缺陷。一方面，当诊断过程中获得的数据是一些高度非线性、非高斯特性的过程数据时，PCA将无法提取包含整体特性的特征集，诊断结果的可靠性将会受到影响。另一方面，当工业过程包含多个变量时，故障诊断难以设定一个符合多个变量特性的控制限，诊断结果的误报率会因此而提高。而且，故障诊断的结果也会受到拖尾效应的影响[30]。因此，该方法在故障是多元变量、4 第一章绪论数据是非高斯特性的工业过程中难以充分地发挥作用。为了能够更好地将PCA运用到故障诊断领域中，国内外学者根据不同的实际情况对PCA进行一些相应的改进。例如：郭明提出的一种基于特征子空间故障检测与识别方法，该方法取得的诊断结果明显要比传统PCA方法的结果好[31]。之后，因为传统的PCA方法大多是针对静态数据，而现代工业过程大多具有动态变化的特性。所以，HuangJ.和YanX.[32]提出一种动态PCA方法（DPCA），利用时滞矩阵有效地解决了动态数据的观测问题。该方法不仅实现了对定动态过程的故障诊断，还能分析数据之间的相关特性。随后，LiL.和YuY[33].等人提出的多元PCA（MPCA）方法，为工业过程提供实时的过程监控，取得了较为理想的诊断结果。而针对实际工业过程中的非线性问题，XuY.、LiuY.和ZhuQ.[34]等人提出了一种基于KPCA的故障诊断方法。该方法能够利用核函数对非线性主元进行巧妙地计算，实现PCA在非线性诊断领域中的应用。虽然这些改进的PCA方法在一定程度上拓展了其在工业工程故障诊断中的应用，但是PCA的固有缺陷仍未被很好的克服。例如，因在故障诊断过程中使用了与变量有关的二阶统计量，在对非高斯的工业过程进行相应的故障诊断时，得到不好的诊断的效果；拖尾效应对诊断结果的影响同样也未被很好的解决。除了基于PCA分析的故障诊断方法，基于重构的故障诊断方法也被广泛应用于故障诊断领域。故障重构方法假设潜在故障方向是已知的，然后利用故障方向进行重构，找出数据样本中发生故障的变量[22]。然而，在实际的工业过程中，故障方向往往难以获得。YueH.H.和QinS.J.[35]提出了一种改进的故障重构方法，该方法可以从大量的历史故障数据中提取故障方向，但是在实际工业过程中，故障数据较少，很难获取足够的历史故障数据，因而限制了其在实际过程中的应用。在假设故障方向已知的条件下，AlcalaC.F.和QinS.J.[36]提出的重构贡献（Reconstruction-basedContribution,RBC）故障诊断方法，通过最值优化的求解方式得到故障幅值，然后计算出重构后的检测统计量，将那些重构统计指标小于其控制限的变量判别为故障变量。但该方法仍然未有效地解决拖尾效应问题。在随后的故障诊断研究中，KariwalaV.etal.[37]和HeB.etal.[38]通过将分支定界（BranchandBound,BAB）与缺失数据分析技术结合方式有效地解决了传统基于重构方法存在的不足，但这种方法的计算复杂度高，尤其是工业过程包含大量变量时，其相应的计算时间呈现指数增长趋势，所以BAB不适合实时的工业过程在线故障诊断。除了上述方法，费舍尔判别分析（FDA）、支持向量机（SVM）等模式分类技术也被用到故障诊断过程中[2,39-41]。这类方法通过比较过程中的故障数据和已知历史过程数据之间相似度的方式来确定过程故障的类型。因而，为了实现准确的故障诊断，这类方法在模型训练时需要足够的故障数据。但在实际过程中，包含完备故障类型的数据集是难以获取的，因而故障分类的方法在实际工业过程中5 温州大学硕士学位论文也难以充分发挥作用。1.4本文的主要研究内容和章节安排本文的主要研究内容如下：第一章介绍本课题研究的背景和意义、过程监控的主要内容与故障诊断方法大致的分类，最后综述国内外工业过程故障诊断领域的研究现状。第二章首先介绍变量选择方法的基本原理，然后简要地介绍了PCA监测统计量、蒙特卡洛方法和核密度估计的相关原理。第三章首先阐述故障重构的基本思路，然后针对工业过程中故障方向未知的问题，提出基于变量选择的故障重构方法。该方法可以有效地利用取得的故障方向进行故障重构，并找出最有可能的故障变量。随后详细地阐述了基于变量选择重构的基本原理与故障诊断的步骤。在实例研究中，简要地介绍了田纳西过程（TennesseeEastman,TE）的基本概况。最后利用数值案例数据和TE过程的实例数据对该方法进行故障诊断准确性以及有效性的验证。第四章介绍基于参数回归的工业过程故障诊断方法，阐述了该方法进行故障诊断的基本原理，并详细地介绍了求解变量故障概率的步骤，然后利用TE过程实例对该方法的故障诊断有效性进行了验证，最后利用交叉验证方法比较两种不同参数回归方法故障诊断的鲁棒性。第五章对论文的研究工作进行简要地总结，并对今后的研究工作进行展望。6 第二章基本原理第二章基本原理2.1变量选择2.1.1变量选择的定义假如提取到的观测数据是{(,yin);x1,2,i,}，则y与xi(,xx,,)xTii12iid分别是其相对应的响应变量与自变量。为了能够更好地进行数据分析，对取得的观测数据进行相应标准化处理：nnn2yixijxnj0,0,d(1,2,)ij(2-1)iii111考虑下面的线性回归模型：yxβξ(2-2)T其中y(,yy)y是响应变量矩阵，x(,,xxx,)是dn预测矩阵，β是d12n12n2维回归系数矩阵，ξ是d维误差矩阵，服从正态分布N(,0I)。若观测值之间nnn互相独立，那么回归模型的最小二乘估计可以通过最小化残差平方和得到：2βargmin||||yxβ(2-3)β形如式（2-3）最小二乘回归解的问题就是非零变量的选择问题。虽然最小二乘估计属于无偏估计，但是当自变量的线性相关度比较高时，得到的预测值方差通常较大。若利用预测值方差的有偏估计近似代替其无偏估计方式，就能使一些回归系数压缩为0，提高模型的预测精度。为了更好地进行回归分析，RobertTbshirani[42]提出一种新的变量选择方法（Leastabsoluteshrinkageandselectionoperator,Lasso）。该方法利用取得模型系数的绝对值函数来压缩模型的系数，使得其中的一些回归系数减小，甚至还可以使一些较小系数的绝对值变为0。一般地，对于式（2-2）的回归模型，其对应的Lasso估计为：d2βargmin||yxβ||,..stjtt(0)(2-4)βj1其中t是调和参数，通过t的调节作用能够使回归系数β的值总体变小。比如当tt00t||,βj时，Lasso就可以使一些较小的回归系数趋于0，甚至直接变为j0。根据拉格朗日乘子法[43]，式（2-4）可以描述成：d2βargminyxβj(2-5)βj1即:7 温州大学硕士学位论文2βargminyxββ1(2-6)β式中χxβ中的表示1-范数。12.1.2变量选择问题求解变量选择问题的实质是一个带约束的二次规划问题。因最小角回归（LeastAngleRegression,LARS）[44]在求解变量选择问题时，计算量小、计算速度快。所以本文在求解式（2-6）的Lasso组合优化问题时，采用最小角回归LARS方法。LARS的核心思想是找出所选变量的新路径。从而，当在这条新的路径上前进时，当前的残差与所选变量之间的相关系数相同，直至找出具有最大相关系数的新变量。从几何意义上讲，在被选入回归集的那些变量构成空间中，当前残差的投影就是该空间的平分线。假如x列向量xx,,x,线性无关，那么定义一个矩阵x(,,)sx，12nAjjjATT11/2其中的sj1，A是1,2,,n的一个子集。定义GxxAAA和A1G1AAA()A，其中的1表示全1向量，长度和A中的元素个数一样。那么在x角平分线方向AA1上的单位向量χ可以表示成：χxω,ωAG1，使χ和每一个xj都有AAAAAAAAAT2相同的角度（小于90度），且xχA1,1χ。假设当前步骤下LARS的AAAAA预测结果是χ，则LARS接下来的具体过程如下：A步骤1：初始化χ0,0β；A步骤2：计算当前的相关系数c：Tcc()χ()xyχA(2-7)步骤3：计算最大相关度和对应的下标集合Cmax{|cjj|},A{:|jc|C}，其中集合A是其中拥有最大（绝对值）相关系数的维度的标号集合。步骤4：计算：xAjjjA(,,)sx(2-8)T11/2TA1G1AAAA(A)A,GAxx(2-9)1χAxAωA,ωAAG1AAA(2-10)TCcCcjj步骤5：更新：χAAχγχA，其中axχA，γmin,jAC，AAaaAjAjCcCcjj新的激活集AAj{},=argminjjAC,，最大相关系数AAaaAjAjCCγA，“min+”表示取最小的正值。A8 第二章基本原理步骤6：更新β：ββγδ，其中δ是ω通过填充0AAA(,j()0)Ajδ，将其扩展成和数据个数一样维度得到的。A步骤7：如果β满足约束，转到步骤4。2.2PCA监测统计量主成分分析（PrincipalComponentAnalysis,PCA）[45]的实质是通过正交变换将系统方差集中分布到少数的主元方向上，变换前后的系统方差并未改变，而只是进行了重新分布。在对数据进行相应主成分分析时，PCA将原本的测量空间分解成为2个正交的子空间，其中的一个是主元子空间（PCS），另一个则是残差子空间（RS）。假设x是dn的观测数据，其中的n是变量个数，d是采样点数。那么，利用PCA可以将x分解成PCS中的x和RS中的x两个部分。即：xxx(2-11)T对得到的x相关系数矩阵Rxx/(1)d进行相应的奇异值分解，可以取得得分矩阵T与载荷矩阵P。R的奇异值分解如下：TRPΛPTxP,(2-12)其中Λdiag{,,,}，是x载荷矩阵P方向上的方差，k是主元的个数。12ki同样地，将载荷矩阵P和得分矩阵T分别在PCS和RS中进行分解，可以得到TTT[]和PPP[]。所以最终可以将x分解为：TTxPPx()IPPx(2-13)TT其中PPx和IPP分别是x在PCS和RS中的投影。基于PCA的T2统计量与SPE统计量常用于判断工业过程是否发生故障。其中，T2用来衡量PCS中样本点投影的变化，SPE用来衡量RS中样本点投影的变化。T2统计量为T21T2TxPΛPx(2-14)T222式中，是T的控制限。若x的分布是多元正态分布，那么控制限为TT2kn(1)Fknka(,)(2-15)T()nk其中，a是置信水平，F[46]a(k,n-k)是自由度为k，n-k的F分布临界值。SPE统计量为2T2SPEIPP(2-16)2[47]式中的是SPE的控制限，可以由下式计算得到9 温州大学硕士学位论文1Ch2Chh(1)h021a02a00(2-17)1211式中：niij,=1,2,3j(2-18)jk1213h012(2-19)32式（2-17）中的Ca是置信水平为α下正态分布的临界值。在对工业过程进行故障检测时，若计算得到的统计量T2或SPE超过了其对应的控制限，则表示工业过程中存在异常或故障。2.3蒙特卡洛方法的概述理论上，贝叶斯推理是比较容易实现的。即，对于任何的先验分布，只需要计算其相应的后验分布，如后验分布的矩和后验概率密度函数等就能够实现贝叶斯分析。但是实际上，未知参数的后验分布大多是高维且复杂的罕见分布，所以对其进行相应的积分计算则变得非常困难，这在一定程度上限制了贝叶斯分析法在实际工业过程中的应用。马尔科夫链蒙特卡洛（MarkovChainMonteCarlo,MCMC）[48]方法将一个数值积分转换成对应迭代抽样的求和，并通过模拟抽样实现高维积分的计算。在抽样次数足够大时，MCMC的估计值则会收敛于其对应的真实值。所以，MCMC能够有效地解决复杂且高维积分的计算问题，并为贝叶斯分析法在实际工业过程中的应用提供了一种新的思路。2.3.1蒙特卡洛数值积分对于任意的函数g(x)，蒙特卡洛积分的表达式为：bgx()qxdz()(2-20)aqx()gx()q(x)是x在区间[a,b]中的概率分布。令fx()，取q(x)中n个样本点，m为qx()平稳时的样本数，当n足够大时，式（2-20）可以用均值近似代替：n1(t)Efxfx(2-21)nmtm12.3.2蒙特卡洛方法的基本思路MCMC是使用蒙特卡洛数值积分实现积分运算的。其基本的思路是：首先构造一条Markov链，使它的平稳分布与待估参数的后验分布一致；然后，通过10 第二章基本原理获得的Markov链生成后验分布的样本；最后利用符合平稳分布的样本进行积分计算。假如R是一空间，n是产生的样本总数，m是链条达到平稳状态时的样本数，那么MCMC的执行步骤可以概括为：（1）构造Markov链。构造一条收敛于平稳分布π(x)的Markov链；（2）产生后验分布的样本。从R中的某一点x(0)起，用（1）中取得的Markov链进行模拟的抽样，生产一个点序列：x(1)，…，x(n)；（3）计算蒙特卡洛数值积分。任一个函数f(x)期望的估计形式为：n1(t)Efxfx。nmtm12.3.3Gibbs抽样Gibbs抽样[49]是在多变量的联合分布未知，而单变量对应的条件分布已知时，根据单变量的条件分布抽取单变量的样本序列，组成一个多变量样本的抽样方法。那么取得的模拟抽样序列符合马尔科夫链特性。在细致平稳状态条件[50]下，构造得到的Markov平稳分布就是所求的多变量联合后验分布。假设给定参数向量(0)(0)(0)(0)x=,,xx12,xn的初始向量xxx=x,,,12n，d时刻的Gibbs抽样值为xx()dxd()dx=d,(),,()，那么d+1时刻的Gibbs抽样值x(d+1)可以通过以下步骤得12n到：()dd()d()(1)d（1）从xx|,,,xx中抽取x；123n1(1)dd()d()(1)d（2）从xx|,,,xx中抽取x；213n2(1)dd(d1)(1)(1)d（n）从xx|,,,xx中抽取x。nn131n最后，利用抽样得到的后验样本来计算后验分布的矩以及概率密度函数，实现预期的贝叶斯分析目的。2.4核密度估计核密度估计（KernelDensityEstimation,KDE）[51]是一种非参数检验方法，常用于估计概率论中未知的概率密度函数。KDE通常使用平滑的峰值函数(“内核”)来拟合测量数据点并模拟真实的概率分布曲线。假如xx,,x,是一个独立12n同分布的F序列，它的概率密度函数为f(x)，所以其核密度估计为：nn11xxifh()xKxxh(i)Kh(2-22)nii11nhh其中，K()是核函数（非负、积分是1，并且符合概率密度的性质、均值是0）。h>0是平滑参数，称作带宽(bandwidth)或窗口[52]。由于核密度估计是一种从数据11 温州大学硕士学位论文本身研究其分布特征的方法，因此它不需要利用数据分布的先验知识，也不需要对数据分布进行任何额外的假设。这样KDE就能够解决预测的参数模型与实际物理模型之间存在较大差距的问题。从式（2-22）能够直观地发现K()与h的选取决定了核密度估计结果的好坏。一般地，核函数对核密度估计的影响要远小于h。而且不同带宽对核密度估计结果的影响差异很大。因此，在选择h时，h的值应能够尽可能地最小化误差。若h太大，用于计算的点太多，虽可以减小方差，但偏倚会比较大；h太小，偏倚减小了，但是用于计算的点太少，方差又很大。因此，在理论层面上，存在一个可以最小化均方误差的h。h的选择取决于n的大小，当n比较大的时候，可以选用一个较小的h。因为即使是一个较小的h，一个较大的n可以保证区间内有足够的数据样本点用于估计概率密度。所以，通常当N→∞时，应使h→0。例1/5如，可以使最优的h赋值为n的-1/5次方乘以一个常数，即hcn。若核函数选择高斯分布，可以计算出c=1.05×标准差。43.532.52probability1.510.50-4-3-2-10123图2-1x的直方图概率密度Figure2-1Histogramprobabilitydensityofx对于一个直方图概率密度（如图2-1），可以利用像直方图一样的计算方式估计x=处的密度函数值：用这个区间里点的个数除以点的总个数得到的值近似代替x=处的密度函数值。即：xxhxhnnxx1#i[,]111if(x)1(xhxixh)11(2-23)2hn2nhii11nh2h1式中分子表示的是落在区间[x-h,x+h]里点数。令Kx()1(x1)，则f(x)的核02估计可以表示为：n1()xxifh()xK0(2-24)nhi1h12 第二章基本原理相应的密度函数积分为：nn11xxifxdx()()dtKdx()KtKtdt000(2-25)nhhiin11因此，只需K()的积分等于1，就可以使得到的密度函数积分为1。常用的四个核函数K()是：1UniformKx:()x1(2-26)2Triangular:()=1-Kxxx1(2-27)ax1GammaKx:(|a,b)x(a1)eb(x0)(2-28)ab(a)2(xxc)12GaussianKxxc:(,)e2(2-29)2其中，Uniform、Triangular、Gamma、Gaussian分别代表均值、三角、伽马以及高斯核函数，xc是核函数的中心，为宽度参数。2.5本章小结本章首先介绍变量选择方法的基本原理以及求解变量选择问题的方法，然后介绍了PCA检测时常用T2和SPE统计量的基本定义、随后概述了蒙特卡洛方法和核密度估计的相关原理，为后续的章节做准备。13 温州大学硕士学位论文14 第三章基于变量选择的故障重构方法第三章基于变量选择的故障重构方法3.1引言故障重构技术能定位故障发生的部位及大小，指导操作人员对生产过程进行干预，使过程尽快恢复到正常状态，提高生产过程的安全性和可靠性。但是，故障重构方法需要知道故障方向，而在实际工业过程中，故障方向往往难以获得。虽然YueHH.与QinSJ.[35]提出的改进重构方法能够从历史故障数据中提取故障方向，但是该方法需要大量历史故障数据，而在实际工业过程通常难以获取足够的历史故障数据；AlcalaCF.和QinSJ.[36]提出的重构贡献（RBC）故障诊断方法虽可以一定程度上弥补传统重构方法的不足，但RBC仍未有效地解决拖尾效应问题；KariwalaV.etal.[37]和HeB.etal.[38]通过将分支定界（BranchandBound,BAB）与缺失数据分析技结合方式解决了故障方向求解的问题，但这种方法的计算复杂度高，尤其当工业过程包含大量变量时，BAB的计算时间呈现指数增长趋势，所以BAB仍然不适合实时的在线故障诊断。因工业过程中的故障通常与一些故障变量有关，所以若能够找出对应的故障变量就能实现准确的故障诊断。最小绝对值收敛和选择算子（Leastabsoluteshrinkageandselectionoperator,Lasso）[42]属于一种变量选择方法，该方法将模型系数的绝对值函数作为惩罚项来压缩模型系数，使其中的一些回归系数变小，甚至能够使一些较小系数的绝对值变为0，以实现变量选择的目的。而且，Lasso在变量选择时，可以通过正则化子项的调节作用避免过拟合现象的出现，从而Lasso可以有效地解决最小二乘回归的稀疏解问题。因此Lasso具有变量选择速度快，计算复杂度小的优点，用Lasso来提取故障方向将会大幅地减少计算的复杂度。受此启发，本章提出一种基于变量选择的故障重构方法，将难以求解的组合优化问题转化成为易于出求解的二次型问题；在此基础上增加一项约束项，将其等价转化为Lasso回归问题，并通过最小角回归算法（LeastAngleRegression，LARS）[44]求解Lasso回归问题以获取未知的故障方向，随后利用取得的故障方向进行故障重构，快速地定位出故障变量。3.2故障重构故障重构的基本思路是当工业过程中发生故障时，根据故障样本数据中的测量变量推断其对应的正常测量值，以此来确定故障的幅度，找出故障变量。假如样本x是故障样本，那么x中的一部分变量值偏离了原来的正常状态。利用故障15 温州大学硕士学位论文重构可以估计出其相应的正常样本x*：xxΞe(3-1)其中，Ξ表示故障方向矢量，Ξ中的元素是0或1，0代表相应的变量是正常变量，1表示对应的变量是故障变量；e是故障幅值；是两个向量的对应元素相乘。故障重构的目标是使得x*尽可能地接近初始正常状态，即x*的统计量T2尽可能小，得到：22minminTTx*xΞe(3-2)Ξ,,eΞe当故障方向Ξ已知时，式（3-2）求解相对比较容易。但在实际工业过程中，难以获取故障方向。当故障方向未知时，式（3-2）中有两个待求变量，待求变量Ξ是离散值0或者1，待求变量e是连续变量。此时，式（3-2）是一个组合优化问题，求解非常困难。3.3Lasso重构当故障方向未知时，为了同时求解故障方向和幅值，令f=Ξe，式（3-2）变为：2minTxf(3-3)ff中的非零元素对应的变量即为故障变量，非零元素值即为故障幅值，因此求解出f即可同时得到故障方向和幅值，但是不能直接根据（3-3）求解。很显*然，式（3-3）的最优解是fx，即认为所有变量都发生了故障，这显然不符合实际情况。在故障重构中，有一个隐含的假设，即如果存在多种变量组合都能够重构回到正常状态，那么被重构故障变量的个数越少越好。故障变量个数可以表示为f，从而式（3-3）可以转化为：02minTxff(3-4)minf0可以发现式（3-4）是一个多目标优化问题，为了便于求解将其转化为单目标优化问题，用一个较小的f值近似代替minf，可以得到：002minTxff(3-5)s.t.fc0其中，fΞe，是0范数，c是指定常数，用以控制故障变量的个数。因00范数是高度非线性且非光滑的，所以求解式（3-5）的仍然非常困难。为了提高求16 第三章基于变量选择的故障重构方法解的效率，用1范数近似替代0范数有：2min(T)xff(3-6)st..df1其中表示1范数，d是一个与c有关的常数。式（3-6）可以利用拉格朗日乘1法表示成：2min()Txff(3-7)f1式中是惩罚系数。将T2的计算公式代入式(3-7)有：T1Tmin()xfP()ΛxfPf(3-8)f1对Λ进行Cholesky[53]分解：TΛLL(3-9)其中L表示一个下三角矩阵。所以，将式(3-9)代入到式(3-8)中有：TTTTTmin()PLx()()()PLfPLxPLff1(3-10)fTT假设令yPLx()、ZPL()以及βf。则，式(3-10)可以表示成：TminyZβyZββ1(3-11)βNMN1其中，ZR表示预测矩阵，yR是由响应测量的量组成的向量矩阵，β是回归系数矩阵，M是测量变量的数目，N是测量样本的总数。可以得出（3-11）是一个典型的Lasso回归问题，所以重构目标函数式（3-2）最终能够用求解Lasso问题的算法进行求解。3.3.1求解故障方向Ξ通过求解形如式（3-11）的Lasso回归问题，可以取得故障方向Ξ。因LARS的计算效率最高，且它的计算代价和最小二乘法的相接近。所以，为减少计算的复杂度，选用最小角回归（LARS）[44]求解式（3-11）。具体步骤如下：步骤1：初始化，f0，初始激活集Γ；011T步骤2：计算相关系数cPΛPxf()和最大的相关系数值Ccmaxi，ji并更新激活集ΓjjΓ1ic:iC；步骤3：计算等角矢量ujSwjj；其中，Szj()siiiΓj，sisignci，iΓjTTT1/2T1zi是zPL()中的第i列，wj(1ss1jjjj)(ssjj)1j，其中的1j是一个全1的矩阵，它的长度和Γ的大小相同；jCciiCcc步骤4：计算步长：min,，Γ是jcTT1/2TT1/2jiΓj(1ss1)aa(1ss1)jjjjijjjji17 温州大学硕士学位论文ccT+j的补集，ΓΓj{1,,}M，ΓΓj，aZui，“min”表示取最小的正jj值；TT步骤5：更新fj：()PLf()PLfjjjju111，如果β满足约束，转到步骤3。那么，通过LARS求解得到f中的非零元素对应的变量即为故障变量，从而可以求得故障方向Ξ。3.3.2参数μ的选择ff1.f2μkμ2μ1μ0μf3图3-1μ与f变化的示意图Figure3-1conceptualillustrationofthechange:μandf从式（3-10）中可以得出，μ值的选择对计算结果十分重要。μ越小，f中的非零元素，即故障变量的数目越多，而根据前面的讨论，被重构的故障变量越少越好。所以在满足重构要求的前提下，μ的值应尽可能大，f随μ变化的示意图如图3-1所示。因此，本文在进行故障重构时，选取一个较大μ值对应的故障方2向开始重构，求解式（3-11），如果重构后的T()x*高于控制限，则不断减小μ值；2并重新求解式（3-11），直到重构后的T()x*低于控制限。3.3.3计算故障幅值e如果已检测出某一故障工况的出现，为了能准确地诊断其原因，将当前工况下的故障数据按照已经取得的故障方向Ξ进行故障重构。故障重构的目标就是利用PCA和故障方向Ξ，获得重构向量x*对应正常运行时样本的最优估计，即监测指标最小。用D表示PCA的监测指标通用形式，则故障检测的指标可以表示成：T2indexxxDxxxΞe(3-12)DD1T其中DP=ΛP，利用index()x对e求导，并使其一阶导数为0。可以得出使最的小的值：18 第三章基于变量选择的故障重构方法dindexxT2xΞeDe(3-13)de令式（3-13）等于零，可以解出：1TTe=DΞDΞΞx(3-14)3.3.4故障诊断将Ξ和e式分别代入到统计量的计算公式，可以得到指数的重构值：2TTxxΞeDxΞe(3-15)如果利用Lasso得出故障方向是实际故障方向，那么重构后的指标应该小于控制限，整个故障诊断的流程如图3-2所示，最后结合变量轨迹分析找出最可能的故障变量。j=0f0=01TcPΛPxf()j更新激活集uSwjjjCciiCcmin,jcTT1/2TT1/2iΓj(1ss1)aa(1ss1)jjjjijjjjij=j+1TT()PLf()PLfujjjj111fjjej2计算重构后的TjT21TT()xxjjΞjePjΛPxjΞjej22TT(x)(控制限)jCLNoYes结束图3-2基于Lasso重构的故障诊断流程图Figure3-2FlowdiagramoffaultdiagnosisbasedonLassoreconstruction19 温州大学硕士学位论文3.4实例研究为了说明基于变量选择故障重构方法故障诊断的有效性，本节利用数值案例和TEBenchmark实例进行验证。通过交叉验证得到数值案例中的主元个数设2，此时训练集数据并未引入偏移量（正常数据）。故障诊断时，在测试集数据第二个过程变量中引入偏移量（故障数据）。TEBenchmark实例研究时，一方面，工业过程中常用4、5、7这三个预定义故障的数据来论证故障诊断方法有效性；另一方面，通过故障诊断方法得出4、5、7这三个预定义故障中的故障变量个数都为一个。因此，为了得到一个比较清晰的故障诊断结果，本文在进行TEBenchmark的实例研究时，选取4、5、7这三个常用预定义故障中的数据进行方法的有效性例证，重构的误差限定为10%。3.4.1数值案例假如一个多变量系统的数值模型是：xGte(3-15)式中，x是Mn过程数据向量，M是过程变量个数，n是采样次数，t是两个独立过程状态，G是过程状态传递矩阵，e是测量噪声。在PCA建模中，百分比控制限为99%，训练集是正常运行状态下，n=500采样得到的数据，通过交叉验证方法得出主元个数为2个。为了模拟故障，在第二个过程变量中引入偏移量，幅度为0.2，在之后的故障检测中，主元个数当成已知量。过程变量个数M会影响故障诊断的结果，所以在M分别为5、20、30、100和200时，将BAB和Lasso重构的诊断结果进行比较（如表3-1），每一个案例中的检测样本包含100个故障采样值。为了能够直观的比较计算时间，将表3-1的计算时间表示成图3-3。表3-1BAB与Lasso重构的故障诊断结果Table3-1FaultdiagnosisresultsusingBABandLassoreconstruction过程变量准确度（%）计算时间（s）个数M分支定界Lasso重构分支定界Lasso重构51001001.852.412010010013.945.145010010073.4411.1410010099353.6060.5920099981657.68141.05根据表3-1和图3-3可以得出，当过程中的变量个数M为5时，BAB和Lasso重构的所用时间基本相同；但是随着过程变量个数不断增加，同等条件下20 第三章基于变量选择的故障重构方法的BAB故障诊断所用时间程指数增长趋势，如：M为20时，BAB诊断所用时间大约是Lasso重构所用时间的3倍；M为200时，BAB诊断所用时间大约是Lasso重构所用时间的12倍。从而可以得出结论，Lasso重构的故障诊断所需计算时间明显少于BAB方法所需时间，即Lasso重构的诊断效率明显高于BAB；并且，即使在故障变量较多的条件下，Lasso重构同样可以进行有效地故障诊断。20001800BABLassoreconsitution1600140012001000800600computationtime(s)4002000050100150200variablenumber图3-3数值案例中两种方法进行故障诊断的计算时间Figure3-3Computationtimeoftwomethodsforfaultdiagnosisinthenumericalexample3.4.2TEBenchmark实例3-4TE过程流程图Figure3-4FlowchartofTEprocess21 温州大学硕士学位论文表3-2TE测量变量Table3-2ThemeasuredvariablesofTennesseeEastman变量描述变量描述x1A进料（流1）x27成分E（流6）x2D进料（流2）x28成分F（流6）x3E进料（流3）x29成分A（流9）x4总进料（流4）x30成分B（流9）x5再循环流量（流8）x31成分C（流9）x6反应器的进料速度（流6）x32成分D（流9）x7反应器压力x33成分E（流9）x8反应器液位x34成分F（流9）x9反应器温度x35成分G（流9）x10排放速度（流9）x36成分H（流9）x11产品分离器的温度x37成分D（流11）x12产品分离器的液位x38成分E（流11）x13产品分离器压力x39成分F（流11）x14分离器的出口流量（流10）x40成分G（流11）x15汽提塔的液位x41成分H（流11）x16汽提塔的压力x42D物料的流量（流2）x17汽提塔的出口流量（流11）x43E物料的流量（流3）x18汽提塔的温度x44A物料的流量（流1）x19汽提塔蒸的汽流量x45总进料量（流4）x20压缩机的功率x46压缩机的再循环阀x21反应器的冷却水出口温度x47排放阀（流9）x22分离器的冷却水出口温度x48分离器的罐液流量（流10）x23成分A（流6）x49汽提塔液体产品流量（流11）x24成分B（流6）x50汽提塔的水流阀x25成分C（流6）x51反应器的冷水流量x26成分D（流6）x52冷凝器的冷却水流量田纳西伊斯曼(TennesseeEastman,TE)[54]是一个标准的测试（Benchmarking,Process）过程，它的原型是一个实际的工业过程，常用于评价过程控制策略和多元统计过程中的监控(MultivariateStatisticalProcessMonitoring,MSPM)方法的优劣，具体流程如图3-4。TE过程共有五个主要单元:反应器、冷凝器、压缩机、分离器和汽提塔，包括A、B、C、D、E、F、G、H八种成分。气体反应物A、C、D、E与惰性气体B在被送入反应器后，产生系统产物G和H，F是产生的副产物。反应器中的主要反应包括：22 第三章基于变量选择的故障重构方法Ag()Cg()Dg()Ggliq()()Ag()Cg()Eg()()(Hgliq)Ag()Eg()Fgliq()()3()DgFgliq2()()其中，g代表气体，liq代表液体，该反应过程是放热且不可逆的。反应器中生成的产物和剩余物以气体或液体的形式进入冷凝器，随后产物流经过冷凝器降温后被送入分离器。其中的大部分的产物会被分离器分离出来，而未分离的部分产物经过压缩机的增压后被重新送入反应器中。同时，为了防止B与F的积累，排放再循环流中的一部分。经冷凝器冷凝后形成的流10进入汽提塔，而流10中的剩余反应物被去除掉。汽提塔中生成的产物G与H经过下游过程的进一步提纯就可以转化为最终产物。惰性气体B和副产物F在分离器中以气体的形式排出[55]。表3-3TE过程的故障描述Table3-3FaultdescriptionofTEprocess序号描述故障类型1A/C进料的比率波动，B的进料稳定（流4）阶跃2B的进料波动，A/C的进料比率稳定（流4）阶跃3D进料温度扰动（流2）阶跃4反应器的冷却水入口温度扰动阶跃5冷凝器的冷却水入口温度扰动阶跃6A进料损失（流1）阶跃7C处存在压力的损失—可用性的降低（流4）阶跃8A，B，C进料成分波动随机变量9D进料温度波动随机变量10C进料温度扰动随机变量11反应器的冷却水入口温度扰动进料温度扰动随机变量12冷凝器的冷却水入口温度扰动随机变量13反应器的动力性能变化慢偏移14反应器的冷却水阀门粘滞15冷凝器的冷却水阀门粘滞16-20未知21流4阀门固定在一个稳态位置恒定位置表3-2给出的是TE过程中的52个测量变量，其中22个是连续的测量变量、19个是成分测量变量、11个是控制变量。TE过程中的仿真数据可以从http://brahms.scs.uiuc.ed得到，过程数据由训练集和测试集组成，分别是500组和960组数据。所有的测试集数据每三分钟进行一次采样，总的仿真时间是48h，故障发生在8小时之后，即从第160组数据之后引入故障。如表3-3所示，仿真23 温州大学硕士学位论文过程包括了21个预定义的故障，其中的16个是已知故障，其余的5个故障则是未知故障。实例研究时，训练集(正常的过程数据)被用于建立监控模型，对应的测试集数据用来监控过程。由于稳态故障特征有助于诊断发生的根本原因[56]，因此在本文的实例研究中，用稳态时的测试集数据进行所提方法诊断效果的验证。最初的TE过程是开环、不稳定的。因而，研究者们为了更好地模拟真实工业过程，对TE提出一些不同的控制算法。文献[57]提出的分散控制算法常被应用于工业过程中的故障诊断。而且，大规模的连续过程控制通常采用分散控制系统，如石化过程。所以本文在进行故障诊断时，采用的是文献[57]提出的分散控制算法。510152025variable30354045500.060.10.30.50.8图3-5故障4中的变量故障方向Figure3-5FaultdirectionsofFault450BeforereconstructionControllimits40Afterreconstruction302T20100360460560660760860960sampleintervals图3-6故障4T2的控制图Figure3-6T2ControlchartforFault4第一个例子是故障4，通过LARS得出的故障方向矢量如图3-5所示，图中的横纵轴分别表示值变化和变量，白色代表对应变量是故障变量，黑色表示对应变量为正常变量。从图3-5中可以得出，随着逐渐增大，被诊断为故障的变量24 第三章基于变量选择的故障重构方法越来越少。比如0.1时，除了x1，x5，x7，x8，x15，x16，x18，x29，x31，x34，x36，x39，x46，x49，x52之外，其余变量都被识别为故障变量；1.0时，只有x51被识别为故障变量。因此，为使被诊断为故障的变量最少，取1.0所对应的故障方向进行故障重构。重构时选取变量轨迹变化稳定的样本数据，即第360至960时的数据，重构前故障4的T2有75%的样本点数据分布在控制限内，重构后故障4的T2有96%的样本点数据分布在控制限内，具体如图3-6所示，因此变量x51是故障4中的故障变量。4123102159051xx-10-2-3-4-501002003004005006007008009000100200300400500600700800900sampleintervalssampleintervals(a)(b)图3-7故障4变量(a):x9,(b):x51的变化轨迹Figure3-7VariabletrajectoriesforFault4:(a):x9,(b):x51故障4是由于反应器冷却水入口温度（x9）变化导致的。由于系统控制器的存在，当x9变化时，控制器通过调节反应器冷却水流量（x51），使x9回到正常值，x51出现偏离，系统恢复到初始正常状态。图3-7给出了x9和x51的变化轨迹，由图3-7可以明显地看出，x51从第161个采样时刻开始偏离初始状态，之后达到一个新的稳态，而x9在整个采样时刻内一直稳定在初始稳态左右。因此可以得出，故障4中的故障变量是x51，这与Lasso重构得出的故障诊断结果一致，进一步说明Lasso重构取得的诊断结果具有较高的准确度。510152025variable30354045500.20.40.60.81.01.2图3-8故障5中的变量故障方向Figure3-8FaultdirectionsofFault525 温州大学硕士学位论文第二个例子是故障5，通过LARS取得的故障方向矢量如图3-8所示，从图中可以得出，随着逐渐增大，被诊断为故障的变量越来越少。比如0.4时，除了x1，x10，x16，x26，x27，x29，x39，x44，x45，x50，x51之外，其余变量都被识别为故障变量；1.4时，只有x52被识别为故障变量。因此，为了使得被诊断为故障的变量最少，取所1.4对应的故障方向进行故障重构。重构时，同样选取第360至960时的采样数据。在重构前故障5的T2有88.5%的样本点数据分布在控制限内，重构后故障5的T2有91%的样本点数据分布在控制限内，具体如图3-9所示，因此变量x52是故障5中的故障变量。40Beforereconstruction35ControllimitsAfterreconstruction30252T20151050360460560660760860960sampleintervals图3-9故障5的T2控制图Figure3-9T2ControlchartforFault5故障5是由于冷凝器冷却水入口的温度变化引起的。在发生故障时，下游的过程变量也会受故障因素的影响。因冷凝器直接与分离器相连，分离器的温度（x11）首先受到影响。然后，故障影响继续传播，其他过程单元的变量也会受到影响，如汽提塔的出口流量（x17）、温度（x18）、水流阀（x50）和分离器冷却水出口温度（x22）等发生改变。为了减少故障对系统的影响，系统反馈控制器通过调节冷凝器使冷却水流量（x52）使过程系统进入一个新的稳态。在工业过程达到一个新的稳定状态后，x52并未恢复到初始的状态，而是稳定在一个新的稳态，其他的变量则再次回到初始的状态。图3-10给出了x11，x17，x18，x22，x50，和x52的轨迹。从图3-10可以明显地看出，x52从第161个采样时刻开始偏离初始状态，之后达到一个新的稳态。而其他的变量在整个采样时刻内，出现短暂偏离后，一直稳定在初始稳态左右。因此可以得出，故障5中的故障变量是x52，这与Lasso重构得出的故障诊断结果一致，同样验证了Lasso重构取得的诊断结果有较高的准确度。26 第三章基于变量选择的故障重构方法1515101055110170xx-5-5-10-10-15-1501002003004005006007008009000100200300400500600700800900sampleintervalssampleintervals(a)(b)10151055018200xx-5-5-10-10-15-1501002003004005006007008009000100200300400500600700800900sampleintervalssampleintervals(c)(d)1561045250052xx0-5-2-10-15-401002003004005006007008009000100200300400500600700800900sampleintervalssampleintervals(e)(f)图3-10故障5变量(a):x11,(b):x17,(c):x18,(d)x22,(e):x50和(f):x51的变化轨迹Figure3-10VariabletrajectoriesforFault5:(a):x11,(b):x17,(c):x18,(d):x22,(e):x50,(f):x52第三个例子是故障7，通过LARS得出的故障方向矢量如图3-11所示，随着逐渐增大，被诊断为故障的变量越来越少。比如0.1时，除了x8，x10，x14，x18，x19，x33，x34，x39，x44，x47，x52之外，其余变量都被识别为故障变量；1.1时，只有x45被识别为故障变量。所以，为了使得到故障的变量最少，取所1.127 温州大学硕士学位论文对应的故障方向进行故障重构。重构时选取600至960时刻的稳态样本数据，重构前故障7的T2只有5%的样本点数据在控制限。重构后，故障7的T2控制限内的有95%的样本点数据，具体如图3-12所示，因此变量x45是故障7中的故障变量。510152025variable30354045500.080.10.30.50.70.9图3-11故障7中的变量故障方向Figure3-11FaultdirectionsofFault790Beforereconstruction80ControllimitsAfterreconstruction7060502T403020100600650700750800850900950sampleintervals图3-12故障7的T2控制图Figure3-12T2ControlchartforFault7故障7与物料C压力损失流4（x4）有关，即流4的可用性降低。当故障发生时流4的供料流急剧下降。反应器进料速率的将会因此而降低，反应器的水平将低于正常水平，过程中的许多变量x7（反应器压强）和x29（流9）等都会因此而受到影响。在此过程中，系统控制器为了使过程恢复正常稳态，将调节流4的进料流量（x45）以弥补物料C压力损失产生的影响，使反应器恢复正常的作用水平。因此，经过一段时间的系统调节后，除x45之外的变量，在之后的采样时28 第三章基于变量选择的故障重构方法刻里恢复正常。图3-13给出了x4、x7、x29和x51的变化轨迹，由图3-13中可以明显地看出，x45从第161采样时刻开始偏离初始状态，之后达到一个新的稳态，而x4、x7以及x29在第161采样发生偏离后，在系统的调节作用下回到初始稳态。由此可以得出，故障7中的故障变量是x45，这与Lasso重构得出的故障诊断结果一致，所以验证了Lasso重构诊断结果的准确性。153025102015510475x0x0-5-5-10-10-15-20-15-2501002003004005006007008009000100200300400500600700800900sampleintervalssampleintervals(a)(b)122081510329455xx-20-7-5-12-1001002003004005006007008009000100200300400500600700800900sampleintervalssampleintervals(c)(d)图3-13故障7变量(a):x4,(b):x7,(c):x29,(d):x45的变化轨迹Figure3-13VariabletrajectoriesforFault7:(a):x4,(b):x7,(c):x29,(d):x453.5本章小结本章首先介绍了故障重构的基本原理，指出故障重构方法存在的问题。随后，为解决实际工业过程中故障方向未知条件的故障诊断问题，以Lasso变量选择技术基础，提出一种基于变量选择的故障重构方法。该方法将难以解决的组合优化问题转换成易于求解的最小二乘回归问题。在此基础上增加一个约束项，将其等价转化为Lasso回归问题，并通过最小角回归算法（LARS）求解Lasso回归问题，获取故障方向，然后利用取得故障方向进行故障重构，快速地找出故障变量。最后分别用数值案例与TEBenchmark实例研究验证了基于变量选择故障重构方法故障诊断结果的有效性与准确性。29 温州大学硕士学位论文30 第四章基于参数回归的工业过程故障诊断方法第四章基于参数回归的工业过程故障诊断方法4.1引言由于噪声、外界扰动等因素，工业过程的实际运行数据通常具有随机性和不确定性，故障诊断的结果也会具有不确定性。本章引入贝叶斯理论，用概率来表达这种不确定性。RobertTbshirani[42]指出，当回归系数具有独立、相同的拉普拉斯（Laplace）先验时，Lasso估计可以等价转化为一个Bayesian的后验估计（BayesianLasso）。同时，FusedLasso[58]可以通过惩罚系数的调节作用，自适应地调整不同变量的回归系数值，减少相互之间的影响。因此，在故障诊断建模时，将变量选择问题转化为FusedLasso回归问题，可以有效地解决BayesianLasso故障诊断的误判问题。而且，正态指数分布（normal-exponential-gamma，NEG）[59]在0处具有较好的尖峰特性。因而在贝叶斯分析时，将不同变量回归系数差jj-1的先验分布设定为NEG，可以更好地将正常变量的回归系数的值涵盖在控制限内，降低正常变量诊断为故障变量的机率，提高故障诊断的鲁棒性。随后的ParkT.和CasellaG.[60]利用Gibbs采样有效地解决了BayesianLasso估计问题。在此基础上，本章提出一种基于参数回归的故障诊断方法，将故障诊断问题转换成为变量选择中的Lasso回归问题。一方面，将Lasso回归问题转化为BayesianLasso回归问题；另一方面，将故障诊断问题转化为一个BayesianFusedLasso回归问题。在贝叶斯分析时，变量回归系数值连续变化差j-j-1的先验分布采用NEG分布、回归系数β的先验分布采用Laplace分布。随后，利用吉布斯采样方法实现Bayesian回归分析，得出回归系数的概率密度估计。在得到的回归系数概率密度的基础上，结合蒙特卡洛方法计算出变量发生故障的概率，最后通过比较找出最有可能的故障变量。实验研究表明，BayesianFusedLasso得到故障诊断结果的准确性明显好于BayesianLasso方法得出的故障诊断结果。4.2参数回归故障诊断的目的是找出造成工业过程中故障发生的故障变量，即在包含故障和正常的数据中找出故障的变量。因此，这等同于在一个包含两类不同数据的判别分析问题，其中一类是正常数据，另一类是故障数据。所以，多变量的故障诊断问题可以转化成判别分析中的变量选择问题。判别分析的目标是找到一个投影31 温州大学硕士学位论文方向，在这个投影方向上能够最大限度地减少类均值，并最大限度地减小类间方差，通常将其转化为一个优化问题的求解来实现。假设ΞΞΞ{,xx,}是121n2122一个d维样本，包含两个不同的类Ξxx,与Ξxx,。ΞΞΞ11n111n212中的n，n1，n2分别代表Ξ,,ΞΞ中观测值的个数，且n=n1+n2。对样本集数据进12行相应的标准化处理后，为了更好地将最小二乘回归问题和判别分析联系起来，定义一个预测矩阵：111xx(4-1)1x-22其中，1i是一个全1列向量，xi是Ξ中一个n1d向量；响应向量：in11n1y(4-2)n12n2其中，常数n/ni用来补偿不平衡样本量的影响。因而对于形如xβy这样的最小二乘回归问题，回归系数向量β可以通过最小化残差平方和（RSS）来获得。目标函数为：Tminyxβyxβ(4-3)βw0其中β，w是预测矩阵x中第一列对应的系数，ω则表示其余列对应的系0ω数。所以，通过求解形如的参数回归问题就能得到回归系数。其中，故障变量对应的值不为0，而正常变量对应的β值则为0，因此可以实现预期的故障诊断目标。但求得的值中有许多不为零的项（故障变量），而实际故障变量的数目是有限个的。因此在进行故障诊断时，将式（4-3）增加一项零范数约束项，以减少得到的故障变量个数。4.3BayesianLasso故障诊断通过增加约束项的方式，式（4-3）转换成一个带约束的优化问题：Tminyxβyxββ0(4-4)β其中β是关于β的零范数，是惩罚系数，调节回归系数β的稀疏度。因为零0范数L0的非凸、非光滑的，在优化问题中不易求解，因此用L1近似代替L0有：32 第四章基于参数回归的工业过程故障诊断方法Tminyxβyxββ1(4-5)β那么通过求解式（4-5）这样一个Lasso回归问题，可以取得连续变化的β值。其中的一些β会变成或接近于零，那些β等于0对应的变量判别为正常变量，而剩余β不为0所对应的变量便是要找的故障变量，从而快速地实现故障诊断。但若只将β为0的系数对应变量判别为正常变量，而将那些非常接近于0的β（比如一些β）系数对应的变量判别为故障变量，容易得出错误的诊断结果，同时又不能给出各变量发生故障时的概率。因实际工业数据通常具有随机性和不确定性，所以将贝叶斯理论引入到故障诊断，以得到更加符合实际的故障诊断结果。在数理统计中，贝叶斯定理[61]常用来描述条件概率和之间的关系，其具体形式为：()AB()A(|)BAB()AB(|)(4-6)式(4-6)也可以用除法表示成：()BAB(|)(|)BA()A(4-7)()ABA(|)(|)AB()B当样本数目接近总数时，事件A或者B在样本中发生的概率接近在A或B在整个实验中发生的概率。如果A与B是非独立事件，那么乘法公式(4-6)可以计算出事件A与B同时发生的概率。同样地，如果已知事件A与B同时发生的概率以及事件A（或B）发生的概率，那么在事件A（或B）出现概率的条件下，可以根据贝叶斯定理计算另一事件发生的概率。对于连续的随机变量，概率通常用函数的概率密度形式来表示。假如一个参数是，样本数据是y的统计模型，和y满足：y~(|)~()f(4-8)那么根据贝叶斯定理可以有：()(|)yy()(|)(|)y()y()(|)yd(4-9)式(4-9)中的分母是一个常数，所以(4-9)表示成：(|)yy()(|)(4-10)即后验概率密度正比于先验概率密度和似然函数的乘积。这样就能在给定观测样本数据的条件下，从样本数据和参数的先验概率密度中推出参数的后验概率密度估计形式。将参数的先验信息与样本信息相结合，通过导出的后验概率密度来估计参数，大致的过程如图4-1所示。33 温州大学硕士学位论文观测样本信息先验信息贝叶斯定理后验信息统计推理图4-1贝叶斯估计过程Figure4-1DiagramofBayesianestimation因此，假定β的先验[60]是：d2β|expj(4-11)j1222其中，22222,2)，2是形状参数，的先验为π()=1/或者π()=IG(0002是尺度参数。逆伽马分布IG(v,)的具体形式为：0v(v1)IG|,xvexpx(4-12)vx那么在给定样本数据时，根据贝叶斯的估计公式，参数的联合后验概率密度函数能够表示成：2221(,β|)(|,yy)|ββ(4-13)()y取对数有：222(1)11dTβ,|()yexpyx2βyxββ1(4-14)()y2后验概率密度可以表示为正比于先验概率密度与似然函数的乘积形式：222(d1)1Tβ,|expyyx2βyxββ1(4-15)2那么极大化后验密度则等于极小化1Tminyxβyxββ(4-16)212所以，可以将Lasso估计等价转化为一个后验的最值问题（BayesianLasso）。然而，后验概率密度函数很复杂，很难从后验分布直接抽样。Gibbs抽样能够在联合分布密度函数比较复杂，但每个参数的条件后验分布较为简单时进行抽样，基本思路是根据每个参数的满条件后验密度函数进行抽样。尽管Gibbs抽样不是实际的随机抽样，但当抽样的次数足够多时，其可以快速收敛于后验联合概率密度函数，因此需要推导出各个参数的满条件后验密度函数。4.3.1参数后验密度函数的推导的先验式(4-11)可以利用AndrewsD.F.和MallowsC.L.混合尺度准则[62]表示34 第四章基于参数回归的工业过程故障诊断方法成：2221β2expexpexpβd(4-17)222222222022在式(4-11)的基础上，Park和Casella[60]假定以下先验形式：d222221jβ|,,,,exp12d22(4-18)j12222jjd22222212,,,expdj(4-19)j122当σ2先验服从逆伽马分布π(2)=IG(2,2)时：00v022202(v21)002exp2(4-20)v022222参数β,,,,的满条件的后验概率密度函数可以表示成：12d22221T21β|,,yxxy,,,,~N(A,A),12dd(4-21)T1222A=xx,,,,DDdiag12d222211|,,,yxβ,,,~IG(,),12d22(4-22)2T1,,ndDyxβββ1010212|,,~IGauss(,),2jj(4-23)222,,1,2,,jd2j其中IGauss(,)的表示逆高斯分布：232xxxexp(0)(4-24)222()x回归系数的后验概率直方图可以通过对获得的后验概率密度函数进行Gibbs来获得。4.3.2Gibbs抽样的构建Gibbs抽样具体的步骤包括：0(0)(0)22（1）根据已知的参数先验分布，确定参数的初始值β,,；（2）在第k+1kk221T2k1(k1)从满条件分布β|,,yx,~NAxy,()A中抽取β，d35 温州大学硕士学位论文2(1)2kk112(1)k从满条件分布|,,yx,()~IGβ,中抽取，22（k+1）11(1)kk2(1)从满条件分布|,(β),~IGauss(,)中抽取；22jj（3）重复执行步骤（2），直到抽样次数满足设定抽样的总次数为止。4.3.3回归系数β控制限的计算根据得到抽样结果，可以画出各参数的频率直方图。为了能够更好计算变量发生故障的概率，需要获取回归系数β的后验概率密度分布。因为从直方图中直接求解概率密度分布较为困难，所以通过核密度估计方法来计算β的后验概率密度。核函数K()设定为高斯核函数（Gaussian）：2()x12Kx(,)e2(4-25)21/5其中，hcn，c1.5。为满足故障诊断的误差要求（误差低于10%），根据拉依达准则概率为0.95（95%）的要求，控制限设置为（-2，+2）。因在Lasso选择得到的回归系数中，正常变量的回归系数非常接近于0，那么β的值越接近于0，对应变量诊断为正常变量的概率就越大。因此高斯核函数的均值设定为0，这样可以更好地统计变量故障概率。同时，由于不同变量在未发生故障时，对应回归系数值的大小不同，因此为了更好地将所有正常回归系数的值包含在控制限中，将设定成变量正常时对应的回归系数方差最大值()β。从max而可以得到控制限是（02()β，02()β）。那么，在抽样的总数足够大，maxmax且分组足够细的条件下，选取频率最高一组的中间值作为参数的后验估计值，这样可以得到的β后验概率密度分布p()β。4.3.4变量故障概率的计算因吉布斯采样的总次数远大于数据集中的变量个数，所以在计算变量发生故障的故障概率时，可以用频率去估计变量发生故障的概率。假定在吉布斯采样过程中，第j个变量对应回归系数j的值有nj次落在控制限中，吉布斯采样的总次数为N。那么第j个变量发故障的概率为P：jnjP1*100%(4-26)jN比较P与95%的大小，若P95%，则第j个变量是故障变量。利用同样的方式，jj可以快速地判别是否还存在其他的故障变量，整个故障诊断的过程如图4-2。36 第四章基于参数回归的工业过程故障诊断方法开始从工业过程中采集要关注的重要变量的数据对采集的变量数据进行标准化处理故障诊断模型构建:利用标准化后的变量数据进行诊断模型的建立，将故障诊断问题转化为判别分析中的变量选择问题，得到一LassoT回归模型:minyxβyxββ1β利用吉布斯采样法求解Lasso回归模型，计算出回归系数β的概率密度函数p()β通过核密度估计方法计算出β的95%控制限j=1利用蒙特卡洛方法计算出第j个变量Pjj(>)的值，152jPj>95%是否第j个变量是故障变量第j个变量是正常变量j=52j=j+1是否输出故障诊断的结果结束图4-2基于参数回归的故障诊断流程图Figure4-2Theflowchartoffaultdiagnosisbasedonparameterregression4.3.5实例研究本节在进行TEBenchmark的实例研究时，同样选取4、5、7这三个常用预定义故障中的数据进行方法的有效性例证。实验研究中，利用交叉验证方法来确定22参数、的初始值以及惩罚系数的值，吉布斯采样的总次数为N设定为1000。22通过交叉验证得到、和为0.1。4.3.5.1案例1：故障4以测试集中第945时刻的数据为例，此时利用交叉验证得到惩罚系数=16.0，37 温州大学硕士学位论文计算得出该时刻下的()1.4β，从而可以得到正常变量对应回归系数的控制max限(2.8,2.8)。然后，依次利用吉布斯采样和核密度估计，可以快速地得出故障4中各变量回归系数的后验概率密度分布。53.553442.53321.522probability(%)probability(%)1probability(%)110.5000-100100510-10010(a)(b)(c)4443.53.53.53332.52.52.52221.51.51.5probability(%)probability(%)probability(%)1110.50.50.5000-505-505-15-10-5(d)(e)(f)图4-3故障4变量(a):x1,(b):x9,(c):x20,(d):x29,(e):x36,(f):x51的后验β概率密度分布图Figure4-3ProbabilitydensitychartofposteriorβforFault4:(a):x1,(b):x9,(c):x20,(d):x29,(e):x36,(f):x51图4-3给出了故障4部分变量的后验β概率密度分布，其中(a)、(b)、(c)、(d)、(e)、(f)分别表示变量x1，x9，x20，x29，x36，x51对应的后验概率密度分布。从图4-3中可以发现x1，x20，x29，x36的后验β主要集中在控制限内；而变量x9对应的后验β主要集中在6附近，变量x51对应的β主要集中在-11附近。这说明变量x1，x20，x29，x36是正常变量的可能性较大，而变量x9与x51是故障变量的可能性较大。38 第四章基于参数回归的工业过程故障诊断方法100%90%80%70%60%50%40%faultprobability30%20%10%0051015202530354045505560variablex图4-4故障4中各变量发生故障的概率Figure4-4FaultyprobabilityofthevariablesforFault4为准确地找出故障4中的故障变量，利用式（4-26）计算故障4中各变量发生故障的概率，直观地表示成图4-4。从图4-4可以得到变量x9的故障概率为99%，x36的故障概率是26%，以及x51的故障概率是100%，其余变量的故障概率都小于10%，这说明其余变量发生故障可能性很小。从而，可以确定故障4中最有可能的故障变量是变量x9和x51。4.3.5.2案例2：故障5同样地，以测试集中第945时刻的数据为例，此时利用交叉验证得到惩罚系数=8.0，计算得出该时刻下的()1.3β，从而可以得到正常变量对应回归系max数的控制限(2.6,2.6)。然后，依次利用吉布斯采样与核密度估计，可以得到故障5中每个变量的回归系数后验概率密度分布。图4-5给出了故障5部分变量的后验β概率密度分布，其中(a)、(b)、(c)、(d)、(e)、(f)分别表示变量x11，x17，x18，x22，x50，x52对应的后验β概率密度分布。从图4-5中可以发现x11，x22的后验β主要集中在控制限内；而变量x17对应的后验β主要集中在-7附近，变量x18对应的后验β主要集中在-15附近，变量x50对应的后验β主要集中在15附近，以及变量x52对应的后验β主要集中在-200附近。这表明变量x11，x22正常的可能性较大，而变量x17，x18，x50和x52出现故障的可能性较大。39 温州大学硕士学位论文56655444333222probability(%)probability(%)probability(%)111000-505-10010-20020(a)(b)(c)45303.52543202.5321521.510probability(%)probability(%)probability(%)1150.5000-505-20020-400-2000(d)(e)(f)图4-5故障5变量(a):x11,(b):x17,(c):x18,(d):x22,(e):x50,(f):x52的后验β概率密度分布图Figure4-5ProbabilitydensitychartofposteriorβforFault5:(a):x11,(b):x17,(c):x18,(d):x22,(e):x50,(f):x52100%90%80%70%60%50%40%faultprobability30%20%10%0051015202530354045505560variablex图4-6故障5中各变量发生故障的概率Figure4-6FaultyprobabilityofthevariablesforFault540 第四章基于参数回归的工业过程故障诊断方法为准确地找出故障5中的故障变量，利用式（4-26）计算故障5中各变量发生故障的概率，具体如图4-6所示。从图4-6可以得出变量x17，x18，x50和x52发生故障的概率为100%，x19的故障概率为58%，x20的故障概率为18%，以及x22发生故障的概率为22%，其余变量的故障概率都小于10%，这说明其余变量发生故障可能性很小。从而，可以确定故障5中最有可能的故障变量是变量x17，x18，x50和x52。4.3.5.3案例3：故障7463.53.53532.542.52231.51.52probability(%)probability(%)probability(%)1110.50.5000-10010-10010-10010(a)(b)(c)453.53.53432.52.53221.521.5probability(%)probability(%)probability(%)1110.50.5000-505-15-10-5-505(d)(e)(f)图4-7故障7变量(a):x4,(b):x7,(c):x20,(d):x29,(e):x45,(f):x50的后验β概率密度分布图Figure4-7ProbabilitydensitychartofposteriorβforFault7:(a):x4,(b):x7,(c):x20,(d):x29,(e):x45,(f):x50与故障4和故障5一样，利用测试集中第945时刻的数据，此时利用交叉验证得到惩罚系数=10.0，计算得出该时刻下的()1.0β，从而可以得到正常变max量对应回归系数的控制限(2.0,2.0)。随后，依次利用吉布斯采样与核密度估41 温州大学硕士学位论文计，就能得到故障7中每个变量的回归系数后验概率密度分布。图4-7给出的是故障7部分变量的后验β概率密度分布，其中(a)、(b)、(c)、(d)、(e)、(f)分别表示变量x4，x7，x20，x29，x45，x50对应的后验β概率密度分布。从图4-7中可以发现x4，x7和x29的后验β主要集中在控制限内；而变量x20对应的后验β主要集中在3附近，变量x45对应的后验β主要集中在-10附近，以及变量x50对应的后验β主要集中在4附近。这说明变量x4，x7与x29是正常变量的可能性较大，而变量x20，x45，与x50是故障变量的可能性较大。100%90%80%70%60%50%40%faultprobability30%20%10%0051015202530354045505560variablex图4-8故障7中各变量发生故障的概率Figure4-8FaultyprobabilityofthevariablesforFault7为准确地找出故障7中的故障变量，利用式（4-26）计算故障7中各变量发生故障的概率，具体如图4-8所示。从图4-8可以得出变量x4发生故障的概率为68%，x20发生故障的概率是57%，x45发生故障的概率为100%，其余变量的故障概率都小于50%，这说明其余变量发生故障可能性较小。从而，可以得出故障7中最有可能的故障变量是变量x4，x20和x44.4BayesianFusedLasso故障诊断对于由式(4-6）的预测矩阵x和式(4-7）的响应矩阵y确定的系统模型，FusedLasso[58]回归可以表示成：dTminyxβyxβ1β11|jj1|(4-27)βj2其中，正则化系数用来调节回归系数的稀疏度，用以调节回归系数变化差12异的光滑度。当时0，FusedLasso模型将转化成为Lasso回归模型。因为142 第四章基于参数回归的工业过程故障诊断方法FusedLasso可以利用和分别调整回归系数的稀疏度与回归系数变化差异的12光滑度，所以其既能更好地调整变量的稀疏性，使回归系数的值分布变得更加集中，又能有效地减少不同变量回归系数之间的相互干扰作用。与Lasso回归模型一样，FusedLasso也可以等价转化为贝叶斯估计的形式[63]。FusedLasso回归系数的β先验分布可以表示成：21dd22212β|expβ1jj1(4-28)j2根据AndrewsD.F.和MallowsC.L.的混合尺度准则[62]，式（4-28）可以转化为：21dd222222221j11β|expexpd22jjj122222jj2(4-29)d1jj1222222expexpd22jjj222222jj2d111dd222T12221jjexp2βββjj122(4-30)dddd2222jjjjddjjjj1212其中111-00022211211111--002222222233111110-222200(4-31)β33341111000++-2222dd11dd111000-22+2ddd因此，将变量选择问题转化为FusedLasso回归问题后，同样可以将其等价转化贝叶斯分析问题，并利用吉布斯采样实现贝叶斯分析，计算出变量发生故障的概率。β的先验概率密度函数表达式为：2jjNEGj|,=kDexp221(4-32)4其中，k2Γ1/2是归一化常数，D−2λ−1是抛物柱面函数，正则43 温州大学硕士学位论文化参数和调节回归系数的稀疏程度，的值越大或的值越小，回归系数的值分布则越稀疏。抛物柱面函数可以通过二阶线性常微分方程的解求得：22dz1a0(4-33)2dz42式（4-32）的积分表达式为：2β11β22βDd21exp2exp(4-34)Γ21402NEG密度函数式（4-32）可以表示为分层表达式的形式：221j2122NEG|,jexpjjjjjexpjjexp22dd2j2jΓ(4-35)2222=N|0,EXP(jjjjjj|j)G|a(,)ddLaplace0.6NEG0.5y0.4tisn0.3ed0.20.10.0-4-3-2-101234β图4-9拉普拉斯分布和NEG分布的示意图Figure4-9TheschematicdiagramoftheLaplaceandtheNEGdistribution图4-9给出了拉普拉斯分布和NEG分布的示意图，从图中可以明显发现NEG在0处有更好的尖峰特性。正常变量的回归系数通过FusedLasso调节作用，其值的连续变化差异将变小，从而正常变量回归系数的值将更多的集中在0附近。所以，将FusedLasso回归模型中回归系数之间的差异变化βj-βj-1先验分布设定为NEG，有利于提高诊断结果的准确度。在BayesianFusedLasso诊断方法中，式（4-3）的目标函数被等价转化为：44 第四章基于参数回归的工业过程故障诊断方法dTminyxβyxβ121β1jj(4-36)βj2因在将式（4-36）FusedLasso回归问题转化为贝叶斯分析时，变量回归系数值连续变化差βj-βj-1的先验分布采用NEG分布，回归系数β的先验分布采用Laplace分布。因而可以将模型中回归系数的先验定义成：212ddd22jjj1β|Laplace|NEG|,122(4-37)22jj12同样，式（4-37）可以表示成分层的形式：2d12dd22jjj1β|Laplace|NEG|,12222jj12dd2221jj11=expexp22jj11222222jj2(4-38)d1jj1d2exp22jjexpj222j22jjj222dddd21222exp2dddjjjjj2jjjj212Γ22222222从而，参数β,,,,,,,,,,,,的满后验概率密度函数为：12d23d23d22222221T21β|,,yxxy,,,,,,,,,,,~(A,A)N12232dddd1A=xxT1(4-39)β2222222|,,,yxβ1,22,3211,dd,n,,,,,,~IG2,22101nd(4-40)TT11yxβyxββββ0122221|,,~IGauss22j,1(4-41)1jj212|,,2,~IGauss,2j(4-42)2jjjj12jjj1222j|jj,2,222~Ga1,(4-43)1其中，β由式（4-31）确定。这样，就可以根据求得的参数满后验概率密度分布进行吉布斯采样。Gibbs抽样具体的步骤为：(0)2(0)2(0)2(0)(0)（1）利用已知参数的先验分布来确定参数的初始值β,(),(),(),；45 温州大学硕士学位论文（2）在第k+1222k1kTkk21从满条件分布(|,,(βyxxy),(),(),)~N(A,A)中抽取d(k1)β，2(1)22kkkk从满条件分布(()|,,yx,(),(β),)~IG(2,2)中抽取j112(k1)()，(1)k122(()1)k1从满条件分布|,(β)(1)kk,~IGauss2(1)2,1中抽取，22(1)1k12jj()j2(kk1)()12()从满条件分布|,(β(k1))kk,k2(~IGauss1)()(),2j中2jj(kk1)(1)2j()(jj)1(1)k1抽取，2j(1)kkk22(1)2从满条件分布|,,~IGauss1,()中抽取jjj2222(k1)j；（3）重复执行步骤（2），直到抽样次数满足设定抽样的总次数为止。其他的故障诊断步骤与BayesianLasso故障诊断方法的步骤相同。4.4.1实例研究本节在进行TEBenchmark的实例研究时，利用4、5、7这三个常用预定义故障中的数据进行两种参数回归方法的交叉验证。同样地，实验研究中，利用交叉验证方法来确定BayesianFusedLasso中的参数222、、、的初始值以及惩罚系数1和2的值，吉布斯采样的总次数为N设定为1000。通过交叉验证得22216.0、到和为0.1，的初始值为0.2，的初始值为0.3；故障4中的119.0，故障5中的8.0、10.0以及故障7中的10.0、12.0。21212这其中，最常用的交叉验证方法是K折交叉验证[64]，本次实例研究中的K取10。实例验证时，首先以测试集中第945时刻数据得出的诊断结果为例，比较BayesianLasso诊断和BayesianFusedLasso方法故障诊断的准确性；然后，用稳态时测试集数据（第936至960采样时刻的数据）的诊断结果进行进一步的验证。为便于图表说明，用BLasso和BFNLasso分别代表BayesianLasso故障诊断和BayesianFusedLasso故障诊断，Lasso重构表示基于变量选择的故障重构方法。4.4.1.1案例1：故障4故障4中正常变量回归系数β的控制限(2.8,2.8)。在第945采样时刻，通过BFNLasso方法得出部分变量的后验β概率密度分布如图4-10，其中(a)、(b)、46 第四章基于参数回归的工业过程故障诊断方法(c)、(d)、(e)、(f)分别表示变量x1，x9，x20，x29，x36，x51对应的后验β概率密度分布。53.553442.53321.522probability(%)probability(%)1probability(%)110.5000-10010-10010-10010(a)(b)(c)4543.53.54332.52.532221.51.5probability(%)probability(%)probability(%)1110.50.5000-505-505-15-10-5(d)(e)(f)图4-10故障4变量(a):x1,(b):x9,(c):x20,(d):x29,(e):x36,(f):x51的后验β概率密度分布图Figure4-10ProbabilitydensitychartofposteriorβforFault4:(a):x1,(b):x9,(c):x20,(d):x29,(e):x36,(f):x51100%90%80%70%60%50%40%faultprobability30%20%10%0051015202530354045505560variablex图4-11故障4中各变量发生故障的概率Figure4-11FaultyprobabilityofthevariablesforFault447 温州大学硕士学位论文从该图中可以明显地发现x1，x9，x20，x29，x36的后验β主要分布在控制限内；只有变量x51对应的后验β主要分布在控制限外，因而可以将变量x51判定为故障4中的故障变量。然后，利用式（4-26）可以计算出故障4中各变量发生故障的概率，将其直观地表示成图4-11。从图4-11可以得到变量x51发生故障的概率为100%，x9的故障概率为2%，其余变量的故障概率都为0，这说明其余变量是正常变量。从而能够以100%的概率确定故障4中的故障变量是变量x51，这与Lasso重构得到故障4中的诊断结果一致。100%100%990%51xx80%90%70%BLasso60%80%BFNLasso50%40%BLasso70%30%BFNLasso20%60%faultprobabilityof10%faultprobabilityof050%935940945950955960935940945950955960sampleintervalssampleintervals(a)(b)图4-12故障4中变量(a):x9,(b):x51发生故障的概率Figure4-12FaultprobabilityforFault4:(a):x9,(b):x51为进一步地验证BFNLasso故障诊断的准确性，以稳态时的测试集数据（第936至960采样时刻的数据）得到的诊断结果为例，比较两种方法故障诊断的准确性。图4-12给出了故障4中，两种方法得出的变量x9和x51发生故障的概率。图4-12中（a）是x9发生故障的概率，其中BLasso得出的x9故障概率几乎超过50%，能够以86%的平均概率将x9判定为故障变量，但实际上x9是正常变量；相比之下，BLasso得出的x9故障概率均小于10%，可以直接将x9判定为正常变量。图中（b）给出两种方法得到x51发生故障的概率，可以发现两种方法均以100%的概率将x51判定为故障变量。因此可以得出结论，BLasso易将正常变量判定为故障变量，该方法的故障诊断误判率高于BFNLasso。表4-1故障4中10折交叉验证结果的误报率Table4-1Thefalsepositiverateby10-foldcross-validationinFault4误报率（%）K12345678910BLasso21.217.417.421.215.426.921.215.428.819.2BFNLasso0000000000为进一步验证该结论的正确性，利用10折交叉验证方法，计算出故障4中48 第四章基于参数回归的工业过程故障诊断方法BLasso和BFNLasso故障诊断时的误判率，具体如表4-1。从表中的数据可以发现，在每一次交叉验证中，BLasso诊断的误判都要明显地高于BFNLasso，BLasso在十次交叉中的平均误判率为20.4%；而同等条件下，BFNLasso在每次交叉验证中的误判率都为0，这足以说明BLasso故障诊断的误判率明显地高于BFNLasso故障诊断的误判率。4.4.1.2案例2：故障5555444333222probability(%)probability(%)probability(%)111000-505-10010-10010(a)(b)(c)3.5653542.543231.522probability(%)1probability(%)probability(%)10.51000-505-10010-200-1000(d)(e)(f)图4-13故障5变量(a):x11,(b):x17,(c):x18,(d):x22,(e):x50,(f):x52的后验β概率密度分布图Figure4-13ProbabilitydensitychartofposteriorβforFault5:(a):x11,(b):x17,(c):x18,(d):x22,(e):x50,(f):x52故障5正常变量回归系数的控制限为(2.6,2.6)。在第945采样时刻，通过BFNLasso方法得出部分变量的后验β概率密度分布如图4-13，其中(a)、(b)、(c)、(d)、(e)、(f)分别表示变量x11，x17，x18，x22，x50，x52对应的后验β概率密度分布。从该图中可以明显地发现x11，x17，x18，x22，x50的后验β主要集中在控制49 温州大学硕士学位论文限内；只有变量x52对应的后验β主要集中在-100附近，因而可以将变量x52判定为故障5中的故障变量。然后，利用式（4-26）计算出故障5中各变量发生故障的概率，将其直观地表示成图4-14。从图4-14可以发现，变量x52发生故障的概率为100%，x18的故障概率为20%，其余变量的故障概率都小于10%，这说明其余变量发生故障的可能性很小。因此可以以100%的概率将故障5中的变量x52判定为故障变量，这与Lasso重构得到故障5中的诊断结果一致。100%90%80%70%60%50%40%faultprobability30%20%10%0051015202530354045505560variablex图4-14故障5中各变量发生故障的概率Figure4-14FaultyprobabilityofthevariablesforFault5图4-15给出的是两种方法在936至960采样时刻得出变量x17、x18、x50和x52发生故障的概率。从图中可以发现，BLasso得出x17、x18、x50发生故障的概率为100%，但实际上x17、x18、x50这三个变量都是正常变量；BFNLasso得出x17发生故障的平均概率是9.5%、x18发生故障的平均概率是18%、x50发生故障的平均概率是1%，均明显小于BLasso得出对应变量的故障概率。两种方法将变量x52判别为故障变量的概率均是100%，因此能够以100%的概率认为变量x52发生了故障，其是故障5中最有可能的故障变量。所以，同样可以得出结论，BLasso易将正常变量判定为故障变量，其故障诊断的误判率高于BFNLasso。为进一步验证该结论，利用10折交叉验证方法，计算出故障5中BLasso和BFNLasso故障诊断时的误判率，具体如表4-2。从表中的数据可以发现，在每一次交叉验证中，BLasso诊断的误判都要明显地高于BFNLasso，BLasso在十次交叉中的平均误判率为19.8%；而同等条件下，BFNLasso在每次交叉验证中的误判率都为0，这同样说明BLasso的故障诊断误判率明显地高于BFNLasso故障诊的误判率。50 第四章基于参数回归的工业过程故障诊断方法100%100%1790%BLasso1890%BLassox80%BFNLassox80%BFNLasso70%70%60%60%50%50%40%40%30%30%20%20%10%10%faultprobabilityoffaultprobabilityof00935940945950955960935940945950955960sampleintervalssampleintervals(a)(b)100%100%5090%BLasso52x80%BFNLassox90%70%BLasso60%80%BFNLasso50%40%70%30%20%60%10%faultprobabilityoffaultprobabilityof050%935940945950955960935940945950955960sampleintervalssampleintervals(c)(d)图4-15故障5中变量(a):x17,(b):x18,(c):x50,(d):x52发生故障的概率Figure4-15FaultprobabilityforFault5:(a):x17,(b):x18,(c):x50,(d):x52表4-2故障5中10折交叉验证结果的误报率Table4-2Thefalsepositiverateby10-foldcross-validationinFault5误报率（%）K12345678910BLasso21.125.023.121.136.519.219.219.219.215.4BFNLasso00000000004.4.1.3案例3：故障7故障7中正常变量回归系数的控制限(2.0,2.0)。在第945采样时刻，通过BFNLasso方法得出部分变量的后验β概率密度分布如图4-16，其中(a)、(b)、(c)、(d)、(e)、(f)分别表示变量x4，x7，x20，x29，x45，x50对应的后验β概率密度分布。从该图中可以明显地发现x4，x7，x20，x29，x50的后验β主要分布在控制限(2.0,2.0)内；只有变量x45对应的后验β主要集中在-7.5附近，因而可以将变量x45判定为故障7中的故障变量。然后，利用式（4-26）计算出故障7中各变量发生故障的概率，将其直观地表示成图4-17。从图4-17可以发现，变量x45发生故障的概率为100%，其余变量的故障概率都小于10%，这说明其余变量发生故障的可能性很小。因此能够以100%的概率将变量x45判定为故障7中的故障51 温州大学硕士学位论文变量，这与Lasso重构得到故障7中的诊断结果一致。4473.53.563352.52.542231.51.5probability(%)probability(%)probability(%)2110.50.51000-10010-10010-10010(a)(b)(c)4443.53.53.53332.52.52.52221.51.51.5probability(%)probability(%)probability(%)1110.50.50.5000-505-15-10-5-10010(d)(e)(f)图4-16故障7变量(a):x4,(b):x7,(c):x20,(d):x29,(e):x45,(f):x50的后验β概率密度分布图Figure4-16ProbabilitydensitychartofposteriorβforFault7:(a):x4,(b):x7,(c):x20,(d):x29,(e):x45,(f):x50100%90%80%70%60%50%40%faultprobability30%20%10%0051015202530354045505560variablex图4-17故障7中各变量发生故障的概率Figure4-17FaultyprobabilityofthevariablesforFault752 第四章基于参数回归的工业过程故障诊断方法图4-18给出的是两种方法在936至960采样时刻得出变量x4和x45发生故障的概率。从图中可以发现，BLasso得出x4发生故障的概率都超过50%，可以以88%的平均概率认为故障4发了故障，但实际上x4是正常变量；BFNLasso得出变量x4的故障概率都小于50%，平均的故障概率是15%，明显地小于BLasso给出变量x4的故障概率。两种方法得出变量x45发生故障的概率均为100%，因此能够以100%的概率认为变量x45发生了故障，其是故障7中最有可能的故障变量。所以，同样可以得出结论，BLasso易将正常变量判定为故障变量，其故障诊断的误判率高于BFNLasso。100%100%490%45xx80%90%70%BLassoBLasso60%BFNLasso80%BFNLasso50%40%70%30%20%60%faultprobabilityof10%faultprobabilityof050%935940945950955960935940945950955960sampleintervalssampleintervals(a)(b)图4-18故障7中变量(a):x4,(b):x45发生故障的概率Figure4-18FaultprobabilityforFault7:(a):x4,(b):x45表4-3故障7中10折交叉验证结果的误报率Table4-3Thefalsepositiverateby10-foldcross-validationinFault7误报率（%）K12345678910BLasso46.242.348.155.853.844.244.25044.242.3BFNLasso0000000000为进一步验证该结论，利用10折交叉验证方法，计算出故障7中BLasso和BFNLasso故障诊断时的误判率，具体如表4-3。从表中的数据可以发现，在每一次交叉验证中，BLasso诊断的误判都要明显地高于BFNLasso，BLasso在十次交叉中的平均误判率为47.1%；而同等条件下，BFNLasso在每次交叉验证中的误判率都为0，这同样说明BLasso的故障诊断误判率明显地高于BFNLasso故障诊的误判率。4.5本章小结针对实际工业运行数据具有的随机性和不确定性，本章提出了一种基于参数53 温州大学硕士学位论文回归的工业过程故障诊断方法。该方法将故障诊断问题转化为一个变量选择中的Lasso回归问题。在此基础上，一方面，将Lasso回归问题转化为BayesianLasso回归问题；另一方面，将故障诊断问题转化为一个BayesianFusedLasso回归问题，在贝叶斯分析时，变量回归系数值连续变化差j-j-1的先验分布采用NEG分布和回归系数β的先验分布采用Laplace分布。进一步等价，将其转化为Bayesian回归问题，利用吉布斯采样方法实现Bayesian回归分析，得出回归系数的概率密度估计。然后结合核密度估计方法计算出正常变量回归系数的控制限，并根据蒙特卡洛方法得出变量发生故障的概率，通过比较找出最有可能发生故障的变量，得到的诊断结果与实际工况相符，最后通过TE过程的3个故障案例验证了该方法进行故障诊断的有效性。实例研究表明BayesianFusedLasso故障诊断的鲁棒性明显好于BayesianLasso故障诊断。54 第五章总结与展望第五章总结与展望5.1总结本文针对实际工业过程中故障方向未知、数据具有随机性和不确定性的问题，分别提出对应的解决方案：针对实际工业过程中故障方向未知的问题，本文提出一种基于变量选择的故障重构方法，将故障诊断问题转化为变量选择中的Lasso回归问题，并通过最小角回归算法（LARS）求解Lasso回归问题取得故障方向，然后利用求得故障方向进行故障重构，快速地找出故障变量。数值案例和TEBenchmark实例研究的结果表明基于变量选择的故障重构方法能够有效地提取故障方向、找出故障变量。数据驱动的故障诊断方法在进行故障诊断时，通常得出一个确定结论，而实际工业数据通常具有随机性和不确定性，且数据驱动的故障诊断方法不能准确地给出变量发生故障的概率。因此，为了增强故障诊断的鲁棒性，准确地计算出变量发生故障的概率，本文提出一种基于参数回归的故障诊断方法，将故障诊断问题转化为一个变量选择中Bayesian回归问题，利用吉布斯采样方法实现Bayesian回归分析，得到回归系数的概率密度估计，然后结合蒙特卡洛方法计算出变量发生故障的概率，最后通过比较找出最有可能的故障变量。5.2展望本文针对实际工业过程中故障方向未知、数据具有随机性和不确定性的问题，进行大量的研究，提出相应的解决方法。接下来有待深入研究的工作包括：1本文提出的诊断方法主要集中在解决故障方向提取和变量故障概率获取上，而对故障诊断的实时性没有进行深入的研究，因此接下来的研究将对方法进行改进和完善，提高本文方法进行故障诊断的实时性。2本文在进行故障诊断时，将故障诊断问题转化为Lasso变量选择问题，利用求解Lasso问题的方法来解决故障诊断问题。而现代工业过程因其自动化与智能化程度越来越高，产生的数据规模巨大，因而在求解Lasso问题面临计算复杂度与时间复杂度的问题。一方面计算复杂度的减少会在一定程度上增加时间复杂度，而另一方面时间复杂度降低也会增加计算复杂度，因此如何实现计算复杂度和时间复杂度的平衡是接下来的研究需要重点关注的问题。55 温州大学硕士学位论文56 参考文献参考文献[1]GeZ,SongZ,GaoF.Reviewofrecentresearchondata-basedprocessmonitoring[J].Industrial&EngineeringChemistryResearch,2013,52(10):3543-3562.[2]ChiangLH,KotanchekME,KordonAK.FaultdiagnosisbasedonFisherdiscriminantanalysisandsupportvectormachines[J].Computers&ChemicalEngineering,2004,28(8):1389-1401.[3]贾庆贤,张迎春,管宇,等.基于解析模型的非线性系统故障诊断方法综述[J].信息与控制,2012(3):356-364.[4]KesavanP,LeeJH.Diagnostictoolsformultivariablemodel-basedcontrolsystems[J].Industrial&EngineeringChemistryResearch,1997,36(7):2725-2738.[5]FrankPM.Faultdiagnosisindynamicsystemsusinganalyticalandknowledge-basedredundancy:asurveyandsomenewresults[J].Automatic,1990,26(3):459-474.[6]张妮,车立志,吴小进.基于数据驱动的故障诊断技术研究现状及展望[J].计算机科学,2017(S1):37-42.[7]GaoZ,CecatiC,DingSX.ASurveyofFaultDiagnosisandFault-TolerantTechniques—PartII:FaultDiagnosisWithKnowledge-BasedandHybrid/ActiveApproaches[J].IEEETransactionsonIndustrialElectronics,2015,62(6):3768-3774.[8]WenCL,Fei-YaLV,BaoZJ,LiuMQ.AReviewofDataDriven-basedIncipientFaultDiagnosis[J].ActaAutomaticaSinica,2016,42(9):1285-1299.[9]FrankPM,DingX.Surveyofrobustresidualgenerationandevaluationmethodsinobserver-basedfaultdetectionsystems[J].JournalofProcessControl,1997,7(6):403-424.[10]吴今培.智能故障诊断技术的发展和展望[J].振动、测试与诊断,1999,19(2):79-86,147.[11]KoiwanitJ,SupapT,ChanC,etal.Anexpertsystemformonitoringanddiagnosisofammoniaemissionsfromthepost-combustioncarbondioxidecaptureprocesssystem[J].InternationalJournalofGreenhouseGasControl,2014,26(7):158-168.[12]马纪明,万蔚,王法岩.基于模型诊断技术的神经网络实现方法[J].北京航空航天大学学报,2013(2):178-183.[13]YuM,WangD,LuoM,etal.Faultdetection,isolationandidentificationforhybridsystemswithunknownmodechangesandfaultpatterns[J].ExpertsystemswithApplications,2012,39(11):9955-9965.[14]王文辉,周东华.基于定性和半定性方法的故障检测与诊断技术[J].控制理论与应用,2002,19(5):653-659.[15]李忠海,刘颖,杨燕.基于模糊推理的故障诊断专家知识库设计[J].沈阳航空航天大学学报,2013,30(02):54-59.[16]王文晶,刘锐,张丽梅.基于故障树的故障诊断专家系统推理机的设计[J].电脑知识与技术,2014,10(02):404-406.57 温州大学硕士学位论文[17]GrbovicM,LiW,XuP,etal.DecentralizedfaultdetectionanddiagnosisviasparsePCAbaseddecompositionandMaximumEntropydecisionfusion[J].JournalofProcessControl,2012,22(4):738-750.[18]GuoY,LiG,ChenH,etal.AnenhancedPCAmethodwithSavitzky-GolaymethodforVRFsystemsensorfaultdetectionanddiagnosis[J].Energy&Buildings,2017,142(5):167-178.[19]易辉.基于支持向量机的故障诊断及应用研究[D].南京航空航天大学,2011.[20]张培林,李胜.基于小波包变换和GA-PLS算法的故障特征选择方法[J].振动、测试与诊断,2014,34(2):385-391.[21]张峰.基于费舍尔判别分析的故障诊断方法研究[D].哈尔滨理工大学,2011.[22]DuniaR,JoeQinS.Subspaceapproachtomultidimensionalfaultidentificationandreconstruction[J].AIChEJournal,1998,44(8):1813-1831.[23]YanZ,YaoY.Variableselectionmethodforfaultisolationusingleastabsoluteshrinkageandselectionoperator(LASSO)[J].ChemometricsandIntelligentLaboratorySystems,2015,146:136-146.[24]Garcia-AlvarezD,FuenteMJ,SainzGI.Faultdetectionandisolationintransientstatesusingprincipalcomponentanalysis[J].JournalofProcessControl,2012,22(3):551-563.[25]PengX,LiQ,WangK.FaultdetectionandisolationforselfpoweredneutrondetectorsbasedonPrincipalComponentAnalysis[J].AnnalsofNuclearEnergy,2015,85:213-219.[26]LiuK,FeiZ,YueB,etal.Adaptivesparseprincipalcomponentanalysisforenhancedprocessmonitoringandfaultisolation[J].ChemometricsandIntelligentLaboratorySystems,2015,146:426-436.[27]ChiangL.H.,RussellE.L.,BraatzR.D.Faultdetectionanddiagnosisinindustrialsystems[M].SpringerScience&BusinessMedia,2000.[28]RussellE.L.,ChiangL.H.,BraatzR.D.Data-drivenMethodsforFaultDetectionandDiagnosisinChemicalProcesses[M].SpringerLondon,2000.[29]MillerP.,SwansonR.E.,HecklerC.E.Contributionplots:Amissinglinkinmultivariatequalitycontrol[J].AppliedMathematicsandComputerScience,1998,8(4):775-792.[30]VandenKerkhofP.,VanlaerJ.,GinsG.,&VanImpeJ.F.Analysisofsmearing-outincontributionplotbasedfaultisolationforstatisticalprocesscontrol[J].ChemicalEngineeringScience,2013,104(50):285-293.[31]郭明,王树青.基于特征子空间的系统性能监控与工况识别[J].化工学报,2004,55(1):151-154.[32]HuangJ,YanX.Dynamicprocessfaultdetectionanddiagnosisbasedondynamicprincipalcomponentanalysis,dynamicindependentcomponentanalysisandBayesianinference[J].ChemometricsandIntelligentLaboratorySystems,2015,148:115-127.[33]LiL,YuY,YangJ,etal.Voltammetricelectronictongueforthequalitativeanalysisofmilkadulteratedwithureacombinedwithmulti-waydataanalysis[J].Int.J.Electrochem.Sci,2015,10:5970-5980.[34]XuY,LiuY,ZhuQ.MultivariatetimedelayanalysisbasedlocalKPCAfaultprognosisapproachfornonlinearprocesses[J].ChineseJournalofChemicalEngineering,2016,24(10):1413-1422.[35]YueHH,QinSJ.Reconstruction-basedfaultidentificationusingacombinedindex[J].Industrial&EngineeringChemistryResearch,2001,40(20):4403-4414.[36]AlcalaCF,QinSJ.Reconstruction-basedcontributionforprocessmonitoring[J].Automatica,2009,45(7):58 参考文献1593-1600.[37]KariwalaV,OdioweiPE,CaoY,etal.Abranchandboundmethodforisolationoffaultyvariablesthroughmissingvariableanalysis[J].JournalofProcessControl,2010,20(10):1198-1206.[38]HeB,YangX,ChenT,etal.Reconstruction-basedmultivariatecontributionanalysisforfaultisolation:Abranchandboundapproach[J].JournalofProcessControl,2012,22(7):1228-1236.[39]OhYS,MoKJ,YoonES,etal.Faultdiagnosisbasedonweightedsymptomtreeandpatternmatching[J].Industrial&EngineeringChemistryResearch,1997,36(7):2672-2678.[40]ChuYH,QinSJ,HanC.Faultdetectionandoperationmodeidentificationbasedonpatternclassificationwithvariableselection[J].Industrial&EngineeringChemistryResearch,2004,43(7):1701-1710.[41]MusulinE,YélamosI,PuigjanerL.Integrationofprincipalcomponentanalysisandfuzzylogicsystemsforcomprehensiveprocessfaultdetectionanddiagnosis[J].Industrial&EngineeringChemistryResearch,2006,45(5):1739-1750.[42]TibshiraniR.Regressionshrinkageandselectionviathelasso[J].JournaloftheRoyalStatisticalSociety.SeriesB(Methodological),1996,(58):267-288.[43]BertsekasDP.Nonlinearprogramming[M].Belmont:Athenascientific,1999.[44]EfronB,HastieT,JohnstoneI,etal.Leastangleregression[J].TheAnnalsofstatistics,2004,32(2):407-499.[45]JolliffeIT.PrincipalComponentAnalysisandFactorAnalysis[M]//Principalcomponentanalysis.SpringerNewYork,1986:115-128.[46]单鸿亮,孟濬.基于主成分分析的一类工业复杂控制系统变量降维选取[J].硅酸盐学报,2002(z1):131-134.[47]JacksonJE,MudholkarGS.Controlproceduresforresidualsassociatedwithprincipalcomponentanalysis[J].Technometrics,1979,21(3):341-349.[48]SpallJC.Estimationviamarkovchainmontecarlo[J].IEEEControlSystems,2003,23(2):34-45.[49]GemanS,GemanD.Stochasticrelaxation,Gibbsdistributions,andtheBayesianrestorationofimages[J].IEEETransactionsonpatternanalysisandmachineintelligence,1984(6):721-741.[50]钱敏平.马链的细致平衡(可逆性)与环流分解[J].科学通报,1979,24(9):389-390.[51]WuWB,MielniczukJ.Kerneldensityestimationforlinearprocesses[J].TheAnnalsofstatistics,2002,30(5):1441-1459.[52]张学工.模式识别.第3版[M].清华大学出版社,2010.[53]陈建平.Cholesky分解递归算法与改进[J].计算机研究与发展,2001,38(8):923-926.[54]DownsJJ,VogelEF.Aplant-wideindustrialprocesscontrolproblem[J].Computers&ChemicalEngineering,1993,17(3):245-255.[55]田文德,孙素莉,汪海.田纳西-伊斯曼过程故障特性仿真[J].青岛科技大学学报(自然科学版),2009,30(1):74-79.[56]YoonS,MacGregorJF.FaultdiagnosiswithmultivariatestatisticalmodelspartI:usingsteadystatefaultsignatures[J].JournalofProcessControl,2001,11(4):387-400.[57]RickerNL.DecentralizedcontroloftheTennesseeEastmanchallengeprocess[J].JournalofProcessControl,59 温州大学硕士学位论文1996,6(4):205-221.[58]TibshiraniR,SaundersM,RossetS,etal.Sparsityandsmoothnessviathefusedlasso[J].JournaloftheRoyalStatisticalSociety:SeriesB(StatisticalMethodology),2005,67(1):91-108.[59]GriffinJE,BrownPJ,GriffinJ.Alternativepriordistributionsforvariableselectionwithverymanymorevariablesthanobservations[J].WorkingPapers,2005:05-10.[60]ParkT,CasellaG.Thebayesianlasso[J].JournaloftheAmericanStatisticalAssociation,2008,103(482):681-686.[61]DelampadyM.ThomasBayes(1702–1761)[J].Resonance,2003,8(4):2-5.[62]AndrewsDF,MallowsCL.ScalemixturesofnormalDistributions[J].JournaloftheRoyalStatisticalSociety.SeriesB(Methodological),1974,36(1):99-102.[63]KyungM,GillJ,GhoshM,etal.Penalizedregression,standarderrors,andBayesianlassos[J].BayesianAnalysis,2010,5(2):369-411.[64]GeisserS.Thepredictivesamplereusemethodwithapplications[J].JournaloftheAmericanStatisticalAssociation,1975,70(350):320-328.60 致谢致谢研究生阶段的学习生活即将结束，回首在温州大学度过的点点滴滴，我得到了许多帮助和关怀，值此论文收稿之际，我谨向各位老师、同学及家人表示衷心的感谢。首先，我要感谢我的论文指导老师舒亮与闫正兵老师。从最初的论文选题与结构设计，到接下来的资料收集、算法设计和仿真以及论文的修改定稿，两位恩师都倾注了大量的心血。并且，导师们深厚的理论素养，渊博的学识，开拓了我的视野，使我的科研能力有了极大的提高。其次，我要感谢数电学院的钱祥忠老师、曾国强老师、张正江老师、戴瑜兴老师以及其他授课老师对我的悉心指导和帮助。另外，我还要感谢实验室的张海洲、张建师兄，曹婷婷、陈倩师姐以及我的室友褚金金、祝贵君同学在生活和学习上对我关怀和帮助。正是因为他们的陪伴，我研究生阶段的学习与生活才变得丰富多彩。最后我还要感谢家人对我的无私关爱与支持。61 温州大学硕士学位论文62 攻读硕士学位期间发表的论文攻读硕士学位期间发表的论文[1]张申波,闫正兵,吴平,张正江,戴瑜兴.基于LASSO的故障重构方法[J].计算机与应用化学,2016,33(11):1227-1230.[2]张申波,闫正兵,吴平,张正江,戴瑜兴.基于贝叶斯Lasso的工业过程故障诊断方法[C].昆明:中国过程系统工程年会,2017.[3]ShenboZhang,ZhengbingYan,PingWu,ZhengjiangZhang.FaultisolationbasedonBayesianfusedlasso[C]//ChineseAutomationCongress(CAC),2017.IEEE,2017:2778-2783.[4]ShenboZhang,ZhengbingYan,ZhengjiangZhang,Xiukai,Ruan.FaultdiagnosisinindustrialprocessesusingBayesiantheory[J].JournalofProcessControl,2017.（投稿）[5]闫正兵,张申波,吴平,张正江,张佳跃,一种基于贝叶斯理论的工业过程故障诊断方法.2017105670361,2017-10-20.（发明专利，实质审查）[6]戴瑜兴,朱志亮,谢晓青,张申波,曾国强,张正江,王环.基于GA-BP网络的液压制动系统多源融合故障预示方法.2016108978628,2017-03-08.（发明专利，实质审查）63 温州大学硕士学位论文64