基于多層自助最大熵法的可靠性評估

夏新濤，葉亮，李云飛，常振

（河南科技大學機電工程學院，河南洛陽471003）

基于多層自助最大熵法的可靠性評估

夏新濤，葉亮，李云飛，常振

（河南科技大學機電工程學院，河南洛陽471003）

在小樣本無任何先驗信息的條件下，提出多層自助最大熵評估模型分析機械產品壽命可靠性。用自助法對當前無失效數據樣本進行再抽樣，獲得足夠的樣本數據?；谧畲箪胤?，改變抽樣個數得到不同的拉格朗日乘子。再次運用自助法對拉格朗日乘子的小樣本數據進行再抽樣，基于最大熵法獲得拉格朗日乘子的區間估計。對每個拉格朗日乘子的上下限進行排列組合，得到多個概率密度函數和可靠性函數，運用最小不確定性原理得到可靠性函數的區間估計。試驗研究表明，多層自助最大熵評估模型可以有效地解決概率分布已知或未知的小樣本無失效數據的可靠性評估問題。

系統評估與可行性分析；可靠性評估；多層自助最大熵法；乏信息；無失效數據；拉格朗日乘子

0 引言

目前，無失效數據的可靠性評估方法主要有經典統計方法和貝葉斯統計方法。經典統計理論認為，通過實驗獲得的失效數據越多，對產品進行可靠性評估的結果越準確。但是，對于許多高可靠性或危險的試驗來說，很難獲得失效數據，例如高端武器和航天航空飛行器的試驗等。因此，無失效數據的可靠性評估問題日益引起學術界的關注。在機械產品可靠性試驗中，獲得的數據通常有3類，即失效數據、不完全失效數據、無失效數據。目前，對失效數據的研究已經達到比較高的水平。然而對無失效數據的研究還不夠深入，因此這也是長期以來國內外研究的熱點與難點。

關于無失效數據的可靠性分析，很多研究工作已經取得了相應的成果。如傅惠民等［1］提出一種Weibull分布定時無失效數據可靠性分析方法，在形狀參數下限已知的情況下，給出了可靠度和使用壽命的單側置信下限，從而能夠根據定時無失效數據對產品進行高置信水平的可靠性評定；王憑慧等［2］在高置信水平下提出了一種衛星推力器無失效數據的可靠性分析方法；李寶盛等［3］求出了工作壽命可靠度置信區間的解析公式，計算和分析了不同試驗樣本工作壽命可靠度的置信限，得出了對可靠性評估有指導意義的結論；蒲星［4］運用Bayes方法對無失效數據問題進行了有效的可靠性研究，分析了無失效數據產生的原因，并提出了運用黃金分割系數求失效概率的方法。但召江等［5］在形狀參數先驗分布分別為均勻分布與擬合出的概率分布時，利用無失效試驗數據得出失效率和形狀參數的Bayes估計，進而計算出Weibull分布的特征壽命；Jiang等［6］提出了運用修整極大似然方法估計威布爾分布的參數。在國外，Nam等［7］以深溝球軸承為例，研究了軸承在高溫條件下的壽命可靠性；Bailey［8］基于無失效數據，建立了二項概率分布的預測模型；Kwon［9-10］基于Bayes方法，分別針對壽命服從Weibull分布和對數正態分布的產品，提出了無失效數據樣本的可靠性驗證方法。

對于未知概率分布的小樣本數據的可靠性評估問題，現有的可靠性評估方法難以解決，到目前為止還是一個重要的科學技術難題。鑒于此，夏新濤等［11］以灰色系統的基本原理為依據，構建出特殊的灰相似空間，允許概率分布未知和數據個數很少（n≥4）。賈波等［12］在仿真實現過程中，發現最大熵方法存在溢出的問題，通過變量變換法成功解決了此問題。夏新濤等［13］將灰自助原理融入泊松過程，提出灰自助泊松方法，以預測滾動軸承振動性能可靠性的變異過程?；谧畲箪胤ǎ蓚ィ?4］提出了火工品類關鍵單元可靠性驗證的LQB型方案。在國外，Yeom等［15］提出改進的最大熵抽樣方法，從而使大多數抽樣點落在一個用復雜的非線性函數表示的可行域內；Young等［16］基于最大熵原理，求解第一結尾時間序列的可靠性邊界。

