基于时序关联矩阵的高炉冶炼过程多重关联时延估计方法

蒋 珂 蒋朝辉 ,2 谢永芳 潘 冬 桂卫华

高炉炼铁是一个连续鼓风、周期性加料和间歇性出铁、在高温高压高粉尘等恶劣环境下发生复杂物理化学反应、剧烈的物质与相变转化、高强度转移与传递的生产过程[1].全过程是在炉料自上而下、煤气自下而上的相互接触过程中完成的,高炉炼铁过程工艺如图1 所示[2].固体燃料(焦炭、煤粉等)、含铁原料(烧结矿、球团矿和块矿)和溶剂(白云石、石灰石和锰矿等)按一定配比从炉顶装入炉内,同时高炉下部喷煤系统喷吹的煤粉与热风炉加热高温热风中的氧气在风口平台发生燃烧反应产生一氧化碳和氢气等高温还原性气体.上升煤气流与下降炉料间进行的一系列传热、传质以及干燥、蒸发、挥发、分解、还原、软熔、造渣、渗碳、脱硫等物理化学变化后生成熔融铁水,当炉缸内的铁水到达一定容量后,使用开孔机钻开出铁口间歇性地排出熔融铁水,在撇渣器处实现渣铁分离后经鱼雷罐车运至后续生产工序[3-5].

图1 高炉三维仿真模拟图Fig.1 Three-dimensional simulation diagram of the blast furnace cast field

时延特性是所有流程工业,包括但不限于高炉冶炼过程的一个共有属性,主要是由炉料传输反应时间和冶炼单元在空间和时间分布上的差异所造成的.由于炉料无法被标记跟踪,导致很难确定高炉上部布料制度改变后多久才能影响到运行状态和炉缸内的铁水质量.此外,在高炉不同的空间位置上分布着不同的冶炼单元,比如高炉中下部连通的热风炉系统和喷煤系统,与喷煤和热风相关的控制变量对运行状态和铁水质量的影响也存在一定的时延,主要是因为上升的高温煤气流和下降的炉料发生物理化学反应需要一定的时间.

高炉冶炼过程时延的存在,不仅增加了运行状态在线监测和铁水质量在线建模的难度,也加大了运行状态和铁水质量调控的难度.当矿源改变或者外界干扰引起高炉冶炼过程稳态工作点改变时,一般需要依靠经验丰富的现场操作者基于经验进行调节,逐渐过渡到新的稳态工作点,这种动态调节较为粗放,主观性较强,并且高炉冶炼过程的滞后性使得动态调整的时间往往过长.长时间的炉况不稳定造成铁水成分波动大、铁水质量不合格,会给钢铁企业带来大量的能源和资源的浪费.因此,分析和准确估计高炉冶炼过程时延特性和参数,对运行状态和铁水质量的建模和优化控制具有重要的意义[6].

高炉的运行状态需要从原料、操作、出铁等方面考虑多个指标综合衡量,铁水质量就是其中一个重要的指标[7-10],因此本文主要聚焦于估计高炉冶炼过程的关键过程变量相对于铁水质量变量的时延参数.铁水温度和硅含量是表征铁水质量的重要参数,考虑到前期工作研发的红外视觉检测系统能实现出铁口铁水温度的实时在线检测[11],为此,本文工作主要是估计关键过程变量相对于铁水硅含量的时延参数,并为铁水硅含量的实时建模和优化控制提供重要的参考信息.关键过程变量与铁水硅含量的时延特性具有明显的多重性和关联性.其中第一重是由于炉料传输反应的时间和冶炼单元在空间和时间分布上的差异所导致的时延.第二重来自于高炉现场人工取样、冷却再化验过程和化验信息录入与审核下达存在不确定的时延.除此以外,高炉内部发生多达数百种物理化学反应使得各过程变量之间又互相影响,导致过程变量的多重时延又存在关联.

现有时延分析是根据专家经验结合相关系数得到最大相关性的某一确定滞后时间.比较常见的做法是单独地计算不同时延步长下的过程变量与质量变量的皮尔逊相关系数,最大相关系数所对应的步长即为时延时间[12-14],尽管这种方法有着便于实现、计算速度快等优点,但是,在计算时忽略了其他辅助过程变量对铁水硅含量的影响以及过程变量之间的相互影响,即忽略了时延之间的关联性,并且皮尔逊相关系数只能考虑两个变量之间的线性关系,并不适合非线性强耦合的高炉冶炼过程.在其他工业领域,如石油化工过程,文献[15]考虑了过程变量间的联合互信息来确定时延参数,但是基于互信息的时延联合估计在输入维数高、样本数量大的情况下往往运算时间较长,不利于实际的工程实践[15].本质上来说,无论是基于相关系数还是互信息的时延参数估计方法,都是在寻找与待估计变量相关性最大的过程变量.

