一种相关向量机优化强重构边缘降噪自动编码器的轴承亚健康识别方法与流程

本发明涉及一种相关向量机优化强重构边缘降噪自动编码器的轴承亚健康识别方法，属于滚动轴承故障诊断领域。

背景技术：

现代工业水平的发展随着科技水平的提高得到长足发展，生产设备也由智能化代替了最初依靠人工的情况。在这个过程中，对于设备本身的要求也在不断提高，任意零部件的损坏都会对生产流程产生重要的影响，滚动轴承作为旋转机械元件的核心部件，也是机械设备中常用的元件之一。如果能准确识别轴承的当前状态，能在很大程度上避免造成的各种损失，由于轴承的故障并不是瞬间造成的，而是一个不断积累的过程，所以对轴承亚健康识别的研究受到了专家学者的广泛关注。

技术实现要素：

为了解决上述存在的问题，提供一种相关向量机优化强重构边缘降噪自动编码器的轴承亚健康识别方法。

本发明创造采用的技术方案为：其步骤为：

1)信号采样：将加速度传感器分别放置于电机壳和风扇的十二点钟方向，matlab处理16通道数据记录器采集的振动信号，12khz的采样频率下收集信号，驱动端信号的采样频率是48000s/s；

2)降噪处理：采用傅里叶变换对信号做预处理，变换之后用最小值–最大值归一化处理，使得数据在保持原有相对关系的基础上，数值波动范围稳定在[0,1]；

3)特征提取：确定边缘降噪自动编码器的隐含层数为三层，将预处理后的数据输入到强重构边缘降噪自动编码器中；在改进约束条件的作用下映射后提取特征，隐含层的输入为上一隐含层的输出加上预处理后的数据，计算隐含层实际输出；在输入输出基础上计算损失函数对权重和偏置的偏导数，更新权重和偏置；对比设置的损失函数和最大设定次数，到达期望值时停止，否则重新计算隐含层的输出；

4)故障诊断划分：使用改进的相关向量机进行故障诊断的划分；

5)机械设备状态诊断：利用实测的轴承数据，识别出轴承数据的状态。

所述的步骤3)中，具体方法如下：

3.1)预处理输入数据，设处理后的输入值及输出期望值为x；

3.2)随机初始化权值w和偏置b，使用交叉熵来调节参数w,b，定义损失函数：

其中：d是输入维数；

a表示批量划分后的大小；

ik表示每批次数据中第i次输入的第k个元素；

先求和输入向量再求和批量数据的每个样本，自动编码器约束条件如式(2)：

其中：n表示输入个数；

m表示隐含层神经元个数；

z表示隐含层特征表达式；

l表示交叉熵函数；

表示输入数据x的第d维的噪声干扰项的方差；

β表示权重系数；

计算方式如式(3):

其中：表示平均激活度；

ρ表示每层神经元的激活度；

3.3)定义前一隐含层的输出为x＇，则下一隐含层的输入为(x＇+x)；

3.4)计算通过隐含层处理后的输出层的实际输出z；

3.5)在输入输出的基础上，损失函数对w,b的偏导数，重建过程中使用标量形式，计算方式为：

其中，wqp表示节点q与p的权值；

wkp表示节点k与p的权值；

byp表示节点p的偏置；

bzk表示节点k的偏置；

对于每批量数据的第i个样本来说，表示第p隐层的输入，表示第k输出层的输入；

f'(x)是sigmoid的一阶导数：f'(x)＝f(x)[1-f(x)]

对损失函数式(1)求偏导得式(6)：

其中，wrs表示节点r与s的权值；

byr表示节点r的偏置；

3.6)按照式(7)规则更新参数值，η是学习率，w是权值，by是偏置：

3.7)对比设置的损失函数和最大设定次数，达到我们的期望值时就停止，否则返回3.6)继续训练。

3.8)将改进边缘降噪自动编码器的输出，输入到改进的相关向量机中分类。

所述的步骤4中，具体方法如下：