根據測量理論和統計理論，任何參數估計都伴隨不確定性，可靠性評估中的不確定性可以用可靠性估計區間來表示。現有可靠性理論把拉格朗日乘子當作常數，而Bayes理論認為統計學上的參數可當作變量，可以對參數進行區間估計，因而也能對可靠性函數做出區間估計。因此本文對于只有無失效數據而沒有概率分布任何先驗信息條件下機械產品的可靠性評估問題，提出了多層自助最大熵評估模型來分析機械產品壽命可靠性。具體過程如下:首先，運用自助法對小數據樣本進行等概率可放回再抽樣，獲得大量樣本數據。改變每次抽樣個數，基于最大熵法，可得到多個不同的拉格朗日乘子。然后，再次利用自助法對所有的拉格朗日乘子進行再抽樣，多次使用最大熵法獲得每個拉格朗日乘子的區間估計。最后，對各個拉格朗日乘子的上下限進行排列組合，得到多個概率密度函數和可靠性函數，根據最小不確定性原理獲得可靠性函數的區間估計。通過試驗證明運用該方法能夠對可靠性函數做出合理的區間估計，分析結果真實可信。

1 多層自助最大熵法原理

1.1 自助法原理

假設通過試驗獲得一組無失效數據序列X，用向量表示為

式中:xl為第l個無失效數據；n為無失效數據的個數。

從無失效數據樣本X中等概率可放回地抽樣，每次抽取t個數據，連續重復抽取B次，可以得到自助樣本Xt為

式中:t為每次抽取原始樣本數據的個數；k為生成自助樣本的數據序號；xt（k）為每次抽取t個數據生成自助樣本的第k個數據。

1.2 基于最大熵方法求解概率密度函數

最大熵方法能夠對未知的概率分布做出主觀偏見為最小的最佳估計。在求解過程中，引入拉格朗日乘子，從而把概率分布問題轉化為拉格朗日乘子的求解問題。

為了敘述方便，將離散的無失效數據序列X連續化，定義最大熵的表達式為

式中:f（x）為連續化后的數據序列的概率密度函數；Inf（x）為概率密度函數的對數；S為積分區間。

通過調整f（x）可以使熵達到最大值，同時采取拉格朗日乘子法求解。設為拉格朗日函數，可得

式中:v為原點矩階數，常用v=5；mλ為第λ階原點矩；xλ為求解第λ階原點矩時f（x）的系數；cλ為第λ+1個拉格朗日乘子，λ=0，1，…，v.

其他m個拉格朗日乘子應滿足

1.3 積分區間的映射

為了使求解收斂，將無失效數據按遞增的順序排列并分成j組，畫出直方圖，同時可得到組中值zq和頻數pq，q=2，3，…，j+1.然后將直方圖擴展成j+ 2組，并令p1=pj+2.將原始數據區間S映射到區間［-e，e］中，具體如下:

令

式中:a、b為映射參數；w為所要變換的自變量；x∈［-e，e］，

由dx=dw/a可得

因此，最大熵概率分布密度函數由（6）式變換為

1.4 求解概率密度函數的上下限

現有可靠性理論把拉格朗日乘子當做常數，不能獲得概率密度函數的上下限。因此，本文根據Bayes原理，把拉格朗日乘子當做變量，運用多層自助最大熵方法求出概率密度函數的上下限。

改變抽樣個數k，基于自助最大熵法可得到不同的拉格朗日乘子Cλ和映射參數a、b.

式中:cλ（t）為每次抽取t個數據時的第λ+1個拉格朗日乘子系數；a（t）、b（t）為每次抽取t個數據時的區間映射參數。

1.4.1 各個拉格朗日乘子的區間估計

將各個拉格朗日乘子作為源信息樣本，再進行等概率可放回地抽樣?；谧畲箪胤椒ǎ汕蠼獾酶鱾€拉格朗日乘子的區間估計。

假設顯著性水平為α，α∈［0，1］，則置信水平P為

各個拉格朗日乘子cλ在置信水平P下的下邊界值設為cλL，則有

式中:cλ0為積分變量初始值。

因此，參數的估計區間為

1.4.2 映射參數的點估計

1.4.3 概率密度的上下界函數求解

分別將拉格朗日乘子的上限或下限值，映射參數的點估計值代入（13）式中，由排列組合原理，可得2v+1條概率密度函數曲線的表達式為

式中:i為概率密度函數的序號；c0i為第i個概率密度函數的第1個拉格朗日系數，c0i=c0L或c0i=c0U；cλi為第i個概率密度函數的第λ+1個拉格朗日系數，cλi=cλL或cλi=cλU，cλU、cλL分別為拉格朗日乘子cλ的上下限。

定義概率密度估計真值函數f0（x）為

式中:f0（x）為每次抽取原始樣本數據個數n時的概率密度函數；c00、cλ0分別為每次抽取n個數據時的第1個和第λ+1個拉格朗日系數；a0、b0為每次抽取n個數據時的區間映射參數。

根據最小不確定性原理，可以從概率密度函數曲線中得到距離概率密度估計真值函數f0（x）最近的上下兩條曲線，即最大熵概率密度分布的上下界函數fU（x），fL（x）.

1.5 可靠性的上下界函數求解

1.5.1 可靠性估計真值函數

將每次抽取原始無失效數據樣本個數n時計算所得的可靠性函數作為其估計真值函數R0（x），具體計算如下:對最大熵概率密度估計真值函數f0（x）積分，得到最大熵概率分布估計真值函數F0（x）為