现有的方法为过程变量之间的时延参数估计提供了很好的思路,但考虑到过程变量与硅含量之间的时延存在多重关联性,导致传统的单变量时延估计方法难以准确地估计出复杂的高炉冶炼过程真实时延.随着集散控制系统和工业互联网成功运用以来,高炉冶炼系统已积累了海量能反映冶炼过程的知识和数据,这些数据中隐含着过程的多重关联时延特性.因此,本文从数据驱动的角度出发,提出了一种基于时序关联矩阵的过程变量多重关联时延参数估计方法,采用多元时间序列描述不同时延系统的多重关联特性,运用时基标识不同时延的反馈信息,改变多元时间序列和时基,来寻找隐藏在数据中的真实时延参数.具体而言,以变量在某段时间稳定运行的数据为对象,根据过程变量的时延序列重构对应的时序关联矩阵,并引入灰色关联度来量化时序关联矩阵的多重关联相关性.

此外,考虑到单个变量的时延参数的取值ci为预估时延区间内采样周期整数倍个数,不失一般性,若估计N个过程变量相对于硅含量的时延参数,可能的时序关联矩阵将会有种可能.穷举计算所有可能时延组合下的时序关联矩阵的相关性是非常耗时的,因此在保证精度的前提下,以尽量小的计算开销,快速估计出最优的时延参数是十分必要的.基于此,本文提出了一种双尺度协同搜索策略的动态多群粒子群算法用来快速搜寻具有最大相关性的时序关联矩阵,将多重关联时延参数估计问题转化为求解一系列时序关联矩阵的最大相关性的优化问题.最后,根据估计的时延参数在时序上重新配准过程变量和硅含量时间序列来恢复真实的数据分布,提高了建模数据的准确性和一致性.数值仿真的结果和高炉冶炼过程定性的分析结果和定量的硅含量在线预测结果也进一步验证本文所提方法的可行性和有效性.

1 多重关联时延估计问题描述

高炉冶炼过程测量数据的准确性和有效性是过程建模、优化控制与性能评估的核心,作为一种典型的具有多重关联时延特性的生产过程,时延特性反映了冶炼过程的动态因果关系.因此,若不考虑高炉冶炼过程的多重关联时延问题,则关键过程变量与铁水硅含量之间实际的因果关系将被打乱.本文的主要工作是充分挖掘大量历史数据的分布规律,从中估计出过程变量相对于硅含量的多重关联时延参数.

为了实时监测高炉冶炼过程的运行状态,假设在高炉四周安装了如图2(a)所示的多个传感器,其中x0表示硅含量,xi表示第i个过程变量.假设样本实时且均匀地在时间轴上采样M次,若不考虑过程变量与硅含量之间的多重关联时延,那么在第t时刻与硅含量对应的过程变量为如图2(b)所示时间序列数据,即:

图2 高炉炼铁过程中变时滞问题描述Fig.2 Illustration of variable time-delay problem in the blast furnace ironmaking process

但是,由炉料传输反应时间和冶炼单元在空间和时间分布上的差异可知,不同的过程变量与铁水硅含量之间有着不同的时延,假设铁水硅含量与各过程变量之间真实的多重关联时延为Γ＇,即:

因此本文的目的就是基于高炉冶炼过程保存的大量的生产运行数据估计出最优的多重关联时延序列 Γ,使 得根据估计出来的时延参数,在时序上重构过程变量与硅含量之间的对应关系,消除时延后的时间序列数据如图2(c)所示,即:

为了准确地估计各辅助过程变量的时延参数,本文提出了一种基于时序关联矩阵的过程变量多重关联时延估计方法.引入灰色关联分析来定量描述时序关联矩阵的多重关联相关性,并提出一种基于双尺度协同搜索策略的动态多群粒子群算法用来快速寻找最优时延参数.

2 基于时序关联矩阵和灰色关联分析的多重关联时延参数估计

2.1 基于时基序列的时序关联矩阵重构

高炉冶炼过程中保存着大量的生产运行数据,这些数据隐含着过程变量与硅含量之间时延特性.为了分析时延的多重关联性,需要将待估计的所有过程变量归纳到一起进行综合考虑.定义铁水硅含量为零时延参考变量,基于时延参数在时序上重构过程变量数据形成对应的多元时间序列,重构的过程变量与硅含量组成能反映不同时延特性的时序关联矩阵.

假设在高炉冶炼过程中,从上到下在竖式容器四周安装了N个传感器,在时间轴上采样M次可以得到过程变量形成的原始时序数据矩阵X:

原始时序数据矩阵中元素的下标代表过程变量,上标代表采样时刻.其中x0为零时延参考变量铁水硅含量的时间序列数据,xi为第i个过程变量的时间序列数据.

为了避免多重关联时延估计过程中量纲对结果的影响,需要将原始时序数据矩阵X进行归一化处理,公式如下:

假设估计的零时延参考变量(铁水硅含量)与各过程变量的多重关联时延参数用时延序列表示,即:

原始时序数据矩阵中各过程变量的采样周期记为T,则第i个过程变量相对于零时延参考变量的时延参数需要满足:

式中,D是时基序列,di是第i个过程变量的时基,是一个无量纲的整数,这样可以保证过程变量的多重关联时延参数为采样周期的整数倍,这是为了确保过程变量在时序上重构的可行性.

