基于马氏链蒙特卡洛方法

1、基于马氏链蒙特卡洛方法的数控系统可靠性评估李斌全，戴怡 （天津职业技术师范大学，机械工程学院，天津 ）摘要：数控系统属于高可靠性产品，应用贝叶斯理论，已经成为其可靠性评估的重要方法。由于数控系统寿命试验数据服从Weibull分布，导致后验分布出现数值大、高维复杂的情况，贝叶斯计算依靠的数值积分方法难以实施。根据马氏链蒙特卡洛（MCMC）方法思想，建立后验平稳分布的马尔科夫链，对分布参数的贝叶斯估计进行求解，解决了数值积分问题，保证可靠性评估的有效实施。通过与BUGS软件结果进行比较，表明提高了模型计算的稳健性、有效性及精度。关键词：数控系统；马氏链蒙特卡洛（MCMC）方法；可靠性评估；Weib

2、ull分布；贝叶斯理论中图分类号：TH 17The reliability evaluation of NC system based on Markov chain Monte Carlo methodLi Bin-quan，Dai Yi(School of Mechanical Engineering, Tianjin University of Technology and Education, Tianjin ，China)Abstract: NC system was high-reliability product, so the application of Bayesian t

3、heory had become an important means of its reliability evaluation. The failure data of NC system obeyed Weibull distribution which had complex forms. So the posterior distribution became more complicated and the numerical integration which Bayesian computing depended on was not available. The proble

4、m of numerical integration was solved by Markov chain Monte Carlo (MCMC) method. It ensures the effective implementation of reliability evaluation. Compared with BUGS package, the result of Bayesian estimation shows that it increases the robustness and effectiveness of the calculation. Key words: NC

5、 system；Markov chain Monte Carlo (MCMC) method；reliability evaluation；Weibull distribution；Bayesian theory0 引言可靠性工程中，对产品进行可靠性研究或可靠性特征量预测，最基本的工作是要确定其寿命分布参数。据多位学者的工作，数控系统寿命服从两参数Weibull分布1,2。其分布参数的贝叶斯估计需要对后验分布进行二重积分。根据贝叶斯公式，后验分布为基金项目：国家自然科学基金（）；国家数控机床科技重大专项（2009ZX04014-013）（1）其中，和为Weibull分布两参数，为样本似然函数，

6、为先验分布。则分布参数的贝叶斯估计为 （2） （3）由于数控系统寿命较长且Weibull分布形式复杂，使得被积函数数值巨大。在实际评估中，计算（2）、（3）式的内层积分时就已经出现浮点数溢出情况，导致依靠被积函数值计算的数值积分方法失效随着马氏链蒙特卡洛（MCMC）方法的不断发展，贝叶斯后验分布计算有望得到彻底解决。通过MCMC方法，复杂的后验分布被直观的模拟出来，这就使得参数求解运用的数值积分方法转化成从简单的分布中抽样并推断。基于抽样所得分布参数样本，后验分布的统计计算简单易行，各种可靠性特征量求解的可操作性显著提高。目前贝叶斯计算多用国外的BUGS软件进行3,4，但是该软件未公布其核心算

7、法，这给掌握算法实质以及进行算法改进带来困难。另一方面BUGS软件采用的Gibbs抽样是MCMC方法的特例，所以自主掌握MCMC方法，并应用于可靠性评估有重要意义。本文利用MCMC方法对参数的贝叶斯估计进行求解，并与BUGS软件结果进行比较，表明其适用于数控系统的可靠性评估。1 基于MCMC方法的贝叶斯可靠性1.1截尾样本的Weibull似然函数构建对数控系统进行有替换定时截尾寿命试验，得到，其中，为完全失效数据样本，为截尾数据样本。两参数Weibull分布概率密度函数为，生存函数为，其中为形状参数，为尺度参数。则截尾样本似然函数为 （4）1.2 后验分布形式的确立取参数的先验分布为伽马分布，