式中:S0為積分區間下限值。

設可靠性系數為rc，由個體無失效數據個數n可得總體的失效概率估計真值函數P0（x）為

總體可靠性估計真值函數也即最大熵可靠性估計真值函數R0（x）為

1.5.2 可靠性函數的區間估計

將最大熵概率密度分布的上下界函數fU（x）、fL（x）分別代入式（26）、（27）、（28）中，可求解出可靠性函數的上下界曲線RU（x）、RL（x）為

給定置信水平P，根據可靠性函數的真值估計曲線R0（x），得到估計真值x0.將其代入到可靠性函數的上下限中，計算出可靠性函數的區間估計為

2 試驗研究

案例1 根據文獻［17］，認為軸承的壽命概率分布為威布爾分布。對20套軸承進行壽命試驗，獲得的無失效數據如下:

X=（422，422，539，539，539，539，602，602，

770，770，770，770，847，847，847，

847，924，924，924，924）.

用自助法每次抽取20個數據，共抽取30 000次，所得數據如圖1所示。

圖1 自助法獲取軸承樣本數據Fig.1 Sample data of bearings obtained by bootstrap method

基于最大熵法計算可得:映射參數a0=0.016 8，b0=-11.974 0；拉格朗日乘子（c00，c10，c20，c30，c40，c50）=（-4.601 9，0.351 5，-1.135 0，-0.022 6，-0.025 2，-0.003 8）.

由（25）式可計算出概率密度估計真值函數f0（x），將f0（x）代入（26）式、（27）式、（28）式中，可得可靠性估計真值函數R0（x），如圖2所示。由圖2可得:R0（677）=0.994 5，與文獻［17］中運用Bayes方法計算出的R（677）=0.9762相差0.0183.這表明自助最大熵方法的可靠性估計真值結果與用Bayes方法的可靠性估計真值結果相差很小。可見自助最大熵方法獲得的可靠性估計真值函數的擬合效果跟Bayes方法的估計結果基本一致。但是自助最大熵方法對概率分布沒有要求，而且能解決小樣本數據的可靠性評估問題。因此，用自助最大熵方法獲得的可靠性估計真值函數評估機械產品的可靠性是有效且可行的。

改變每次抽樣個數t，共抽取30 000次，可得拉格朗日乘子如表1所示。

圖2 可靠性估計真值函數曲線Fig.2 Estimated true value function curve of reliability

表1 改變抽樣個數所得各個拉格朗日乘子Tab.1 Lagrange multipliers obtained by changing the number of samples

設置信水平P=90%，即α=0.10，將所得的17個c0值作為樣本進行自助再抽樣可得:置信區間的下邊界值c0L=-5.156 4；置信區間的上邊界值c0U=-4.869 1.

因此，參數c0的估計區間:［c0L，c0U］=［-5.156 4，-4.8691］.

同樣算得參數c1的估計區間:［c1L，c1U］=［0.4367，0.634 0］；參數c2的估計區間:［c2L，c2U］=［-1.004 8，-0.842 8］；參數c3的估計區間:［c3L，c3U］=［-0.015 7，0.009 9］；參數c4的估計區間:［c4L，c4U］=［-0.044 1，-0.028 4］；參數c5的估計區間:［c5L，c5U］=［-0.007 6，-0.002 3］.

改變每次抽樣個數，共抽取30 000次，可得映射參數如表2所示。

經過計算，參數a、b的估計值分別為

分別將拉格朗日乘子的上限或下限值，映射參數的點估計值代入（24）式中，由排列組合原理和最大熵方法，可得64個概率密度函數。

表2 改變抽樣個數所得映射參數Tab.2 Mapping parameters obtained by changing the number of samples

經過計算，由最小不確定性原理可得:當c0= c0L=-5.156 4，c1=c1U=0.634 0，c2=c2U=-0.842 8，c3=c3U=0.009 9，c4=c4U=-0.028 4，c5=c5U= -0.002 3時，可得到可靠性函數的上限曲線RU（x）；當c0=c0U=-4.869 1，c1=c1L=0.436 7，c2=c2L= -1.004 8，c3=c3L=-0.015 7，c4=c4L=-0.044 1，c5=c5L=-0.007 6時，可得到可靠性函數的下限曲線RL（x）.如圖3所示。

圖3 可靠性函數曲線Fig.3 Curves of reliability functions

由圖3可知可靠性函數取值隨著自變量時間的增大而逐漸減小，這是符合工程實際的。而且在x= 700 h之前，可靠性函數上下限曲線圖均與可靠性估計真值函數曲線圖基本完全重合。在x=700 h之后，隨著時間的增加，可靠性函數上下限曲線均越來越偏離可靠性估計真值函數曲線。這是因為在當前時間段內，對可靠性函數估計的不確定度較小。隨著時間變量的增大，可靠性函數的估計難度逐漸增大，對可靠性函數估計的不確定度也逐漸增大，即上下限值越來越偏離估計真值。