经过归一化处理后的原始时序数据矩阵中隐藏着高炉冶炼过程的时延特性.为了估计出不同过程变量的多重关联时延参数,在时空维度上重构不同时基序列下的过程变量与铁水硅含量组成的时序关联矩阵,运用时基标识不同时延的反馈信息.其重构过程如下,首先从零时延参考变量中选取从第t时刻开始采样频率为T的F个数据构成时间序列

使得序列所含的数据至少能包含从高炉顶部布料到炉缸出铁的一个完整的冶炼周期,并且从采样时刻t开始选取的F个数据需要保证包含零时延参考变量的一个过渡态.

考虑到高炉是一个竖式的冶炼容器,铁水硅含量的采样点在垂直方向上等效于高炉底部,其余的过程变量xi监测点位于高炉底部及上方,根据上文定义的多重关联时延序列可知,相对于铁水硅含量,过程变量xi相对于零时延参考变量的滞后系数为τi,因此取过程变量xi从t-τi时刻开始的F个数据生成新的时间序列x˙i,即从t-diT开始的连续F个采样数据:

按照上述规则就可以根据时延序列重构出新的数据矩阵,即根据过程变量的时延参数在时空上重构的时序关联矩阵

2.2 基于灰色关联分析的时序关联矩阵相关性度量

时序关联矩阵将过程变量在不同时延参数下的时间序列数据归纳在一起,为了有效地描述多元时间序列数据之间的相关程度,必须对时序关联矩阵的相关性进行量化.灰色关联分析是多元时间序列关联分析中常用的方法,其通过灰色关联度来定量描述多组变量之间的变化趋势,多组变量之间的关联度越高,灰色关联度越大,反之亦然[16].此外,灰色关联分析对数据数量和分布规律没有要求[17].因此,本文采用灰色关联度来衡量多个时间序列数据之间的相关程度,即表示相互影响的过程变量时间序列数据与铁水硅含量时间序列数据之间的相关性.

令重构的时序关联矩阵中的零时延铁水硅含量序列为母序列,需要估计时延参数的多元过程变量为子序列,首先计算母序列与子序列对应元素之差,用矩阵A表示:

令a=min{aji},b=max{aji},则灰色关联系数矩阵可用γ表示:

式中,ρ是分辨系数,ρ越小表示分辨能力越大,根据现有参考文献可知ρ=0.5 计算的灰色关联系数矩阵更能反映实际的多元变量之间的相关关系[18].灰色关联系数矩阵描述的是多元子序列(过程变量) 与母序列(零时延硅含量) 在不同采样时刻之间的相关程度,为了描述子序列与母序列之间的总体关联程度,将灰色关联系数矩阵中的元素按列求取均值,即子序列相对于母序列的灰色关联度定义为:

为了定量描述重构的时序关联矩阵各过程变量与铁水硅含量之间的多重关联程度,采用灰色关联度之和R来描述多元过程变量与硅含量之间的总体相关性,其中R表示为:

灰色关联度之和R能在一定程度上度量时序关联矩阵中多元时间序列之间的相关性,R的值越大表明多元时间序列相关性越强,此时对应的多重关联时延 Γ 即为高炉冶炼过程的真实多重关联时延 Γ＇.

2.3 基于专家经验的过程变量时延区间预估

根据式(6)可知,每个过程变量xi相比于零时延参考变量硅含量都有一个时延参数.考虑到高炉四周安装的传感器和控制点位置不同,导致过程变量相对于铁水硅含量的时延参数是不一样的.如高炉顶部上料系统检测数据的时延明显就要比高炉中下部热风炉系统检测数据的时延要长.加上入炉矿源和冶炼参数的变化,每个过程变量的时延参数也处于动态变化的过程中,因此需要对各过程变量时延取值范围的上下界进行适当放宽,这样就可以得到过程变量xi的取值范围

由本文所研究的中国西南某钢铁厂2#2650 m3高炉的冶炼工艺特点可知,其冶炼过程是一个周期性布料,间歇式出铁的过程.根据现场的专家经验知识可知,下降的炉料与上升的热风在炉内发生1～2 h 物理化学反应后会打开炉缸下方的出铁口间歇性地排出铁水,考虑到这种情况下投入的炉料可能没有完全完成一个冶炼周期,因此对基于专家经验确定的时延区间进一步放宽到 1～4 h.再考虑到时延的多重性,即铁水硅含量现场采样、冷却后离线化验和信息录入还存在 1～1.5 h 的滞后,综合考虑认为高炉冶炼过程中各过程变量相对于零时延参考变量硅含量的时延区间预估为1～6 h.根据过程变量xi时延区间的取值范围可以计算出对应的时基区间为考虑到di只能取整数,记di的取值范围为因此过程变量xi时基的可能取值个数为:

需要注意的是,本文所有过程变量的时延区间都设置为一样的范围,过程变量xi时延区间下限为所有过程变量中时延最小值,上限为所有过程变量中时延最大值.尽管这样会放宽高炉中部或者底部检测到的部分过程变量的时延区间范围,增加时序关联矩阵的数量,但能在一定程度上避免由于专家经验的主观性导致部分过程变量时延区间预估不准,进而导致时序关联矩阵不全和估计的多重关联时延序列不准确的问题.