4.1)设有两类输入样本x1,x2，x1(x11,x12,...x1i),x2(x21,x22,...x2i)(i＝1,2,...l1,j＝1,2,...l2)，μ1,μ2是x1,x2映射后的均值向量，由(8),(9)算得：

类内散布矩阵由式(10),(11)解得：

其中，k(,)为核函数；

为k(，)所对应的非线性映射；

k1，k2，k11，k22，k12公式为：

输入样本映射到特征空间中的方差计算方式如式(17)(18)：

在最大熵方法的约束作用下，同样求解输入样本映射到特征空间中的方差本文的提出的权重系数λ的计算方式如式(19)：

4.2)按照改进核函数方式，定义相关向量机结构为：

权重w＝(w0,...,wn)，k(x,xi)表达式如式(21)：

在fisher准则约束下求解输入样本映射到特征空间中的方差在最大熵准则约束下求解输入样本映射到特征空间中的方差根据所提方法计算核函数权重系数λ，方式如式(22)：

4.3)超参数α初始化，权重协方差σ如式(23):

σ＝-h＝(φ^tbφ+a)^-1(23)

其中h为海森矩阵，并且y＝[y1,y2,λ,yn]^t，b＝diag(β1,β2,λ,βn)，βn＝σ[y(xn)]{1-σ[y(xn)]}，a＝diag(α0,α1,λ,αn)；

4.4)重新估值超参数，估算方式如式(24)：

μi是xi映射后的均值向量；

4.5)相关向量机分类时一般采用sigmoid概率函数，其表达形式如式(25)：

其中：y(x,wmp)＝φ(x)wmp，wmp由式(26)算得：

求解得到样本对类别的概率预测值p(t＝1|x)，p(t＝0|x)，值中较大值对应的类别便是x的识别结果。

本发明创造的有益效果为：在原本边缘降噪约束条件基础上加入稀疏性的约束条件，并且将隐含层的输入改为前一隐含层的输出加上预处理的数据，解决边缘降噪自动编码器约束条件有限、数据压缩能力弱、重构误差大，特征数据表达能力弱的问题。在特征样本数据中构造信息不同及分布不平坦时，使用两种核函数权重相加的方式代替原来单一核函数单独作用，利用提出的方法，在fisher准则和最大熵准则共同作用下求解输入样本映射到特征空间中的方差，确定核函数权重系数。由提出的计算方式求得权重系数，然后初始化参数、估计超参数，最后得到亚健康识别分类。利用实测的轴承数据，比较了相关向量机边缘降噪自动编码器(mdae-rvm)、相关向量机强重构边缘降噪自动编码器(imdae-rvm)和相关向量机优化强重构边缘降噪自动编码器(imdae-irvm)三种方法的亚健康识别准确率，相关向量机优化强重构边缘降噪自动编码器(imdae-irvm)相比其他两种方法，有较高的亚健康识别准确率，改进后的方法能正确识别轴承数据的状态，以此验证本发明的有效性。

附图说明：

图1是改进约束条件mdae拓扑结构图。

图2是改进隐含层输入的imdae拓扑结构图。

图3是方法imdae-irvm流程图。

图4a是正常状态驱动端时域振动信号。

图4b是内圈故障驱动端时域振动信号。

图4c是滚珠故障驱动端时域振动信号。

图4d是外圈故障驱动端时域振动信号。

图5a是傅里叶变换后的正常状态驱动端时域振动信号。

图5b是傅里叶变换后的内圈故障驱动端时域振动信号。

图5c是傅里叶变换后的滚珠故障驱动端时域振动信号。

图5d是傅里叶变换后的外圈故障驱动端时域振动信号。

图6a是部分样本的傅里叶变换前滚珠驱动端振动信号图。

图6b是部分样本的傅里叶变换后滚珠驱动端振动信号图。

图7a是部分样本的傅里叶变换前内圈驱动端振动信号图。

图7b是部分样本的傅里叶变换后内圈驱动端振动信号图。

图8a是部分样本的傅里叶变换前外圈驱动端振动信号图。

图8b是部分样本的傅里叶变换后外圈驱动端振动信号图。

图9a是部分样本的傅里叶变换前正常驱动端振动信号图。