設可靠性系數rc=0.1，置信水平P=90%，即顯著性水平α=0.1，計算出估計真值x0:x0=826 h.可靠性函數的區間估計為

案例2 根據文獻［18］，認為導彈壽命的概率分布為指數分布。對19個導彈進行壽命試驗，獲得的無失效數據如下:

用自助法每次抽取19個數據，共抽取30 000次，所得數據如圖4所示。

圖4 自助法獲取導彈樣本數據Fig.4 Sample data of missiles obtained by bootstrap method

基于最大熵法計算可得:映射參數a0=3.575 6，b0=-6.774 8.拉格朗日乘子（c00，c10，c20，c30，c40，c50）=（0.730 4，-0.150 8，-0.939 5，0.009 6，-0.103 9，-0.013 1）.

由（25）式可計算出概率密度估計真值函數f0（x），將f0（x）代入（26）式、（27）式、（28）式中，可得可靠性估計真值函數R0（x），如圖5所示。

文獻［18］中運用Bayes方法計算出R（3）= 0.916 4，運用多層自助最大熵方法計算出的R0（3）=0.898 0，二者相差0.018 4.這表明多層自助最大熵方法的可靠性估計真值結果與用Bayes方法的可靠性估計真值結果誤差很小?？梢娮灾畲箪胤椒ǐ@得的可靠性估計真值函數的擬合效果跟Bayes方法的估計結果基本一致。但是自助最大熵方法對概率分布沒有要求，而且能解決小樣本數據（n=19）的可靠性評估問題。因此，用多層自助最大熵方法評估機械產品的可靠性是有效且可行的。

改變每次抽樣個數t，共抽取30 000次，所得拉格朗日乘子如表3所示。

設置信水平P=90%，即α=0.10，將所得的16個c0值作為樣本進行自助再抽樣可得:置信區間的下邊界值c0L=0.359 4，置信區間的上邊界值c0U=0.592 4.

圖5 可靠性估計真值函數曲線Fig.5 Estimated true value function curve of reliability

表3 改變抽樣個數所得各個拉格朗日乘子Tab.3 Lagrange multipliers obtained by changing the number of samples

因此，參數c0的估計區間:［c0L，c0U］=［0.359 4，0.592 4］.

同樣算得參數c1的估計區間:［c1L，c1U］=［0.1934，0.363 0］；參數c2的估計區間:［c2L，c2U］=［-1.213 1，-0.998 8］；參數c3的估計區間:［c3L，c3U］=［-0.016 3，0.010 3］；參數c4的估計區間:［c4L，c4U］=［-0.043 2，0.011 1］；參數c5的估計區間:［c5L，c5U］=［-0.019 3，0.007 6］.

改變每次抽樣個數，共抽取30 000次，可得映射參數如表4所示。

經過計算，參數a、b的估計值分別為

分別將拉格朗日乘子的上限或下限值，映射參數的點估計值代入（24）式中，由排列組合原理和最大熵方法，可得64個概率密度函數。

經過計算，由最小不確定性原理可得:當c0= c0L=0.359 4，c1=c1U=0.363 1，c2=c2U=-0.998 8，c3=c3U=0.010 3，c4=c4U=0.011 1，c5=c5L= -0.019 3時，可得到可靠性函數的上限曲線RU（x）；

當c0=c0U=0.592 4，c1=c1L=0.193 4，c2=c2L= -1.2131，c3=c3L=-0.016 3，c4=c4L=-0.0432，c5=c5L=-0.019 3時，可得到可靠性函數的下限曲線RL（x）.如圖6所示。

表4 改變抽樣個數所得映射參數Tab.4 Mapping parameters obtained by changing the number of samples

圖6 可靠性函數曲線Fig.6 Curves of reliability functions

由圖6可知:可靠性函數取值隨著自變量時間的增大而逐漸減小，這是符合工程實際的。而且在x=1.9 a之前，可靠性函數上下限曲線圖均與可靠性估計真值函數曲線圖基本完全重合。在x=1.9 a之后，隨著時間的增加，可靠性函數上下限曲線均越來越偏離可靠性估計真值函數曲線。這是因為在當前時間段內，對可靠性函數估計的不確定度較小。隨著時間變量的增大，可靠性函數的估計難度逐漸增大，對可靠性函數估計的不確定度也逐漸增大，即上下限越來越偏離估計真值。

設可靠性系數rc=0.1，置信水平P=90%，即顯著性水平為α=0.1，計算出估計真值x0:x0= 2.375 a.可靠性函數的區間估計為

案例3 某電子產品壽命的概率分布未知，對其進行模擬實驗，得到的無失效數據如下:

用自助法每次抽取10個數據，共抽取30 000次，所得數據如圖7所示。