3 基于双尺度协同搜索策略的动态多群粒子群算法的最优时基序列搜索

时序关联矩阵的数量与过程变量和时基个数呈爆炸式增长,利用穷举计算的方法虽能保证得到最优的时基序列,但是该方法具有很大的计算开销和时间复杂度.为此,本文提出一种基于双尺度协同搜索策略的动态多群粒子群来快速搜寻高炉冶炼过程真实的多重关联时延,也就是基于时延参数重构的具有最大相关性的时序关联矩阵,将多重关联时延估计问题转化为求解一系列时序关联矩阵的最大灰色关联度之和的优化问题.

3.1 基于双尺度协同搜索策略的动态多群粒子群优化算法

粒子群算法(Particle swarm optimization,PSO)是一种经典的智能优化算法,具有结构简单、调整参数少、收敛速度快等优点[19],因此本文采用粒子群算法来搜索最优的多重关联时延参数.粒子群算法将群体中的鸟抽象为没有质量和体积的粒子,每个粒子包含两个向量: 位置向量和速度向量.位置向量代表了粒子在搜索空间的位置,即优化问题的一个备选解[20],速度向量代表了粒子在搜索过程中的运动状态.备选解的好坏由适应度函数来评价.在本文所研究的多重关联时延参数估计问题中,相对于零时延参考变量硅含量有N个过程变量的时延参数需要估计.因此,粒子群需要在N维空间中搜索最优时基序列,评价时基序列好坏的适应度函数为重构的时序关联矩阵的灰色关联度之和R.假设N维空间中第i个粒子的位置向量表示为其中其对应的速度向量表示为其中第i个粒子的速度和位置更新的数学描述如下:

式中,pbesti是第i个粒子在当前迭代次数之前的个体历史最优位置向量,gbest是种群中所有粒子在当前迭代次数之前的全局最优位置向量.ω是速度的惯性权重,c1和c2分别表示第i个粒子向pbesti和gbest靠近的权值系数,r1和r2表示在区间 [ 0, 1] 内均匀分布的随机数.

由式(20)和式(21)可以看出,经典的PSO 算法中粒子主要通过向全局最优粒子学习来更新自己的速度和位置,当全局最优粒子陷入局部最优时进化停滞,导致算法将过早收敛.因此为了解决PSO算法易于过早收敛的缺点,本文提出了一种基于双尺度协同搜索策略的动态多群粒子群优化算法(Dynamic multi-swarm particle swarm optimization based on double-scale collaborative search strategy,DMS-PSO-CS).首先,为了平衡粒子群的全局搜索能力和局部搜索能力,将粒子群分为两个尺度,即探索性子群和探测性子群.大尺度的探索性子群主要聚焦全局搜索能力,能快速定位到最优解的区域,具有一定的空间勘探能力.小尺度的探测性子群主要聚焦于局部精细搜索能力,能在最优解附近进行精确的搜索,具有较强的开采能力.其中,双尺度的粒子主要是通过粒子的惯性权重和最大速度来描述.较大的惯性权重可以提升粒子的跳跃能力,增强全局探索能力,而较小的惯性权重则可以精确搜索当前区域,增强局部探测能力.另一方面,粒子的最大速度vmax决定了粒子在迭代过程中的最大移动距离,较大vmin能获得较强的全局搜索能力但容易错过最优解,反之能获得较强的局部精细搜索能力但容易陷入局部最优解.

大尺度探索性子群的权重调整策略和最大速度分别为:

式中,ωl和vl,max分别为探索性粒子的惯性权重和最大速度,ωmax和ωmin分别为权重的上限和下限,kmax为最大的迭代次数.需要注意的是探索性子群的惯性权重是跟迭代次数有关的周期性振荡函数.大尺度探索性子群的速度向量更新公式如下:

式中,lbest是探索性子群的全局最优位置,位置向量通过式(21)更新.

小尺度探测性子群的权重调整策略和最大速度分别为:

式中ωs和vs,max分别为探测性粒子的惯性权重和最大速度,m ax(vmax) 和 m in(vmax) 分别为最大速度的上限和下限.需要注意的是探测性子群的惯性权重是跟迭代次数有关的单调递减函数,最大速度是跟迭代次数有关的线性下降函数.

其次,为了增加探测性子群的多样性,将探测性子群分为多个动态小子群,避免子群最优粒子陷入局部最优导致整个子群提前收敛的问题.并且,为了能快速地搜索到全局最优解,提出了一种探索性子群和探测性子群协同搜索策略,小尺度探测性子群粒子的位置更新不仅要根据上一时刻粒子本身的速度和粒子历史时刻最优位置决定,还要考虑到探测性子群的局部最优位置与探索性子群局部最优位置之间的关系.具体来说,当探索性子群局部最优位置的适应度函数优于探测性子群的局部最优位置对应的适应度函数值,探测性子群的粒子在更新的时候就会有选择地向探索性子群局部最优粒子方向靠近,精细搜索该区域提高收敛精度.小尺度探测性子群的速度向量更新公式如下:

式中,sbestj为探测性子群中第j个小子群的全局最优位置,f(·) 为粒子的适应度函数.将探测性子群分为多个小子群能增加空间搜索的广泛性,避免跟随全局最优粒子的位置而提前收敛,位置向量通过式(21)更新.

最后,迭代一定次数后对粒子群进行随机打乱和重组,这种随机重组策略可以看作粒子间周期性的信息交换,重组后能改变粒子的邻域拓扑结构,使得原本处于探索状态的粒子对搜索空间进行探测,处于探测状态的粒子对搜索空间进行探索.该策略能够增强粒子群的搜索空间,很大程度上避免算法陷入局部最优.综上所述,基于双尺度协同搜索策略的动态多群粒子群优化算法流程如下:

算法 1.双尺度协同搜索策略的动态多群粒子群优化算法

3.2 基于双尺度协同搜索策略的动态多群粒子群算法的最优时基序列估计

将本文所提的基于双尺度协同搜索策略的动态多群粒子群算法用来寻找最优时基序列,需要注意的是,考虑到时基di只能取整数,因此每个粒子的位置向量都被限制为正整数.基于双尺度协同搜索策略的动态多群粒子群算法寻找最优时基序列流程如图3 所示,具体步骤如下:

图3 基于DMS-PSO-CS 算法的时延参数估计框架Fig.3 Time-delay parameter estimation framework based on DMS-PSO-CS algorithm

步骤 1.根据专家经验和高炉冶炼过程机理分析,从高炉历史数据库中挑选出与铁水硅含量相关的过程变量组成原始数据集X=[x0,x1,···,xi,···,xN],其中为零时延参考变量铁水硅含量,xi为第i个对硅含量有影响的过程变量;

步骤 2.对采集的数据进行相关数据预处理,包括异常值剔除、缺失值处理、归一化处理,预处理之后的数据集记为

步骤 3.分析各过程变量相对于铁水硅含量的时延预估区间,确定时延序列τN]中各时延参数的取值范围

步骤 4.随机初始化双尺度协同搜索粒子群中粒子的参数,根据粒子的位置参数构造对应的时序关联矩阵并计算各粒子对应的适应度函数更新个体最优和子群最优;

步骤 5.判断算法是否达到最大迭代次数,如达到,则停止寻优并输出全局最优时基序列D,并基于此计算出时延序列 Γ;否则,执行步骤6;

步骤 6.判断算法运行次数是否为重组周期的整数倍,如果是,则对双尺度协同搜索粒子群进行随机打乱和重组,更新个体最优和子群最优并转到步骤7;否则,直接执行步骤7;

步骤 7.利用式(21)和(23)、式(21)和(25)分别更新探索性子群粒子和探测性子群粒子的速度和位置参数,并重复执行步骤 4～6 直到达到终止条件,输出全局最优时基序列D,并基于此计算出时延序列 Γ.

4 多重关联时延估计结果分析与验证

本文提出的多重关联时延估计方法的有效性和可行性首先在一个数值仿真的实例上进行了验证,接着在中国西南某钢铁厂2# 高炉上采集的数据进行了工业试验,基于本文方法估计出来的时延结果与现场的专家经验较为符合,并且时序配准后的数据也进一步提高了后续硅含量软测量模型的性能.

4.1 数值仿真

为了体现过程的非线性和时序性,数值仿真中数据的生成规则为:

式中,A∈R4×3是一个随机生成的矩阵,Γ∈R3×1是服从正态分布 N (0, 2) 生成的一个向量,f(·) 是一个矩阵函数,其基本元素为非线性的sin、cos、exp 和平方函数,ε∈R4×1和β∈R 是服从正态分布N(0,0.01) 生成的一个干扰项.对于生成的1000组数据,时延向量被设置为Γ＇=[τ′1,τ′2,τ′3,τ′4]=[1, 2, 3, 4],即将变量x1、x2、x3、x4分别相对于x0往前推 1、2、3、4 采样周期.变量时延区间设置为[0,10],因此对应的时基序列和重构的时序关联矩阵共有个.为了比较本文所提双尺度协同搜索动态多群粒子群优化算法的性能,经典的粒子群算法也被用来寻优比较.其中,两种算法粒子群规模都设置为30 个.但DMS-PSO-CS 算法中,大尺度的探索性子群中粒子个数设置为15 个,小尺度的探测性动态子群个数设置为3 个,其中粒子个数设置为5 个.根据文献[21]的分析可知,PSO算法权重线性下降时能取得较好的寻优效果[21],其调整策略为:

式中,ωmax=0.9,ωmin=0.4,c1=c2=2,最大迭代次数kmax设置为500.双尺度子群的ωmax、ωmin、c1、c2与PSO 保持一致,小尺度的探测性子群的速度上下限设置为 1 和-1,即 m ax(vmax)=1,min(vmin)=-1,重组周期Re设置为10.