图9b是部分样本的傅里叶变换后正常驱动端振动信号图。

图10是隐含层数随神经元个数的增加误差变化图。

图11是均方误差变化趋势图。

具体实施方式

一种相关向量机优化强重构边缘降噪自动编码器的轴承亚健康识别方法，其步骤为：

1)信号采样：将加速度传感器分别放置于电机壳和风扇的十二点钟方向，matlab处理16通道数据记录器采集的振动信号，12khz的采样频率下收集信号，驱动端信号的采样频率是48000s/s。

2)降噪处理：采用傅里叶变换对信号做预处理，变换之后用最小值–最大值归一化处理，使得数据在保持原有相对关系的基础上，数值波动范围稳定在[0,1]。

3)特征提取：确定边缘降噪自动编码器的隐含层数为三层，隐含层神经元个数为400个，将预处理后的数据输入到强重构边缘降噪自动编码器(imdae)中，在改进约束条件的作用下映射后提取特征，隐含层的输入为上一隐含层的输出加上预处理后的数据，计算隐含层实际输出，在输入输出基础上计算损失函数对权重和偏置的偏导数，更新权重和偏置，对比设置的损失函数和最大设定次数，到达期望值时停止，否则重新计算隐含层的输出。其中，隐含层神经元个数是通过测试调整能取得最好效果时所得到的。

具体方法如下：

3.1)预处理输入数据，设处理后的输入值及输出期望值为x。

3.2)随机初始化权值w和偏置b，使用交叉熵来调节参数w,b，定义损失函数：

其中：d是输入维数；

a表示批量划分后的大小；

ik表示每批次数据中第i次输入的第k个元素；

先求和输入向量再求和批量数据的每个样本，自动编码器约束条件如式(2)：

其中：n表示输入个数；

m表示隐含层神经元个数；

z表示隐含层特征表达式；

l表示交叉熵函数；

表示输入数据x的第d维的噪声干扰项的方差；

β表示权重系数；

计算方式如式(3):

其中：表示平均激活度；

ρ表示每层神经元的激活度。

3.3)定义前一隐含层的输出为x＇，则下一隐含层的输入为(x＇+x)。

3.4)计算通过隐含层处理后的输出层的实际输出z。

3.5)在输入输出的基础上，损失函数对w,b的偏导数，重建过程中使用标量形式，计算方式为：

其中，wqp表示节点q与p的权值；

wkp表示节点k与p的权值；

byp表示节点p的偏置；

bzk表示节点k的偏置；

对于每批量数据的第i个样本来说，表示第p隐层的输入，表示第k输出层的输入；

f'(x)是sigmoid的一阶导数：f'(x)＝f(x)[1-f(x)]

对损失函数式(1)求偏导得式(6)：

其中，wrs表示节点r与s的权值；

byr表示节点r的偏置。

3.6)按照式(7)规则更新参数值，η是学习率，w是权值，by是偏置：

3.7)对比设置的损失函数和最大设定次数，达到我们的期望值时就停止，否则返回3.6)继续训练。

3.8)将改进边缘降噪自动编码器的输出，输入到改进的相关向量机中分类。

4)故障诊断划分：使用改进的相关向量机进行故障诊断的划分；在强重构边缘降噪自动编码器的亚健康识别方法(imdae)的基础上，改进相关向量机，使用两种核函数权重相加的方式代替原来单一核函数单独作用，改进后核函数中权重系数在fisher准则和最大熵准则的约束下利用特征空间中样本方差解得，可以改善核函数处理数据的性能，提高识别准确率。

具体方法如下：

4.1)设有两类输入样本x1,x2，x1(x11,x12,...x1i),x2(x21,x22,...x2i)(i＝1,2,...l1,j＝1,2,...l2)，μ1,μ2是x1,x2映射后的均值向量，由(8),(9)算得：

类内散布矩阵由式(10),(11)解得：

其中，k(,)为核函数；

为k(，)所对应的非线性映射；

k1，k2，k11，k22，k12公式为：

输入样本映射到特征空间中的方差计算方式如式(17)(18)：