圖7 自助法獲取樣本數據Fig.7 Sample data obtained by bootstrap method

基于最大熵法可得:映射參數a0=0.150 4，b0=-5.412 8.拉格朗日乘子（c00，c10，c20，c30，c40，c50）=（-2.473 8，-0.321 1，-0.994 8，-0.023 4，-0.031 2，0.004 7）.

由（25）式可計算出概率密度估計真值函數f0（x），將f0（x）代入（26）式、（27）式、（28）式中，可得可靠性估計真值函數R0（x），如圖8所示。

改變每次抽樣個數t，共抽取30 000次，可得拉格朗日乘子如表5所示。

設置信水平P=90%，即α=0.1，將所得的7個c0值作為樣本進行自助再抽樣可得:置信區間的下邊界值c0L=-2.772 4，置信區間的上邊界值c0U= -2.561 2.

因此，參數c0的估計區間:［c0L，c0U］=［-2.772 4，-2.561 2］.

同樣算得參數c1的估計區間:［c1L，c1U］=［-0.224 7，-0.095 6］；參數c2的估計區間:［c2L，c2U］=［-0.985 3，-0.785 4］；參數c3的估計區間:［c3L，c3U］=［0.002 0，0.032 5］；參數c4的估計區間:［c4L，c4U］=［-0.058 1，0.022 5］；參數c5的估計區間:［c5L，c5U］=［-0.002 6，0.003 1］.

表5 改變抽樣個數所得各個拉格朗日乘子Tab.5 Lagrange multipliers obtained by changing the number of samples

圖8 可靠性估計真值函數曲線Fig.8 Estimated true value function curve of reliability

改變每次抽樣個數，共抽取30 000次，可得映射參數如表6所示。

表6 改變抽樣個數所得映射參數Tab.6 Mapping parameters obtained by changing the number of samples

經過計算，參數a、b的估計值分別為

分別將拉格朗日乘子的上限或下限值，映射參數的點估計值代入（24）式中，由排列組合原理和最大熵方法，可得64個概率密度函數。

經過計算，由最小不確定性原理可得:當c0=c0L= -2.772 4，c1=c1U=-0.095 6，c2=c2U=-0.785 4，c3=c3U=0.032 5，c4=c4U=0.022 5，c5=c5L= -0.002 6時，可得到可靠性函數的上限曲線RU（x）；當c0=c0U=-2.561 2，c1=c1L=-0.224 7，c2=c2L=-0.985 3，c3=c3U=0.002 0，c4=c4L= -0.058 1，c5=c5U=0.003 1時，可得到可靠性函數的下限曲線RL（x）.如圖9所示。

圖9 可靠性函數曲線Fig.9 Curves of reliability functions

由圖9可知:可靠性函數取值隨著自變量時間的增大而逐漸減小，這是符合工程實際的。而且在x=33 h之前，可靠性函數上下限曲線均與可靠性估計真值函數曲線基本完全重合。在x=33 h之后，隨著時間的增加，可靠性函數上下限曲線均越來越偏離可靠性估計真值函數曲線。這是因為在當前時間段內，對可靠性函數估計的不確定度較小。隨著時間變量的增大，可靠性函數的估計難度逐漸增大，對可靠性函數估計的不確定度也逐漸增大，即上下限值越來越偏離估計真值。

取可靠性系數rc=0.1，置信水平P=90%，即顯著性水平為α=0.1，計算出估計真值x0:x0= 47.5 h.

可靠性函數的區間估計:

案例4 根據文獻［19］，機電設備壽命T服從對數正態分布LN（μ，σ2）時。4臺電機分別工作到3 782.2 h、4 212.6 h、4 219.2 h和4 476.2 h均未失效。由于試驗數據收集的困難性，希望利用上述數據對該型號電機的可靠性進行評估。

用自助法每次抽取4個數據，共抽取10 000次，所得數據如圖10所示。

圖10 自助法獲取樣本數據Fig.10 Sample data of bearings obtained by bootstrap method

基于最大熵法計算可得:映射參數a0=0.003 0，b0=-12.476 7.拉格朗日乘子（c00，c10，c20，c30，c40，c50）=（-7.272 5，1.176 6，0.038 2，-0.637 4，-0.053 8，0.071 0）.

由（25）式可計算出概率密度估計真值函數f0（x），將f0（x）代入（26）式、（27）式、（28）式中，可得可靠性估計真值函數R0（x），如圖11所示。

圖11 可靠性估計真值函數曲線Fig.11 Estimated true value function curve of reliability