此外,传统的单变量时延估计方法,皮尔逊相关系数(Pearson correlation coefficient,PCC)[22]和互信息相关系数(Mutual information coefficient,MIC)[23]也用来估计时延进行相应的对比,即,单独地计算不同时延步长下的x1～x4与x0的相关系数,最大相关系数所对应的步长即为时延时间.基于不同方法估计的时延结果见表1.由表1 可以看出,在单变量估计方法中,基于互信息系数估计的时延结果明显优于皮尔逊系数估计的结果,这也说明了基于线性关系估计时延的方法并不适合具有非线性特点的数据模式.相比于单变量时延估计方法,基于时序关联矩阵估计的时延参数与真实设置的时延更加一致,这也说明了本文所提方法的有效性.此外,本文所提的双尺度协同搜索粒子群算法在同样的解空间中能找到更准确的时延参数,也验证了本文所提优化算法在寻优能力上的优越性.

表1 数值仿真中基于不同方法估计的过程变量时延值Table 1 The estimated variable time-delay values based on different methods in numerical simulation

4.2 高炉冶炼过程多重关联时延估计

1)过程数据分析与预处理.根据高炉冶炼过程的工艺和现场专家经验确定了对铁水硅含量有影响的过程变量,其详细描述如表2 所示.为了估计选取的过程变量与铁水硅含量的时延参数,选取了2 # 高炉2020 年8 月1 日至2020 年12 月17 日连续生产的、24 个过程变量的采样时间序列数据,其中过程变量的采样周期大概为10 秒一次,铁水硅含量在一个班次(8 h)大概有 6～8 个化验数据,选取时间范围内过程变量有1 160 141 组,铁水硅含量数据有7 282 组.此外,为了定性地分析时延参数估计的结果准确性,与硅含量有着一定正相关关系的铁水温度数据也加入分析,其中铁水温度数据是由前期工作中研发的出铁口红外视觉检测系统实时在线检测的[11].需要说明的是,铁水温度的检测跟硅含量离线取样基本是在同一垂直位置进行的,因此铁水温度跟硅含量,在理论上只存在化验过程的时延,选取时间范围内铁水温度数据有172 352 组.

表2 基于不同方法估计的过程变量时延参数Table 2 The estimated process variable time-delay values based on different methods

为提高多重关联时延估计结果准确性,需要对数据集进行相关预处理,得到标准的和干净的数据提供给后续的流程.对于休风以及设备故障等原因造成的缺失数据直接删除,对于设备故障或者人工录入错误而导致的异常数据通过箱线图直接剔除[24].对于采样频率不一致导致过程变量与铁水硅含量数据时标不对齐的问题,首先将过程变量和硅含量数据集通过时标进行匹配,再以小时为单位进行均值化处理,预处理后的过程变量数据共有111 041 组,带标签的数据共有3 120 组,即为原始时序数据矩阵

2)过程变量多重关联时延估计结果定性分析.根据前面数据预处理分析可知,需要估计出表3 中选取的25 个过程变量相对于硅含量的时延参数,所以时延序列的长度N=25,即粒子的位置向量根据时延区间预估分析可知考虑到原始时序数据矩阵的采样周期T=1 h,时基区间为di ∈即可取值 1、2、3、4、5、6.

为了定性地分析估计的时延参数的准确性,将不同方法估计的结果与现场专家经验进行定性的比较,其中DMS-PSO-CS 和PSO 算法中的参数与数值仿真中的设置保持一致.基于不同方法估计的过程变量时延参数如表3 所示.理论上,铁水温度相对于硅含量只存在化验过程的时延,高炉现场从取样冷却到化验得到硅含量的百分比含量,一般需要1～1.5h 的时间间隔,考虑到原始数据矩阵是按小时采样的,因此过程变量的时基只能取小时的整数倍,理论上铁水温度相对于零时延硅含量只存在1h 滞后.基于时序关联矩阵和单变量估计方法寻找到的铁水温度时延参数都为 1 h,符合高炉冶炼现场的实际情况.根据现场专家经验,高炉中下部的煤粉喷吹系统对调控的滞后时间为 2～3 h,而基于互信息和PSO 算法估计出来的本小时实际喷煤量对铁水硅含量影响滞后时间都为 1 h,显然不符合现场的专家经验,尽管基于皮尔逊相关系数估计出来的本小时实际喷煤量的时延为 2 h,但是其估计的实际喷煤量对铁水硅含量影响滞后时间为 1 h,这显然也是与专家经验相违背的,而基于DMS-PSO-CS算法估计的本小时和上小时实际喷煤量对铁水硅含量影响滞后时间分别为 3 h 和 5 h,这与现场操作者的先验知识更符合.高炉中下部的热风炉系统一般需要 1 h 的滞后时间才能作用到高炉调控上,基于互信息和PSO 算法估计出来的热风温度对硅含量的滞后时间分别为 6 h 和 2 h,而热风压力对硅含量的影响滞后时间分别为 1 h 和 5 h,理论上热风压力与热风温度都是由热风系统控制的,因此与热风相关的变量估计的时延参数不应相差过大.相比而言,基于皮尔逊相关系数和DMS-PSO-CS 算法估计出热风温度和热风压力的时延更加符合现场的经验.总而言之,基于不同的方法估计出来的时延参数对于部分变量都有一定的合理性,但总体来看,本文所提的方法对于大部分变量估计的时延结果都比较符合现场的专家经验,这也进一步说明了本文所提的方法的有效性和可靠性.

