ISSN 
-
重型数控机床对现代军工、航天、船舶、轨道交通、发电设备等重工业发展有重要意义[1]。目前，研究人员大多借鉴中小型数控机床的可靠性模型对重型数控机床进行可靠性分析[2]。重型数控机床具有加工工况复杂、故障溯源困难、维护成本较高、故障样本数据匮乏等特点，这对重型数控机床的可靠性评估及后续研究都造成极大困难[3]。小样本数据下的重型机床可靠性评估，是目前可靠性研究的一个难题[4]。
数控机床的可靠性研究最早开始于20世纪70年代[5]，前苏联专家对建立机床参数模型和机床工艺可靠性方面进行了相关研究[6]。近些年来，一些专家将贝叶斯理论引入小样本条件下的重型数控机床可靠性建模[7]，根据先验信息和主观经验进行模型参数估计[8]。文献[9]采用威布尔双参数和三参数分布对重型数控机床零部件或系统建模[9]。文献[10]采用极大似然估计法和Edgeworth级数法对重型数控机床的故障模型进行参数估计，并解释参数估计结果与数据拟合的关系。文献[11]利用第二类极大似然方法得到折合因子的贝叶斯估计[11]。文献[12]在马尔科夫链蒙特卡洛方法(Markov Chain Monte Carlo, MCMC)方法的基础上针对数据的多源性，建立多源异种数据融合模型，为小样本的可靠性建模提供新思路。上述研究对小样本数据下的可靠性建模问题提出了不同的方法，给重型数控机床可靠性建模及评估奠定了理论基础。
本文基于改进的贝叶斯方法对重型数控机床进行可靠性建模与评估。首先，选取双参数的威布尔分布对机床故障数据进行拟合，采用极大似然参数估计法和贝叶斯参数估计法对可靠性模型求解。其次，对传统的贝叶斯估计法进行改进，得到分层贝叶斯模型，通过对比分析及仿真实验，证明改进后的贝叶斯模型在小样本数据条件下更具有优越性。最后，利用分层贝叶斯模型对该机床进行可靠性评估，得到其理论故障间隔时间。
-
一般地，对系统建立可靠性模型均采用双参数威布尔分布。以重型数控机床故障间隔时间为参数，经过变换得到服从威布尔分布的可靠性模型。此时，机床可靠性服从威布尔分布的概率密度函数和概率分布函数分别为：
$$f(t) = \frac{m}{\eta }{\left( {\frac{t}{\eta }} \right)^{m - 1}}\exp \left[ { - {{( {\frac{t}{\eta }} )}^m}} \right]$$ (1) $$F(t) = 1 - \exp \left[ { - {{( {\frac{t}{\eta }} )}^m}} \right]$$ (2) 式中，m和η分别为威布尔分布的形状参数和尺度参数；t为机床的故障时间。
-
一般地，采用极大似然估计对分布函数进行参数估计。当已知分布模型中有未知参数
$ {\theta _1},{\theta _2}, \cdots ,{\theta _k}$ 时，可以建立似然函数表达式：$$L' = \prod\limits_{i = 1}^n f \left( {{t_i};{\theta _1},{\theta _2}, \cdots, {\theta _k}} \right) $$ (3) 式中，
$ {{t}}_{{i}} $ 为实际采集的数据子样(i=1, 2, ···, n)。在威布尔分布模型中，未知参数为m和η，将式(1)的威布尔分布概率密度函数代入式(3)中，建立的极大似然函数为：
$$L'\left( {{t_i},m,\eta } \right) = \prod\limits_{i = 1}^n f \left( {{t_i}} \right) = \prod\limits_{i = 1}^n {\frac{m}{\eta }} {( {{t_i}} )^{m - 1}}\exp \left[ { - {({\frac{{{t_i}}}{\eta }})^m}} \right]$$ (4) 为简化问题，对上式取对数处理，得到新的似然函数：
$$L'\left( {{t_i},m,\eta } \right) = n\ln ( {\frac{m}{\eta }} ) + \left( {m - 1} \right)\sum\limits_{i = 1}^n {\ln } \left( {{t_i}} \right) - \frac{1}{\eta }\sum\limits_{i = 1}^n {{t_i}^m} $$ (5) 为了求解
$ L{\rm{'}}({t_i},m,\eta )$ 的极值，令$ {\rm{d}}L/y({\theta _i}) = 0$ ，建立方程组：$$ \left\{ {\begin{array}{*{20}{c}} {\dfrac{{{\rm{d}}L}}{{{\rm{d}}\eta }} = - \dfrac{n}{\eta } + \dfrac{1}{{{\eta ^2}}}\displaystyle\sum\limits_{i = 1}^n {{t_i}^m} = 0}\\ {\dfrac{{{\rm{d}}L}}{{{\rm{d}}m}} = \dfrac{n}{m} + \displaystyle\sum\limits_{i = 1}^n {\ln } \left( {{t_i}} \right) - \dfrac{1}{\eta }\displaystyle\sum\limits_{i = 1}^n {\left[ {{t_i}^m\ln \left( {{t_i}} \right)} \right]} = 0} \end{array}} \right. $$ (6) 由式(6)可以求出m和η的最大似然估计值。
通过式(7)与式(8)，对待估参数处于置信水平为
$ {1}-\alpha $ 下的区间估计值进行求解：$$\left( {m - {z_{\frac{\alpha }{2}}}{\rm{se}}(m){\rm{,}}\;\;m + {z_{\frac{\alpha }{2}}}{\rm{se}}(m)} \right)$$ (7) $${\left( {\eta - {z_{\frac{\alpha }{2}}}{\rm{se}}(\eta ){\rm{,}}\;\;\eta + {z_{\frac{\alpha }{2}}}{\rm{se}}(\eta )} \right)}$$ (8) 式中，
$ {{z}}_{\alpha /2} $ 为标准正态分布的$ \alpha $ /2位点；$ {{\rm{se}}(m)} $ 与${{\rm{se}}(\eta )} $ 为参数的渐进分布标准差，具体求解过程不进行赘述。 -
由贝叶斯参数估计理论可知，在对分布函数的参数估计中，所有未知参数都被认为是随机变量。求解待估参数值时，需要先根据观测数据和主观经验设定其模型的先验分布，进而求出参数估计结果。其中，先验分布的选择优劣会直接影响模型参数估计的精度。
基于贝叶斯理论对重型数控机床建立可靠性模型，具体思路为在研究对象样本数据和先验信息的基础上，建立其贝叶斯模型，计算模型参数的后验分布并得出相应的概率分布。
结合贝叶斯定理分析可知，分布模型参数的后验分布可以表示为：
$$P(\left. \theta \right|y) = \frac{{P(\left. y \right|\theta )P(\theta )}}{{P(y)}} \propto P(\left. y \right|\theta )P(\theta )$$ (9) 式中，θ为待估计的参数；y为基于观测数据的信息；
$ {P}{(}{\theta }{)} $ 为参数$ {\theta } $ 的先验概率；$ P(\left. y \right|\theta )$ 是设定参数θ条件下的数据所服从的后验分布；$P(\left. y \right|\theta )/P(y)$ 为事件y对参数设定的支持程度，即似然函数。因为正态分布是威布尔分布的共轭先验分布模型，所以设定威布尔分布的参数m和η服从正态分布，这里用N来表示。
经过推导，得到重型机床基于贝叶斯理论的可靠性模型，如式(10)所示：
$$\begin{split} &f(t|m,\eta ) = \frac{m}{\eta }{( {\frac{t}{\eta }} )^{m - 1}}\exp \left[ { - {{( {\frac{t}{\eta }} )}^m}} \right]\\ &\;\;\;\;\;m \sim {{N}}(\alpha 1{\rm{,}}\;\beta 1){\rm{,}}\;\eta \sim {{N}}(\alpha 2{\rm{,}}\;\beta 2) \end{split}$$ (10) -
为了优化贝叶斯参数估计的方法，将待估参数的先验分布进一步分化，得到层次贝叶斯模型。层次贝叶斯模型理论上可以分为N层，但随着分化层数的增加，模型的计算效率逐渐降低。为了保证计算效率，采用两层贝叶斯模型进行分析。
由一般贝叶斯模型可知，
$m \sim {{N}}(\alpha 1{\rm{,}}\;\beta 1) $ 为参数m的先验分布，$\eta \sim {{N}}(\alpha 2{\rm{,}}\;\beta 2)$ 为参数$ \eta $ 的先验分布。将$ \alpha 2 \sim {\rm{unif}}(\alpha 1'{\rm{,}}\;\beta 1'){\rm{,}}\;\beta 2 \sim {\rm{unif}}(\alpha 2'{\rm{,}}\;\beta 2')$ 分别作为参数m , η的先验分布，这样就得到了更进一层的贝叶斯模型，这里称$ \alpha 2{\rm{,}}\;\beta 2 $ 为超参数，unif表示均匀分布，表达如下：$$ \begin{split} &f(t|m,\eta ) = \frac{m}{\eta }{( {\frac{t}{\eta }} )^{m - 1}}\exp \left[ { - {{( {\frac{t}{\eta }} )}^m}} \right]\\ &\;\;\;m \sim {{N}}(\alpha 1{\rm{,}}\;\beta 1){\rm{,}}\;\eta \sim {{N}}(\alpha 2{\rm{,}}\;\beta 2)\\ &{\alpha 2 \sim {\rm{unif}}(\alpha 1'{\rm{,}}\;\beta 1'){\rm{,}}\;\beta 2 \sim {\rm{unif}}(\alpha 2'{\rm{,}}\;\beta 2')} \end{split} $$ (11) 经过上式推导，即可得到基于层次贝叶斯理论的可靠性模型。
-
求解贝叶斯模型的参数需要计算高维积分的边缘后验分布，而且随着模型组成的复杂会使参数求解变得更加困难。通过引入 MCMC方法对模型进行数据模拟、迭代以提升参数求解的效率。该方法中包含多种抽样方法，包括Metropolis-Hastings 抽样方法和 Gibbs 抽样方法。其中，Metropolis-Hastings抽样方法更为高效，具体操作流程如图1所示。
图 1 蒙特卡洛仿真流程
通过MCMC法处理得到模型待估参数的后验分布，取其期望值作为参数估计的结果，基于排序取分位数法，可以求出待估参数的置信区间。
-
为了分析参数估计结果的优劣，本文采用统计学中的均方根误差值(normalized root mean square error, NRMSE)，NRMSE值越小则结果精度更优。计算方法为：
$${\rm{NRMSE}} = \sqrt {\frac{{\displaystyle\sum\limits_{i = 1}^n {{{(F({{{t}}_{{i}}}) - F(t_{{{\rm{ies}}}}))}^2}} }}{{\displaystyle\sum\limits_{i = 1}^n {{F^2}} ({t_i})}}} $$ (12) 式中：F(ti)是机床的累计故障概率，
$F(t_{{\rm{ies}}}) $ 是将参数估计值代入累积故障概率函数后所得到的累积故障概率值。 -
重型机床故障数据样本量较小、数据跨度较大、故障发生具有不确定性，若采用一般统计学方法无法对此类问题进行有效分析，本文基于贝叶斯理论进行分析，可以一定程度优化分析结果。
-
本文采集某企业近3年内，某型号重型数控机床的停机维修记录共58条。通过数据预处理，得到该机床的故障间隔时间，其中部分结果为22 h，23 h，24 h，···，912 h，936 h，1248 h。
-
将故障间隔时间代入机床的双参数威布尔分布可靠性模型中，利用最大似然参数估计法和贝叶斯参数估计法分别对分布模型的参数进行求解。采用最大似然参数估计法得到的结果如表1所示，用下标MLE来表示。
表 1 最大似然参数估计结果
参数估计结果 置信区间 mMLE=0.9967 0.8222～1.2083 ηMLE=236.9915 180.2095～311.6648 将机床的故障间隔时间从大到小排列，代入式(13)所示的中位秩公式：
$$F({t_i}) = \frac{{i - 0.3}}{{n + 0.4}}\;\;\;\;i = 1,2, \cdots ,n$$ (13) 计算该重型机床实际的累计失效概率密度，得到其点位分布图，如图2所示。
图 2 样本的数据点分布图
将参数估计结果代入机床的可靠性模型，得到基于威布尔分布的机床拟合累计失效概率曲线，如图3所示，图中虚线为95%置信区间的上下限。通过比较两条曲线的重合度来分析模型的拟合精度。
图 3 最大似然估计下的威布尔模型分布函数
由图3可知，通过极大似然估计得到的机床可靠性模型分布曲线，部分与参考数据点位重合，没有重合的部分也基本分布在置信区间内。
基于贝叶斯理论得到的参数估计结果如表2所示，其中，用下标GBR表示一般贝叶斯参数估计的结果，用下标HBR表示分层贝叶斯参数估计的结果。
表 2 一般贝叶斯和层次贝叶斯参数估计结果
模型 参数估计结果 95%置信区间 一般贝叶斯 mGBR=0.9808 mGBR=(0.8411～1.1670) ηGBR=208.80 ηGBR=(124.31～286.30) 层次贝叶斯 mHBR=0.9617 mHBR=(0.8248～1.1470) ηHBR=206.60 ηHBR=(135.57～283.57) 通过对参数进行多层迭代，输入不同的初始值，构造马尔科夫链然后观察其抽样迭代过程是否收敛，运算迭代过程如图4、图5所示。将
$ \alpha 1 = 1{\rm{,}}\; \alpha 2 = 200{\rm{,}}\;\beta 1 = \beta 2 = 0.5 $ 代入式(10)，$\alpha 1 = 1{\rm{,}}\;\beta 1 = 0.5{\rm{,}}\; \alpha 1' = \alpha 2' = 0{\rm{,}}\;\beta 1' = 100{\rm{,}}\;\beta 2' = 1 $ 代入式(11)，得到待估参数的后验分布，如图6、图7所示。
图 4 基于一般贝叶斯模型的仿真迭代轨迹图
图 5 基于分层贝叶斯模型的仿真迭代轨迹图
图 6 基于一般贝叶斯模型的参数后验分布概率密度函数
图 7 基于分层贝叶斯模型的参数后验分布概率密度函数
从待估参数的后验分布图中可以看出，分层贝叶斯参数估计的结果更加集中。将参数估计结果分别代入式(1)、式(13)中，得到图8、图9所示的分布曲线。由此可知，贝叶斯模型与传统的模型相比，其拟合曲线与真实数据的重合度更高。
图 8 一般贝叶斯估计下的曲线分布
图 9 分层贝叶斯估计下的曲线分布
-
将机床的3种可靠性模型代入式(12)，计算得到NRMSE值，结果分别为MLE=0.0971，GBR=0.1210，HBR=0.0900。对比结果可知：分层贝叶斯模型相对于其他两个模型来说误差值更小，数据拟合效果更优；基于极大似然估计法模型的NRMSE值低于一般贝叶斯模型的NRMSE值。虽然此处传统可靠性模型精度比一般贝叶斯模型更优，可是随着经验信息的丰富，优化参数先验信息的设定后，可以逐渐改善参数估计的精确性。
-
改变数据样本的容量，分别使用上述3种参数估计方法，经过计算求出不同样本容量下的标准均方根误差值，结果如图10所示。
图 10 不同样本量下的建模NRMSE值
从图中可以看出，随着样本容量的增加，经过3种参数估计方法得到的可靠性模型NRMSE值均在下降，说明建模误差在降低，精度越来越好，当数据量趋于大数据样本时，贝叶斯方法依然有良好的可靠性。但当数据量较小时，其NRMSE值更小，分层贝叶斯参数估计法有显著优势。特别是在样本量趋于0时，基于极大似然估计建模得到的NRMSE 值近似于1，而基于贝叶斯理论建模得到NRMSE 值依然小于1，说明在样本量很小时，传统的方法几乎失去精确性，而基于贝叶斯理论得到的结果依然有效。综上所述，在数据量不足时，分层贝叶斯模型的拟合效果更优。
-
基于matlab 软件，将服从形状参数Beta=2及尺度参数Eta=1的双参数威布尔分布函数随机筛选10组数据，把产生的小样本量代入到3种参数估计方法中，得到表3的结果。
表 3 仿真估计结果对比
参数 最大自然估计 贝叶斯估计 分层贝叶斯估计 真值 Beta 1.5270 1.5230 1.5360 2 Eta 0.7828 0.7705 0.8586 1 经过仿真分析也可以验证分层贝叶斯在进行小样本参数估计方面的优势，其估计的参数相对来说最接近真实值。
-
平均故障间隔时间(mean time between failures, MTBF)是常用的可靠性评价指标之一，是指产品从发生故障到下一次发生故障的平均工作时间，其表达式为：
$${\rm{MTBF}} = \int_0^\infty {tf} (t){\rm{d}}t = \eta \varGamma ( {1 + \frac{1}{m}} ) $$ (14) 式中，f(t)为概率密度函数。
将层次贝叶斯参数估计结果代入式(14)后，求得机床平均故障间隔时间为
$ \text{MTBF=210.215}\;\text{h} $ 。分析可知，计算得到的MTBF值小于所采集数据中的最大值1248.000 h，210.215 h大于所采集数据中的60%，经过排序分析处于中上位置。由此可见，平均故障间隔时间对于机床预防故障和维修保养有参考意义，说明了贝叶斯理论对于重型机床可靠性建模的有效性。 -
本文在重型数控机床故障数据小样本条件下，基于贝叶斯理论建立服从威布尔分布的可靠性模型，并对传统贝叶斯模型进行改进，得到分层贝叶斯模型。引入MCMC算法求解贝叶斯模型的参数估计结果，对比极大似然估计法与贝叶斯参数估计的优劣，结果显示这两种参数估计结果近似，但贝叶斯参数估计法得到的参数置信区间更窄。采用NRMSE分析参数估计结果精度，随着样本数据容量增大，传统参数估计法与贝叶斯参数估计法趋于一致，当数据样本较小时，分层贝叶斯模型得到的建模精度更高。对贝叶斯理论模型进行仿真实验，结果表明贝叶斯参数估计是有效的。利用贝叶斯可靠性模型对该重型数控机床进行可靠性评估，在机床故障数据小样本条件下，其评估结果更精确有效，对重型数控机床的可靠性研究有重要意义。
Reliability Research of Heavy CNC Machine Tools Based on Improved Bayesian
-
摘要: 重型数控机床在机械加工领域占据重要地位，因此提高其可靠性以及加工精度，对我国工业发展有重要意义。重型数控机床具有结构复杂、故障溯源困难、样本少、数据不足等缺点，因此对其进行可靠性研究比较困难。针对这一问题，采用双参数的威布尔分布建立机床的可靠性模型，引入贝叶斯理论对其进行参数估计，并通过马尔科夫链蒙特卡洛方法(MCMC)计算参数估计结果。对贝叶斯参数估计法中的待估参数进一步分析，得到多层次的贝叶斯模型，并通过参数仿真实验分析其准确性。采用标准均方根误差值及置信区间宽度进行模型优劣的对比，结果表明，改进后的贝叶斯方法参数估计结果精度更优，更有利于建立机床可靠性模型。Abstract: Heavy-duty CNC machine tools occupy an important position in the field of machining, and improving its reliability and machining accuracy is of great significance to the industrial development of China. Compared with conventional machine tools, heavy-duty CNC machine tools have the characteristics of complex structure, difficulty in fault tracing, few samples, and insufficient data, which make it difficult to conduct reliability research on them. Aiming at this problem, this paper uses the Weibull distribution of two parameters to establish the reliability model of the machine tool, introduces Bayesian theory to estimate its parameters, and calculates the parameter estimation results through the Markov Chain Monte Carlo method (MCMC). In order to improve the accuracy of parameter estimation, the traditional Bayesian method is improved, and the standard root mean square error value and confidence interval are used for evaluation and comparison. The results show that the improved Bayesian method has better parameter estimation accuracy and is more conducive to the establishment of machine tool reliability models.
-
表 2 一般贝叶斯和层次贝叶斯参数估计结果
模型 参数估计结果 95%置信区间 一般贝叶斯 mGBR=0.9808 mGBR=(0.8411～1.1670) ηGBR=208.80 ηGBR=(124.31～286.30) 层次贝叶斯 mHBR=0.9617 mHBR=(0.8248～1.1470) ηHBR=206.60 ηHBR=(135.57～283.57) 
下载: 导出CSV
表 3 仿真估计结果对比
参数 最大自然估计 贝叶斯估计 分层贝叶斯估计 真值 Beta 1.5270 1.5230 1.5360 2 Eta 0.7828 0.7705 0.8586 1
下载: 导出CSV
-
[1] 彭卫文. 重型数控机床可靠性建模与评估技术研究[D]. 成都: 电子科技大学, 2016. PENG W W. Research on methods for reliability modeling and assessment of heavy-duty CNC machine tools[D]. Chengdu: University of Electronic Science and Technology of China, 2016. [2] JIANG G J, LI Z Y, QIAO G, et al. Reliability analysis of dynamic fault tree based on binary decision diagrams for explosive vehicle[J]. Mathematical Problems in Engineering, 2021, 2021: 1-13. [3] ASLAM M, BALAMURALI S, PERIYASAMYPANDIAN J, et al. Designing of an attribute control chart based on modified multiple dependent state sampling using accelerated life test under weibull distribution[J]. Communications in Statistics Simulation & Computation, 2021, 50(3): 902-916. [4] 张英芝, 牛序磊, 申桂香, 等. 基于竞争失效的数控机床可靠性建模[J]. 系统工程理论与实践, 2014, 34(8): 2144-2148. doi: 10.12011/1000-6788(2014)8-2144 ZHANG Y Z, NIU X L, SHEN G X, et al. Reliability modeling of CNC machine tools based on competition failure[J]. Systems Engineering-Theory & Practice, 2014, 34(8): 2144-2148. doi: 10.12011/1000-6788(2014)8-2144 [5] 刘以平. 重型数控机床精度可靠性建模的研究[D]. 成都: 电子科技大学, 2015. LIU Y P. Research on the accuracy reliability model of heavy NC machine tool[D]. Chengdu: University of Electronic Science and Technology of China, 2015. [6] 邓超, 吴军, 毛宽民, 等. 面向大型数控机床的工艺可靠性评估[J]. 计算机集成制造系统, 2010, 16(10): 2250-2256. doi: 10.13196/j.cims.2010.10.238.dengch.026 DENG C, WU J, MAO K M, et al. Process reliability evaluation for large-sized NC machine tools[J]. Computer Integrated Manufacturing Systems, 2010, 16(10): 2250-2256. doi: 10.13196/j.cims.2010.10.238.dengch.026 [7] 韩凤霞. 高端数控机床服役过程可靠性评价与预测[D]. 北京: 机械科学研究总院, 2020. HAN F X. Reliability evaluation and prediction of high performance CNC machine tools in service[D]. Beijing: China Academy of Machinery Science and Technology, 2020. [8] YACINE K, ALISTAIR F, YANG Q P. A Bayesian conformity and risk assessment adapted to a form error model[J]. Measurement: Sensors, 2021, 18: 100330. [9] 凌丹. 威布尔分布模型及其在机械可靠性中的应用研究[D]. 成都: 电子科技大学, 2011. LING D. Research on Weibull distribution and its applications in mechanical reliability engineering[D]. Chengdu: University of Electronic Science and Technology of China, 2011. [10] 李晓雷. 高档数控车床的可靠性设计及试验技术研究[D]. 哈尔滨: 哈尔滨工业大学, 2020. LI X L. The research of reliability design and test technology for high grade CNC lathes[D]. Harbin: Harbin Institute of Technology, 2020. [11] 冯静. 小子样复杂系统可靠性信息融合方法与应用研究[D]. 长沙: 国防科学技术大学, 2004. FENG J. Research on methods and applications of reliability information fusion for complex system with small sample test[D]. Changsha: National University of Defense Technology, 2004. [12] 齐乐. 基于MCMC方法的参数生存回归模型的贝叶斯估计[D]. 大连: 大连理工大学, 2018. QI L. Bayesian estimation of parameter survival regression models based on MCMC method[D]. Dalian: Dalian University of Technology, 2018. -
点击查看大图
图(10) / 表(3)
计量
- 文章访问数: 3997
- HTML全文浏览量: 1050
- PDF下载量: 76
- 被引次数: 0