運用多層自助最大熵方法計算出的R0（4 000）= 0.955 8，R0（4 500）=0.925 0.文獻［19］中運用Bayes方法計算出當λ=0.2時，μ的估計值為10.721 8，σ的估計值為1.767 3，從而可得R（4 000）=0.915 2；R（4 500）=0.904 4.可見兩種方法計算結果分別相差0.040 6和0.020 6，多層自助最大熵方法的可靠性估計真值結果與用Bayes方法的可靠性估計真值結果誤差很小。因此自助最大熵方法獲得的可靠性估計真值函數的擬合效果跟Bayes方法的估計結果基本一致。自助最大熵方法可以自動識別樣本數據內部規律，從而計算出小樣本數據（n=4）的概率分布函數。在數據處理過程中，并未利用已知的概率分布先驗信息。因此，用多層自助最大熵方法評估機械產品的可靠性是有效且可行的。

值得注意的是，由于樣本數據較少，無法通過改變抽樣個數獲得拉格朗日乘子的樣本數據，因此無法進行可靠性區間估計。進行可靠性區間估計需要原樣本數據個數n≥7.

案例5 針對案例4，進行模擬試驗，得到的無失效數據如下:

用自助法每次抽取7個數據，共抽取10 000次，所得數據如圖12所示。

圖12 自助法獲取樣本數據Fig.12 Sample data obtained by bootstrap method

基于最大熵法可得:映射參數a0=0.009 2，b0=-38.354 0.拉格朗日乘子（c00，c10，c20，c30，c40，c50）=（-5.653 8，1.294 6，-1.345 5，0.412 5，0.045 0，0.092 4）.

由（25）式可計算出概率密度估計真值函數f0（x），將f0（x）代入（26）式、（27）式、（28）式中，可得可靠性估計真值函數R0（x），如圖13所示。

改變每次抽樣個數t，共抽取30 000次，可得拉格朗日乘子如表7所示。

設置信水平P=90%，即α=0.1，將所得的7個c0值作為樣本進行自助再抽樣可得:置信區間的下邊界值c0L=-5.804 9，置信區間的上邊界值c0U= -5.703 2.

因此，參數c0的估計區間:［c0L，c0U］=［-5.804 9，-5.703 2］.

同樣算得參數c1的估計區間:［c1L，c1U］=［1.248 6，1.624 0］；參數c2的估計區間:［c2L，c2U］=［-1.582 4，-0.704 4］；參數c3的估計區間:［c3L，c3U］=［0.002 0，0.032 5］；參數c4的估計區間:［c4L，c4U］=［-0.098 7，0.082 6］；參數c5的估計區間:［c5L，c5U］=［-0.065 2，0.091 5］.

改變每次抽樣個數，共抽取30 000次，可得映射參數如表8所示。

圖13 可靠性估計真值函數曲線Fig.13 Estimated true value function curve of reliability

表7 改變抽樣個數所得各個拉格朗日乘子Tab.7 Lagrange multipliers obtained by changing the number of samples

表8 改變抽樣個數所得映射參數Tab.8 Mapping parameters obtained by changing thenumber of samples

經過計算，參數a、b的估計值分別為

分別將拉格朗日乘子的上限或下限值，映射參數的點估計值代入（24）式中，由排列組合原理和最大熵方法，可得64個概率密度函數。

經過計算，由最小不確定性原理可得:當c0= c0L=-5.804 9，c1=c1U=1.624 0，c2=c2U=-0.704 4，c3=c3U=0.032 5，c4=c4U=0.082 6，c5=c5U=0.091 5時，可得到可靠性函數的上限曲線RU（x）；當c0= c0U=-5.703 2，c1=c1L=1.248 6，c2=c2U=-0.704 4，c3=c3L=0.002 0，c4=c4L=-0.098 7，c5=c5L= -0.065 2時，可得到可靠性函數的下限曲線RL（x）.如圖14所示。

圖14 可靠性函數曲線Fig.14 Curves of reliability functions

由圖14可知:可靠性函數取值隨著自變量時間的增大而逐漸減小，這是符合工程實際的。而且在x=4 230 h之前，可靠性函數上下限曲線均與可靠性估計真值函數曲線基本完全重合。在x=4 230 h之后，隨著時間的增加，可靠性函數上下限曲線均越來越偏離可靠性估計真值函數曲線。這是因為在當前時間段內，對可靠性函數估計的不確定度較小。隨著時間變量的增大，可靠性函數的估計難度逐漸增大，對可靠性函數估計的不確定度也逐漸增大，即上下限值越來越偏離估計真值。

取可靠性系數rc=0.1，置信水平P=90%，即顯著性水平為α=0.1，計算出估計真值x0:x0= 4 425.4 h.可靠性函數的區間估計:

3 討論