表3 基于不同建模策略下的铁水硅含量软测量模型性能Table 3 Soft-sensor model performance of silicon content in molten iron based on different modeling strategies

3)过程变量多重关联时延估计结果定量分析.为了定量地验证上述多重关联时延参数估计结果的准确性,将未经时序配准的数据、基于上述四种方法估计出时延后配准的数据分别用于铁水硅含量软测量建模分析.本文的软测量模型选用了主流的数据驱动模型: 支持向量回归机(Support vector regression,SVR)[6]、随机权神经网络(Random vector functional-link network,RVFLN)[9]和深度网络(Deep network stacked by denoising autoencoder,SDAE)[25].考虑到本文的主要工作聚焦于铁水硅含量软测量建模前期的多重关联时延参数估计,因此软测量模型都是基本的模型结构并没有根据过程或者数据的特点进行相应改进.由于支持向量机和随机权神经网络都采用有监督的训练方式,因此采用2 820 组带标签的样本训练模型,300 组样本用于测试模型的性能.而深度网络采用半监督的训练方式,即深度模型预训练阶段,先用111 041 组过程变量无监督的训练模型参数,模型微调阶段,在预训练好的深度网络后加一层回归层,并用2 820组带标签的硅含量数据有监督的微调网络参数.计算机配置如下: Windows 10 (64-b)操作系统,32 GB内存,酷睿i7-9700 (3.0 GHz) CPU,Python 编程语言.

为了评价不同软测量模型的性能,本文采用均方根误差(Root mean squared error,RMSE)和平均绝对误差(Mean absolute error,MAE)进行定量衡量,定义如下:

式中,yt代表铁水硅含量实际化验值,代表模型预测的硅含量值,N为样本的个数.统计指标RMSE 和MAE 越小表示模型的性能越好.此外,为了满足工业现场的实时性需求,模型的训练时间也用来评估其性能.

为了验证不同方法估计的时延参数对硅含量预测性能的影响,实验组1、2、3、4 和 7 开展相关的对比实验.不同对照组的详细的实验对比结果如表3 所示,其中Tr 和Ts 分别代表训练集和测试集,深度网络的训练时间由无监督预训练和有监督微调两部分组成.由表3 可知,相比于过程变量无时延估计直接建立硅含量软测量模型,基于单变量法和时序关联矩阵确定的时延都能在一定程度上描述高炉冶炼过程的时延特性,提高硅含量预测模型的精度.在单变量时延估计方面,基于互信息估计的时延参数重构后的样本训练的硅含量模型相比于皮尔逊相关系数方法有一定的提升,说明了基于非线性关系估计出来的时延参数更符合真实的数据分布且有利于提高模型的性能.进一步对比本文所提出的基于时序关联矩阵和单变量的时延估计方法,可以看出,考虑了时延多重关联特性所确定的时延参数更能反映真实的输入-输出因果关系,这一点体现在实验组1、2、3、4 性能评价指标上,也验证了本文所提出的基于时序关联矩阵估计过程变量多重关联时延的可行性和有效性.此外,基于DMSPSO-CS 算法时序配准后的数据训练的硅含量软测量模型,在测试集上统计的RMSE 和MAE,相比于PSO 算法时序配准后的数据训练的模型,分别降低了0.0036 和0.0032,这也说明相比于直接用粒子群算法寻找的最优时延参数,基于本文所提的双尺度协同搜索策略的动态粒子群算法寻找到的最优时延参数,能更好地恢复真实的输入输出时序对应关系.

为了验证不同软测量模型对硅含量预测性能的影响,实验组1、5 和 6 比较了同一批数据下相关的对比实验.由表3 可以看出,相比于传统的机器学习算法,深度学习算法的性能有一定的提升,其中一个主要的原因是深度网络在无监督预训练过程中充分利用了无标签的过程变量数据,此外,这种无监督的训练方式迫使隐含层能学到有用的特征表示在输出层更好的重构输入,有助于挖掘数据中隐含的关系,进而提高模型的性能.尽管支持向量回归模型在小样本数据集上有着出色的性能,但是并没有在硅含量预测任务上表现出明显的优势,主要原因可能是这种浅层的模型结构并不能很好地描述具有强非线性的高炉冶炼过程.相比而言,随机权神经网络的非线性描述能力使得模型的性能有了一定的提升,但网络初始化权重的随机性和数据噪声的影响使得模型性能还有进一步提高的空间.在模型泛化能力方面,深度网络和随机权神经网络在测试集上的性能回落相比于训练集均在可接受的范围内,说明这两种模型能学习到隐藏在数据背后的规律并对新的样本展示出较好的预测精度.在模型训练时间方面,随机权神经网络在这三种算法中有着明显的优势,展示了其在工业过程实时性方面的巨大潜力,但模型的预测精度是其工业应用需要进一步考虑的问题.尽管深度网络两阶段的训练时间远远多于其他两种方法,但在工业应用现场时,无监督的预训练完全可以离线完成,只需使用带标签的样本有监督地在线微调预训练好的网络结构.因此,本文所提的方法在工业现场也具有一定的应用价值.