在最大熵方法的约束作用下，同样求解输入样本映射到特征空间中的方差本文的提出的权重系数λ的计算方式如式(19)：

4.2)按照改进核函数方式，定义相关向量机结构为：

权重w＝(w0,...,wn)，k(x,xi)表达式如式(21)：

在fisher准则约束下求解输入样本映射到特征空间中的方差在最大熵准则约束下求解输入样本映射到特征空间中的方差根据所提方法计算核函数权重系数λ，方式如式(22)：

4.3)超参数α初始化，权重协方差σ如式(23):

σ＝-h＝(φ^tbφ+a)^-1(23)

其中h为海森矩阵，并且y＝[y1,y2,λ,yn]^t，b＝diag(β1,β2,λ,βn)，βn＝σ[y(xn)]{1-σ[y(xn)]}，a＝diag(α0,α1,λ,αn)。

4.4)重新估值超参数，估算方式如式(24)：

μi是xi映射后的均值向量；

4.5)相关向量机分类时一般采用sigmoid概率函数，其表达形式如式(25)：

其中：y(x,wmp)＝φ(x)wmp，wmp由式(26)算得：

求解得到样本对类别的概率预测值p(t＝1|x)，p(t＝0|x)，值中较大值对应的类别便是x的识别结果。

5)机械设备状态诊断：利用实测的轴承数据，识别出轴承数据的状态。

实施例1：

据统计，滚动轴承70％以上的故障都是以振动形式表现出来。通常造成滚动轴承故障的原因有很多，但是只要能在其发生故障前发现问题及时采取措施就为时未晚，常见的故障问题主要有内圈故障、外圈及滚珠故障，对其“健康”状况的监测意义重大。因此，实验数据选用滚动轴承的实测数据，反映工作中的真实状态，通过对比预处理前后的数据，可能更清楚的看到各种故障状态下的数据分布情况。

提出的方法选择振动信号进行故障特征提取。针对滚动轴承信号中包含大量复杂的信号特征，传统的特征提取方法往往依靠人工来处理，同时对振动信号需要特殊的预处理过程，因此提出一种相关向量机改进强重构边缘降噪自动编码器的亚健康识别方法，改进原型自动编码器约束条件和隐含层的输入，将通过强重构边缘降噪自动编码器(imdae)重构后的数据输入到改进的相关向量机进行亚健康识别。能够很好的解决特征样本数据中构造信息不同及分布不平坦，以及使用单一核函数不能有效处理数据结构复杂的样本数据的问题。

一、本发明方案的理论依据：

1.自动编码器的提出

自动编码器作为一种非监督学习的神经网络模型，具有三层神经网络的结构，包含输入层，隐含层和输出层。它可以分为编码器和解码器两个部分。和其他的三层神经网络相比，形式是一样的，但是对输入和输出做了一部分限制，通过编码器部分，可以将高维的输入数据转换为低维的隐含层节点，又通过解码器部分，将低维的隐含层又重构为原来的输入数据，由于隐含层对输入信号进行了重构，所以通过隐含层可以得到数据的另一种表达。

自动编码器的运行步骤如下：

(1)输入x，同时将输出期望值设置定为x。

(2)构建编码器和解码器，随机初始化生成权值w和偏置b，b'。

(3)由给定的输入输出计算各个隐含层，输出层输出。

(4)利用输入和输出，计算误差函数计算各个神经单元的偏导数。

(5)更新权值w和偏置b、b'。

(6)对比设置的损失函数和最大设定次数，当达到期望条件时即停止训练。否则，选取下一个学习样本返回到第(3)步，进入下一轮的学习。

降噪自动编码器作为自动编码器的改进可以有效的去除原始特征的噪声，对特征提取有更好的效果。边缘降噪自动编码器通过边缘化处理噪声干扰的方式，克服降噪自动编码器在计算量、处理时间等方面的不足。最关键的部分在于损失函数的展开方式，用泰勒展开式的方式可达到近似模型的期望损失函数的目的，对于输入向量给隐含层神经元造成的敏感性问题也进行了充分的考虑。