1）對于概率分布已知的產品壽命可靠性評估問題，自助最大熵方法的評估誤差很小。案例1用多層自助最大熵法計算出的R0（577）=0.995 2與文獻［17］中運用Bayes方法計算出的R（577）= 0.984 5相差0.010 7；R0（677）=0.994 5與文獻［17］中運用Bayes方法計算出的R（677）=0.976 2相差0.018 3.文獻［18］中運用Bayes方法計算出R（3）=0.916 4，案例2運用多層自助最大熵法計算出的R0（3）=0.898 03，二者相差0.018 37.案例4運用多層自助最大熵方法計算出的R0（4 000）= 0.955 8，R0（4 500）=0.925 0.文獻［19］中運用Bayes方法計算出R（4 000）=0.915 2，R（4 500）= 0.904 4，兩種方法計算結果分別相差0.040 6和0.020 6.試驗案例1、案例2、案例4表明運用自助最大熵方法得到的可靠性評估結果與運用Bayes方法得到的評估結果基本相同。這是因為最大熵方法能夠對未知的概率分布做出主觀偏見為最小的最佳估計，可以自動識別樣本數據內部規律，從而計算出樣本數據個數極少（n≥4）條件下的概率分布函數。多層自助最大熵法在數據處理過程中，并未利用已知概率分布（Weibull分布、指數分布、對數正態分布）這一信息。

2）對于概率分布未知的產品壽命可靠性評估問題，古典統計理論無法解決。而案例3運用多層自助最大熵方法，計算得到R0（56.97）=0.898 9；案例5運用多層自助最大熵方法計算得到R0（4 500）= 0.897 5.試驗案例表明多層自助最大熵方法不但適用于概率分布已知的產品壽命可靠性評估問題，也適用于概率分布未知的產品壽命可靠性評估問題。

3）對于可靠性函數的區間估計問題，現有方法無法解決。本文中的多層自助最大熵方法把拉格朗日乘子當作變量，對可靠性函數進行區間估計，得出結論:可靠性函數取值隨著自變量時間的增大而逐漸減小，在一定范圍內，可靠性函數上下限曲線圖與可靠性估計真值函數曲線圖基本完全重合。超出該范圍，可靠性函數上下限曲線均越來越偏離可靠性估計真值函數曲線。

4）實際應用中，一般在達到產品無失效數據樣本最大值的一半之前，可靠性函數上下限曲線圖與可靠性估計真值函數曲線圖基本完全重合。如案例1中x=700 h＞1/2×924=462 h，案例2中x= 1.9 a＞1/2×3=1.5 a，案例3中x=33 h＞1/2× 56.97=28.485 h.

4 結論

將拉格朗日乘子當做變量，提出多層自助最大熵方法，以解決機械產品的乏信息可靠性評估問題。該方法對概率分布沒有要求，可以有效地解決乏信息條件下無失效數據的可靠性評估問題與可靠性函數的區間估計問題，試驗證明該方法是對現有可靠性評估方法的有益補充。

（References）

［1］ 傅惠民，張勇波.Weibull分布定時無失效數據可靠性分析方法［J］.航空動力學報，2010，25（12）:2807-2810. FU Hui-min，ZHANG Yong-bo.Method of reliability analysis for time truncated zero-failure data based on Weibull distribution［J］. Journal of Aerospace Power，2010，25（12）:2807-2810.（in Chinese）

［2］ 王憑慧，范本堯，傅惠民.衛星推力器可靠性評估和壽命預測［J］.航空動力學報，2004，19（6）:745-748. WANG Ping-hui，FAN Ben-yao，FU Hui-min.Reliability assessment and life prediction for satellite engine［J］.Journal of Aerospace Power，2004，19（6）:745-748.（in Chinese）

［3］ 李寶盛，何洪慶.“小子樣、零失效”情況下壽命可靠度的置信分析方法［J］.兵工學報，2001，22（2）:234-237. LI Bao-sheng，HE Hong-qing.Method of confidence analysis for life reliability on small and zero failure samples［J］.Acta Armamentarii，2001，22（2）:234-237.（in Chinese）

［4］ 蒲星.一種無失效數據可靠性分析的方法［J］.電子質量，2014（7）:29-31. PU Xing.A method for reliability analysis of zero-failure data［J］. Journal of Electronic Quality，2014（7）:29-31.（in Chinese）

［5］ 但召江，樓洪梁，李興林，等.無失效數據下滾動軸承的可靠性估計［J］.軸承，2013（9）:22-24. DAN Shao-jiang，LOU Hong-liang，LI Xing-lin，et al.Reliability estimation for rolling bearings based on zero-failure data［J］.Journal of Bearing，2013（9）:22-24.（in Chinese）

［6］ Jiang P，Lim J H，Zuo M J，et al.Reliability estimation in a Weibull lifetime distribution with zero-failure field data［J］. Quality and Reliability Engineering International，2010，26（7）: 691-701.

［7］ Nam J S，Kim H E，Kim K U.A new accelerated zero-failure test model for rolling bearings under elevated temperature conditions［J］.Journal of Mechanical Science and Technology，2013，27（6）:1801-1807.

［8］ Bailey R T.Estimation from zero-failure data［J］.Risk Analysis，1997，17（3）:375-380.