8、其中为形状参数，为尺度参数。参数的先验分布为逆伽马分布，其中为形状参数，为尺度参数。则两参数的联合先验分布为 （5）由贝叶斯公式得到后验分布： （6）由于后验分布数值巨大，高维复杂，贝叶斯计算依靠的数值积分方法难以实施。但是可以利用MCMC方法思想解决上述问题。2 后验平稳分布的马尔科夫链建立贝叶斯后验量的计算可归结为关于后验分布的积分5。 MCMC方法思想是通过建立后验平稳分布的马尔科夫链6，并从中抽样，来获取后验分布样本进行统计推断，这就不需要对被积函数进行数值计算，从而保证求解过程的有效实施和稳定性。 令，表示从易于抽样的建议分布产生一个潜在的转移，然后根据概率决定是否转移。如果链在时刻

9、t处于状态x，即，从中抽一个随机数u， （7）则此时，为马氏链一个转移核。在选取后，根据（6）式，概率取为 （8）这时，转移核以为马氏链平稳分布。选取均匀分布分布为建议分布，此时， （9）其中，为转移概率，和是在一步迭代中，从建议分布中抽取的建议值，和是一步迭代的初始值。过程如下:（1）（2）从建议分布中产生候选点（3）令，即让建议值随机游动一个距离（4）根据式（9）计算接受概率（5）以概率接受，否则令通过上述过程，后验分布的马尔科夫链被建立起来，并得到Weibull分布参数和的样本，基于所得样本就可以对分布参数进行各种统计推断。3 仿真分析通过有替换定时截尾寿命试验，获得数控系统寿命的现场样

10、本如表1所示，表1数控系统寿命现场样本（单位：h）完全失效数据176256339455495117612601896236221342133322431853158151120134102438013604251截尾数据11369315232231245625883322102242131021250246614503210310250116011145141522332根据表1数据，用上述方法编程进行计算，得到结果如表2表2MCMC方法计算结果建议分布接受率E(m)E()sd(m)sd()均匀分布0.511.3540990.171.8BUGS是进行贝叶斯推断的软件包，用该软件处理表1数据，得

11、到结果如表3表3应用BUGS软件计算结果nodemeansdstartsamplem1.320.1810015700041114.1100157000通过表2计算结果，表明MCMC方法的有效性与稳定性，与表3进行比较，可以看出其计算精度也有所提高。图1、2分别是m和的全部迭代轨迹和后验分布的密度直方图。（a）参数m（b）参数图1 参数m和的迭代轨迹图（a）参数m（b）参数图2 后验分布的密度直方图4 数控系统可靠性特征量的计算用建议分布为均匀分布时建立的马氏链对数控系统可靠性进行估计，Weibull分布概率密度函数采用传统方法计算MTBF值的区间估计是复杂而且难以实施的。利用MCMC方法可以有

12、效地予以估计。图3是MTBF值的迭代轨迹。图3 MTBF值的迭代轨迹图可以得到，MTBF值的90%置信区间为3666.76，3989.30 。5 结论（1）应用贝叶斯理论进行可靠性评估是一种有前途的方法。由于以往不能很好地解决高维数值积分问题，从而制约了该理论在实际中的应用。MCMC方法解决了高维数值积分问题带来的不便，有利于贝叶斯理论在可靠性评估中的推广。（2）由于后验分布被模拟构造出来，后验分布的各种统计推断都可从简单分布中抽样进行。基于抽样所得分布参数样本，可靠度、失效率等各种可靠性特征量的求解都可以有效实施，不再受数值积分的限制，提高了模型计算的稳定性、可操作性与适应性参考文献1 戴怡

13、，周云飞，陈学东等. 加工中心故障分布规律及其研究方法J. 系统工程与电子技术，2004，26(3)：413-415.2 张海波，贾亚洲，周广文. 数控系统故障间隔时间分布模型的研究J. 哈尔滨工业大学学报，2005，37(2)：198-200.3 林静. 基于MCMC的贝叶斯生存分析理论及其在可靠性评估中的应用D. 南京：南京理工大学，2008.4 林静，韩玉启，朱慧明等. 基于MCMC稳态模拟的Weibull回归模型及其可靠性应用J.系统仿真学报，2006，18(5) ：1161-1163.5 茆诗松，王静龙，濮晓龙. 高等数理统计M.北京：高等教育出版社，1998.6 Geof H. Given