强重构边缘降噪自动编码器改善了边缘降噪自动编码器数据压缩能力弱、约束条件能力有限，重构前后误差较大的问题，一方面在边缘降噪自动编码器的约束条件基础上加上稀疏性的约束条件，解决其在约束条件和数据压缩方面的不足，改进约束条件mdae拓扑结构如附图1所示；另一方面将隐含层的输入改进为前一隐含层的输出加上预处理后的数据，增强特征样本的表达能力，减小重构误差，改进隐含层输入的mdae拓扑结构附图2所示；最后用相关向量机亚健康识别强重构边缘降噪自动编码器(imdae)的输出。

2.改进imdae-rvm训练流程

imdae-rvm的操作流程如附图3所示，具体操作如下：

(1)利用傅里叶变换和归一化方法预处理输入数据。

(2)定义强重构边缘降噪自动编码器(imdae)的期望值，随机初始化权值和偏置，设置imdae的约束条件，预处理后的数据作为强重构边缘降噪自动编码器(imdae)的输入，在约束条件的作用下映射后提取特征，隐含层的输入为上一隐含层的输出加上预处理后的数据，计算隐含层实际输出，在输入输出基础上计算损失函数对权重和偏置的偏导数，更新权重和偏置，对比设置的损失函数和最大设定次数，到达期望值时停止，否则重新计算隐含层的输出。

(3)将经过强重构边缘降噪自动编码器(imdae)重构后的数据输入到改进核函数的相关向量机中进行亚健康识别，首先定义相关向量机的结构，核函数为高斯核函数乘以权重系数加上线性核函数乘以权重系数之和，对于核函数权重系数的确定，选择在fisher准则和最大熵准则的约束下，计算样本映射到样本空间中的方差，由提出的计算方式求得权重系数，然后初始化参数、估计超参数，最后得到亚健康识别分类。

本发明技术方案的实现过程：

1.信号采集：实验数据的采集工作是通过加速度传感器完成的，将传感器分别放置于电机壳和风扇的十二点钟方向，matlab处理数据记录器16通道采集的振动信号，12khz的采样频率下收集信号，驱动端信号的采样频率是48000s/s，数据分为0hp和3hp负载下，轴承处于正常状态，外环故障、内环故障和滚动故障状态，对应的转速情况分别是1797r/min和1730r/min,故障直径选择的0.007mm，列出四种故障的振动信号如附图4a-d所示。实验数据的选取情况如表1所示：

表1八组轴承实验数据参数

2.降噪处理：由于神经网络正常情况下不能处理时间序列，而且滚动轴承的振动信号在多数下是一种有规律的周期信号，所以采用傅里叶变换对滚动轴承的振动信号做预处理，转换后的信号如附图5a-d所示。再做相应的变换以后用最小值–最大值归一化处理，使得数据在保持原有相对关系的基础上，数值波动范围稳定在[0,1]，公式如式(27)：

这里，最小值mina、最大值maxa，v是属性a的原始、v'是归一化后的值。由上述方程可以看出，特征值的最大值和最小值分别对应于1和0。为了更清楚的看到傅里叶处理前后的差距，截取其中一部分实验的数据(1000个)进行对比，结果分别如附图6a-b、7a-b、8a-b、9a-b所示。

3.特征提取：针对实验模型拓扑结构的选择，对比一层、两层、三层隐含层数随神经元个数增加，样本集误差百分比的变化情况，如附图10所示，一定范围内，在隐含层数目不变的情况下随着神经元个数的增加，样本集的误差呈下降趋势，对于两层隐含层来说，当神经元数目达到500时误差最低，约为0.07，相比两层隐含层，三层隐含层在神经元数为400时达到最低点，误差约为0.03，因此，为保证最后实验的准确性，此次试验选择的边缘降噪自动编码器的隐含层数为三层，隐含层神经元个数为400个。

4.故障诊断划分：通过改进的rvm可以对不同故障进行更加有效划分,针对改进后的核函数，需要计算核函数的权重系数，提出一种根据特征空间中输入样本的方差计算权重系数的方法，表2中计算了fisher准则约束下特征空间中输入样本的方差和最大熵方法约束下的特征空间中输入样本的方差算得λ＝0.52，1-λ＝0.48。