考虑到本文的主要目的是估计高炉冶炼过程的多重关联时延参数并提高后续建模数据的准确性和有效性,为了更直观地展示不同方法确定的时延参数对硅含量软测量模型的影响,图4～ 8 分别绘制了实验组1、2、3、4 和 7 的预测细节信息,即不同方法配准后的数据训练的软测量模型预测值曲线和实际化验曲线以及对应的预测误差.从图中可以看出,相比于没有经过时序配准的数据,使用单变量分析和时序关联分析矩阵进行过程变量时延参数估计后,软测量模型的准确度和跟踪能力均要显著优于无时延估计的模型,说明引入过程变量的时延信息能有效提高软测量模型的准确度.但是,相比于基于皮尔逊相关系数和互信息系数分别计算不同时延步长下的最大相关性确定的时延参数,考虑了时延多重关联特性的模型性能明显更好,模型预测误差大都集中在现场可接受的 [-0.1, 0.1] 范围之内.此外,对比图4 和图5 可以看出,经过DMS-PSOCS 算法时序配准后的数据训练的软测量模型数值和趋势跟踪的更好,样本预测误差更小,也进一步说明了本文所提的DMS-PSO-CS 寻优算法找到的时延参数序列,更能真实地反映模型输入跟输出之间的因果关系,通过提高建模数据集的质量进而提高模型的预测性能.

图4 基于DMS-PSO-CS 算法时延估计的铁水硅含量预测结果Fig.4 The prediction details of silicon content in molten iron with time-delay estimation based on DMS-PSO-CS algorithm

图5 基于PSO 算法时延估计的铁水硅含量预测结果Fig.5 The prediction details of silicon content in molten iron with time-delay estimation based on PSO algorithm

4)粒子群算法寻优性能对比分析.为了快速地寻找到最优的时延序列,本文在原始粒子群优化算法基础上提出了一种双尺度协同搜索策略的动态多群粒子群优化算法,为了对比两种算法在寻优性能上的差异,绘制了如图9 所示的寻优迭代曲线图,其适应度函数为本文优化的目标函数,即时序关联矩阵的灰色关联度.由图9 可以看出,经典的PSO算法在第26 次迭代时就能收敛,但是难以跳出局部最优值.而本文所提的DMS-PSO-CS 算法在第18 次迭代的时候就能达到相同的寻优效果,并且能多次跳出局部最优值,迭代到第296 次才收敛.这也说明了本文设计的双尺度粒子群能较好地在全局探索和局部探测能力之间进行很好的平衡,在保证当前最优的情况下通过周期性信息交换策略增加探索更优解的可能.从适应度函数来看,DMS-PSOCS 算法和原始PSO 算法最终得到的最大灰色关联度分别为18.109 和18.442.从寻优结果上也能看出,本文所提优化算法的有效性和可靠性.

图6 基于MIC 算法时延估计的铁水硅含量预测结果Fig.6 The prediction details of silicon content in molten iron with time-delay estimation based on MIC algorithm

图7 基于PCC 算法时延估计的铁水硅含量预测结果Fig.7 The prediction details of silicon content in molten iron with time-delay estimation based on PCC algorithm

图8 无时延估计的铁水硅含量预测结果Fig.8 The prediction details of silicon content in molten iron without time-delay estimation

图9 基于不同算法的寻优迭代曲线Fig.9 Optimization iteration curve based on different algorithms

5 结束语

高炉冶炼过程大时滞特性导致了变量之间存在多重关联时延,为了有效地估计变量间的时延参数,本文从冶炼过程保存的运行生产数据出发,提出了一种基于时序关联矩阵的过程变量多重时延参数估计方法.引入灰色关联分析来描述重构后的时序关联矩阵之间的相关性,将多重关联时延估计问题转化为求解一系列时序关联矩阵最大相关性的优化问题.为了降低穷举所有时序关联矩阵的计算开销和时间复杂度,提出了一种基于双尺度协同搜索策略的动态多群粒子群算法用于快速寻优,提出的双尺度粒子能很好地平衡粒子群的全局探索能力和局部探测能力,粒子之间的周期性重组操作能在一定程度上保证算法多次跳出局部最优.相比于传统的单变量时延估计方法,考虑了多重关联时延特性估计出来的时延参数能提高数据的准确性进而提高后续硅含量预测模型的精度.此外,估计出来的时延参数也能为铁水质量精细化调控提供重要的参考信息,避免现场操作者频繁调整或者过调整.

本文从数据的角度出发提出了一种多重关联时延估计方法,但高炉冶炼过程的动态性也是影响时延结果准确性的一个重要性因素,因此炉况变动下时延参数的动态调整也是下一步需要考虑的问题.此外,尽管本文提出的时延估计方法能提高数据的质量进而提高后续硅含量模型的预测精度,但是相比于最新的硅含量软测量结果还有进一步提升的空间,因此下一步可以考虑将时延参数作为软测量模型的参数,模型的性能作为时延参数准确性的评价指标,进而将后续的软测量建模与时延估计有机结合起来,在提高模型性能的同时提高时延估计结果的准确性和可靠性.