［9］ Kwon Y.Design of Bayesian zero-failure reliability demonstration test for products with Weibull lifetime distribution［J］.Journal of Applied Reliability，2014，14（4）:220-224.

［10］ Kwon Y.The effect of scale parameter in designing reliability demonstration test for lognormal lifetime distribution［J］.Journal of Applied Reliability，2014，14（1）:53-57.

［11］ 夏新濤，陳曉陽，張永振，等.兩個數據序列的灰關系及其工程應用的兩個例子［J］.兵工學報，2008，29（1）:73-77. XIA Xin-tao，CHEN Xiao-yang，ZHANG Yong-zhen，et al.E-valuation for stability of manufacturing process based on grey relation［J］.Acta Armamentarii，2008，29（1）:73-77.（in Chinese）

［12］ 賈波，劉福，雷正偉.基于最大熵方法對測量數據估計的改進方法研究［J］.價值工程，2010（28）:228-230. JIA Bo，LIU Fu，LEI Zheng-wei.Study on evaluation of measurement result based on improved maximum entropy method［J］. Journal of Value Engineering，2010（28）:228-230.（in Chinese）

［13］ 夏新濤，孟艷艷，邱明.用灰自助泊松方法預測滾動軸承振動性能可靠性的變異過程［J］.機械工程學報，2015，51（9）: 97-103. XIA Xin-tao，MENG Yan-yan，QIU Ming.Forecasting for variation process of reliability of rolling bearing vibration performance using grey bootstrap Poisson method［J］.Journal of Mechanical Engineering，2015，51（9）:97-103.（in Chinese）

［14］ 胡松偉.最大熵法在發射裝置可靠性評估中的應用［J］.四川兵工學報，2014，35（5）:31-33. HU Song-wei.Maximum entropy method application on missile launching system reliability evaluation［J］.Journal of Sichuan Ordnance，2014，35（5）:31-33.（in Chinese）

［15］ Yeom D H，Bae S G.A constraint-based maximum entropy sampling method for Kriging models in fuel cell applications［J］. Quality and Reliability Engineering International，2013，29（6）: 819-828.

［16］ Young H S，Byung M K.Reliability bound based on the maximum entropy principle with respect to the first truncated moment［J］.Journal of Mechanical Science and Technology，2010，24（9）:1891-1900.

［17］ 茆詩松，夏劍鋒，管文琪.軸承壽命試驗中無失效數據的處理［J］.應用概率統計，1993，9（3）:326-331. MOU Shi-song，XIA Jian-feng，GUAN Wen-qi.Processing of zero-failure data in bearing life test［J］.Chinese Journal of Applied Probability and Statistics，1993，9（3）:326-331.（in Chinese）

［18］ 孫亮，徐廷學，王冬梅.某型導彈無失效數據的處理方法［J］.戰術導彈技術，2004（3）:29-32. SUN Liang，XU Ting-xue，WANG Dong-mei.Process method for zero-failure data of certain type missile［J］.Tactical Missile Technology，2004（3）:29-32.（in Chinese）

［19］ 孫建中.機電設備無失效數據的可靠性評估［J］.艦船科學技術，2006，28（1）:50-53. SUN Jian-zhong.Reliability estimation of zero-failure data for mechanical and electrical products［J］.Ship Science and Technology，2006，28（1）:50-53.（in Chinese）

Reliability Evaluation Based on Hierarchical Bootstrap Maximum Entropy Method

XIA Xin-tao，YE Liang，LI Yun-fei，CHANG Zhen
（School of Mechatronical Engineering，Henan University of Science and Technology，Luoyang 471003，Henan，China）

A hierarchical bootstrap maximum entropy evaluation model is proposed to analyze the life reliability of mechanical products under the condition of small samples without any prior information.Adequate sample data is obtained by using bootstrap method to re-sample the current zero-failure data samples.Based on maximum entropy method，the different Lagrange multipliers can be obtained by changing the number of samples.In order to get the interval estimation values of Lagrange multipliers，the bootstrap method is used again to re-sample the small sample data of Lagrange multipliers.The probability density functions and reliability functions are achieved by carrying on permutation and combination for the upper and lower limit values of each Lagrange multiplier，so the interval estimation values of reliability functions can be gained using minimum uncertainty principle.Experimental investigation shows that the hierarchical bootstrap maximum entropy evaluation model can effectively solve the reliability evaluation problem for zero-failure data of small samples with known or unknown probability distributions.

system assessment and feasibility analysis；reliability evaluation；hierarchical bootstrap maximum entropy method；poor information；zero-failure data；Lagrange multipliers

TB114.3

A

1000-1093（2016）07-1317-13

10.3969/j.issn.1000-1093.2016.07.022

2016-01-05

國家自然科學基金項目（51475144、51075123）

夏新濤（1957—），男，教授，博士生導師，博士。E-mail:xiaxt1957@163.com；

葉亮（1990—），男，碩士研究生。E-mail:172682823@qq.com