表2不同约束下特征空间样本方差和权重系数

5.机械设备状态诊断：选取四次的实验结果计算均值对比分析了不同方法的识别准确率，表3为方法改进前后准确率对比，mdae-rvm的亚健康识别方法通过边缘降噪自动编码器(mdae)重构输入再使用相关向量(rvm)进行分类，四次实验亚健康识别准确率分别为85.6％、80.3％、82.7％、83.1％，四次均值是82.9％；imdae-rvm的亚健康识别方法先用imdae重构输入信息，再使用相关向量机(rvm)分类，四次实验的准确率分别为90.3％、90.7％、89.9％、91.4％，四次实验均值为90.55％，相比mdae-rvm方法，改进后的imdae-rvm方法在识别准确率上提高了7.65个百分点；imdae-irvm的亚健康识别方法是在强重构边缘降噪自动编码器的基础上改进单一的核函数，并用特征空间的样本方差求解改进后的核函数的权重系数，四次实验的识别准确率分别为94.2％、95.7％、94.6％、94.3％，四次均值为94.7％，与imdae-rvm相比，在准确率上分别提高了11.8和4.15个百分点。均方误差随迭代次数增加的变化趋势如附图11所示。

表3方法改进前后分类准确率对比

表4为几种不同改进方法在轴承亚健康识别的准确率对比，从表中可以看出，提出的rvm优化强重构mdae的亚健康识别方法(imdae-irvm)，较改进粒子群优化bp神经网络(ipso-bp)、改进蚁群方法优化支持向量机(iaco-svm)、改进遗传方法优化bp神经网络(iga-bp)在准确率上分别提高了7.2个百分点、2.2个百分点、3.46个百分点。

表4不同改进方法准确率对比

实验用到的数据为实测的轴承数据，采样频率12khz、轴承负载0hp、实验转速1797r/min的驱动端数据作为正常状态；采样频率12khz、轴承负载3hp、实验转速1730r/min的驱动端数据作为亚健康状态；采样频率12khz、轴承负载0hp、实验转速1797r/min、故障直径选择0.007mm的驱动端故障数据作为故障状态。二进制编码输出的数据，即正常状态为[100]，亚健康状态为[010]，故障状态为[001]。实验结果如表5，6所示：

表5实验结果对比表

从表5可以看出，测试样本的第二组数据，期望输出是[100]的正常状态，使用imdae-rvm的亚健康识别方法输出为[0.20620.65350.0001]，从数据可以看出应该划分为亚健康状态，这不符合期望输出，使用imdae-irvm的亚健康识别方法输出为[0.79030.01340.0988]，从数据可以看出划分为正常状态，符合期望输出，说明改进后的方法识别正确。

表6实验结果对比表

从表6可以看出，测试样本的第七组数据，期望输出是[001]的故障状态，使用imdae-rvm的亚健康识别方法输出为[0.02450.57830.1591]，从数据可以看出应该划分为亚健康状态，这不符合期望输出，使用imdae-irvm的亚健康识别方法输出为[0.03510.27690.8984]，从数据可以看出划分为故障状态，符合期望输出，说明改进后的方法识别正确。此外，观察测试样本的第五组数据，期望输出为[010]的亚健康状态，使用imdae-rvm的亚健康识别方法输出为[0.38260.44270.0032]，使用imdae-irvm的亚健康识别方法输出为[0.00010.63520.0093]，从数据可以看出，两种方法都能成功识别测试样本，但是imdae-rvm方法输出值中0.3826和0.4427相差不大，很容易出现分类不准确的问题，相比imdae-rvm方法，imdae-irvm方法的输出值[0.00010.63520.0093]中数值之间差异大，不容易出现误分类的情况。

通过对比表5和6发现，imdae-irvm方法较imdae-rvm方法能准确识别出测试样本的正确状态，并且数值之间大小差距明显，不会出现误分类的情况，证明了rvm优化强重构mdae的亚健康识别方法(imdae-irvm)的有效性。