技术领域
[0001] 本
发明是一种针对退化数据，基于贝叶斯理论的可靠性综合评估方法，属于可靠性综合评估技术领域。
背景技术
[0002] 通常，在评估产品的可靠性时，通过
加速寿命试验来获得产品的寿命数据，并利用所获得的寿命数据来对产品的可靠性进行评估。但随着科学技术的发展，长寿命、高可靠性产品正逐渐成为当今的主流，因此需要评估可靠性的对象的特点亦向着长寿命、高可靠性方面转移。由于产品的长寿命、高可靠性特点，在通过试验来评估产品的可靠性时，寿命数据越来越难以获得，这就对评估造成了一定的困难。针对这一困难，加速退化试验应运而生，通过试验时实时获得的大量的产品性能退化数据来对产品的可靠性进行评估。目前，应用比较广泛的加速退化可靠性评估方法是北京航空航天大学的姜同敏，李晓阳等人提出的基于漂移
布朗运动的加速退化试验可靠性评估方法。
[0003] 由于产品自身、试验设备、测试设备等多方面的原因，一次试验的结果往往具有一定的随机性、片面性，并不能展现产品可靠性特征的全貌。从评估
精度的
角度上看，如果用于评估的信息越多、越全面，那么所得到的评估结果就越精确，越接近评估对象的实际情况。为了更加准确全面的评估产品的可靠性，就需要将产品在研制过程中的各试验数据进行综合利用，有效的利用各方面的数据并使用合理的方法对其进行融合，从而得到更加精确的评估结果。因此，可靠性综合评估技术就成为产品可靠性评估中重要的一环。
[0004] 对于综合评估的方法，目前研究的比较多的为贝叶斯方法，在20世纪60年代，国际上就己有人将贝叶斯方法用于可靠性统计分析。贝叶斯可靠性评估技术的优势在于能充分利用大量的先验可靠性信息，然后结合少量样本数据进行可靠性评估，从而达到节约试验经费、缩短产品研制周期的目的。贝叶斯方法的这一特点引起了人们的强烈兴趣，许多学者在此领域开展了相关研究工作并取得了大量的研究成果。到了80年代，己有这方面的专著，系统而详尽地总结了这一方向的工作。近些年来，学者们更注重与可靠性工程实践相结合，在现场试验数据较少的情况下，利用贝叶斯方法综合各种先验信息和小样本试验数据，以期得出更加合理和可信的评估结果。Pate-cornel，Coolen，Lichtenstein和Newman等学者对工程中存在的专家经验知识的收集、
整理和合理利用问题进行了相关研究。Erto，Selby和Shoukri，Elperin和Gertsbakh，Arturo等学者利用工程中常见的先验信息类型，分别对不同寿命分布类型的可靠性评估问题进行了分析。在国内，国防科技大学在这方面作了许多的工作，张金槐教授著有《贝叶斯试验分析方法》一书，2000年张士峰对多源验前信息的融合方法进行了研究，2003年张湘平结合武器装备小子样试验分析与评估的特点，以贝叶斯方法为小子样统计推断与融合理论的研究主线，重点研究了适应于贝叶斯方法工程运用的若干共性问题。并在武器装备的精度和可靠性评估两个方面展开了应用研究，2005年段晓君研究了基于先验融合的贝叶斯递归精度评定方法，以及试验鉴定中信息的量化评估模型，满军，张湘平研究了基于可信度的贝叶斯融合方法在可靠性评估中的应用，
2006年刘晗研究了基于贝叶斯理论的小子样可靠性评定方法，刘飞针对固体火箭
发动机可靠性增长试验中存在的理论与实际问题，综合利用可靠性增长试验的多种信息，系统地提出了可靠性增长试验管理，规划和分析方法，2007年吴建业利用组成导弹系统的各分系统，单元的试验信息，运用贝叶斯方法综合利用各验前信息以及不同阶段的试验数据，并结合少量现场试验数据对反舰导弹飞行可靠性进行了评估，2008年李欣欣通过贝叶斯方法在多阶段可靠性增长试验中，来融合历史阶段数据对现阶段可靠性
水平进行综合评估，冯静研究了基于贝叶斯
模糊逻辑算子的小子样可靠性信息融合方法。西北工业大学的刘松林在
2006年研究了以Kullback信息作为融合准则，
串联系统在具有多源验前信息的情形下的可靠性评估问题，沈政等人在2008年研究了多源验前信息下结合专家经验的融合方法。
[0005] 但目前在可靠性综合评估技术方面的研究主要集中在基于寿命数据的综合评估方面，而对于基于退化数据的综合评估方法，虽然国防科技大学的张永强、赵炤、冯静，哈尔滨工业大学的李常有，华中科技大学的卫军，清华大学的王剑等人在将贝叶斯方法应用于退化数据方面有一定的研究，但并没有一套系统的综合利用多源退化数据、建立统计模型、并对不同应
力条件下的退化数据进行综合评估的具体方法。
[0006] 因此在可靠性评估领域中，随着应用加速退化试验及其退化数据来进行可靠性评估的日益广泛，提出针对加速退化试验及其退化数据的可靠性综合评估方法是目前可靠性评估领域的当务之急。
发明内容
[0007] 本发明的目的是为了解决利用退化数据对产品的可靠性进行综合评估的问题，提出了针对加速退化试验及其退化数据的可靠性综合评估方法，以满足目前可靠性评估领域中，日益广泛使用加速退化试验及其退化数据来进行可靠性评估的现状，以及为提高评估精度而使用综合评估方法的实际需求。本发明采用贝叶斯方法，对退化数据进行综合，从而实现了基于多源退化数据的产品可靠性的综合评估。
[0008] 本发明是一种基于多源退化数据的贝叶斯可靠性综合评估方法，包括以下几个步骤：
[0009] 步骤一、多源退化数据的搜集；
[0010] 步骤二、建立退化统计模型并确定可靠度函数；
[0011] 步骤三、退化数据预处理；
[0012] 步骤四、先验分布、后验分布的确定及融合；
[0013] 步骤五、评估产品的可靠性；
[0014] 本发明的优点在于：
[0015] (1)本发明提出了退化统计模型的建立方法。通过对产品性能退化模型进行变换处理，对退化数据的形式进行转化，将并不服从概率分布的退化数据经此服从了某种概率分布。从而使通常处理寿命数据的贝叶斯方法能够处理退化数据，实现了基于退化数据的贝叶斯统计评估。因此也就将贝叶斯方法引入了加速退化试验的评估方法之中，扩展了加速退化试验的评估方法。
[0016] (2)本发明可以处理不同环境
应力，工作应力水平下的退化数据。通过回归方法以及应力水平间的关系对不同应力水平下的退化数据进行折合，将他们折合到同一应力水平下，从而提高了用于评估的信息量，扩展了可用于评估的数据范围。
[0017] (3)本发明提出针对多源退化数据，基于贝叶斯方法对产品的可靠性进行综合评估方法。将来自多个来源的退化数据，经处理统一到相同条件及形式下，通过这种扩大了数据量的手段使评估结果更加全面、准确。且目前的综合评估方法大多针对寿命数据，针对于退化数据的综合评估方法还很少。因此，本发明所提出的方法解决了对退化数据进行可靠性综合评估的问题。
附图说明
具体实施方式
[0020] 下面将结合附图和实施例对本发明作进一步的详细说明。
[0021] 本发明是一种基于多源退化数据的贝叶斯可靠性综合评估方法，流程图如图1所示，包括以下几个步骤：
[0022] 步骤一、多源退化数据的搜集；
[0023] 1.确定需要搜集的退化参数；
[0024] 在对产品的寿命及可靠性进行评估时，由于样本量、试验时间、试验条件等因素的限制，仅仅使用产品某次加速退化试验所获得的退化数据有时较为片面，而在产品的研制阶段、生产阶段、使用阶段存在大量可以利用的信息，通过这些信息来扩大评估的数据量，则可以提高评估的精度。
[0025] 因此，在确定需要搜集的退化参数时，应以加速退化试验以及寿命及可靠性评估中所使用的退化参数为准。
[0026] 2.搜集产品在研制阶段、生产阶段、使用阶段中与已确定的退化参数有关的退化信息；
[0027] 在产品的研制阶段、生产阶段、使用阶段的过程中，会有一些试验项目、测试项目，这些项目中往往存在大量有用的退化信息。
[0028] 3.对搜集到的退化信息进行筛选、整理，得到所需要的退化数据。退化数据应满足以下条件；
[0029] 1)选取同一退化参数的退化信息；
[0030] 以已确定的退化参数为准，选取退化信息。
[0031] 2)获取的退化信息为退化数据，或退化信息能够转化为退化数据；
[0032] 对于形式为退化数据的退化信息，可以直接应用。如果退化信息为退化数据已知的某个函数或是函数表达式，可以通过已知的函数将其转换为退化数据、或是通过仿真的方法得到退化数据，也可以选取。
[0033] 3)不同来源的退化数据的退化机理应该相同；
[0034] 由于产品参数的退化可以由不同的退化机理造成，为了提高评估的精度以及一致性，应以加速退化试验中产品的退化机理为准，来选取退化数据。
[0035] 4.根据退化数据的来源，将这些退化数据进行划分，如数据来自于m个信息源，那么就将退化数据划分为m个组，得到m组信息源数据；
[0036] 步骤二、建立退化统计模型并确定可靠度函数；
[0037] 在贝叶斯统计方法中，样本信息应服从总体分布。因此，首先应将退化数据进行
转化处理，使其能够服从某种概率分布，得到样本信息及其总体分布的形式，从而使其能够应用于贝叶斯方法之中。具体步骤如下：
[0038] 1.确定产品退化参数的退化统计模型；
[0039] 使用回归的方法对产品退化参数随时间变化的退化模型进行拟合建模。在具体对产品退化过程进行拟合建模时，根据需要可以按以下三种方法来建立退化参数随时间变化的退化统计模型：
[0040] 1).假设退化参数随时间变化的确定函数模型中的参数服从某种分布；
[0041] 2).对退化参数随时间变化的确定函数，添加一个不随时间变化服从某分布的误差量；
[0042] 3).对退化参数随时间变化的确定函数，添加随机过程描述；
[0043] 2.建立退化统计模型并确定可靠度函数；
[0044] 根据建立退化模型时的不同情况，对退化模型进行变换，并对退化数据进行处理，得到能够运用到贝叶斯方法中的分布及数据形式；并以此为
基础结合失效域值l及失效判据，建立产品的退化统计模型并确定可靠度函数；
[0045] 针对步骤1中提到的三种建模方法，分别提出了其处理方案：
[0046] 1)假设退化参数随时间变化的确定函数模型中的参数服从某种分布；
[0047] 通过函数中参数服从的某种分布，将函数中的其它参数及变量纳入到这个分布中，从而建立起一个新的分布。并将变换形式后的退化数据作为新的数据来使用，建立了退化统计模型，得到可靠度函数R(t)；
[0048] 举例说明：假设某产品退化参数随时间变化的确定函数模型为：
[0049] y(t)＝β0+β1t (1)
[0050] 其中y(t)为产品退化参数，β0、β1为参数项，t为时间；若β1□N(θ1，σ12)，根据正态分布的性质，可以得到：
[0051] y(t)□N(β0+θ1t，σ12t2) (2)
[0052] 进而得到：
[0053]
[0054] 将 作为新的数据形式，公式(3)作为其分布，从而确定了分布及数据的形式，建立了退化统计模型；
[0055] 如果设l为参数的失效
阈值，即设y(t)-l＜0时产品失效，那么就可以由上面推导的公式得到产品的可靠度函数：
[0056]
[0057] 2)对退化参数随时间变化的确定函数，添加一个不随时间变化服从某分布的误差量；
[0058] 通过函数中误差量服从的某种分布，将函数中的参数及变量纳入到这个分布中，从而建立起一个新的分布。并将变换形式后的退化数据作为新的数据来使用，建立了退化统计模型，得到可靠度函数R(t)；
[0059] 现举例说明：
[0060] 假设某产品退化参数随时间变化的确定函数模型为：
[0061] y(t)＝Ct+ε (5)
[0062] 其中y(t)为产品退化参数，C为常数，ε□N(θ，σ2)，t为时间；则可根据正态分布的性质，得到
[0063] y(t)□N(Ct+θ，σ2) (6)
[0064] 进而得到
[0065] y(Δt)□N(CΔt，0) (7)
[0066] 将y(Δt)作为新的数据形式，公式(7)作为其分布，从而确定了分布及数据的形式，建立了退化统计模型；
[0067] 如果设l为参数的失效阈值，即设y(t)-l＜0时产品失效，那么就可以由上面推导的公式得到产品的可靠度函数：
[0068]
[0069] 3)对退化参数随时间变化的确定函数，添加随机过程描述；
[0070] 对于添加随机过程描述的性能退化参数随时间变化的确定函数，根据随机过程的性质对函数进行处理，建立了退化统计模型，得到可靠度函数R(t)；
[0071] 如对于漂移布朗运动模型：
[0072] Y(t)＝σB(t)+g(t，s)·t+y0 (9)[0073] 其中Y(t)为产品参数的退化过程；B(t)为均值为0，方差为时间t的标准布朗运动B(t)□N(0，t)；σ为扩散系数，不随应力和时间而改变，为常数；g(t，s)为漂移系数；y0为产品性能的初始值；
[0074] 因此由漂移布朗运动的性质可知，单位时间Δt的退化增量ΔY服从均值为g(t，s)·Δt，方差为σ2Δt的正态分布，即
[0075] ΔY□N(g(t，s)·Δt，σ2Δt) (10)[0076] 将ΔY作为新的数据形式，公式(10)作为其分布，从而确定了分布及数据的形式，建立了退化统计模型；
[0077] 如果设l为参数的失效阈值，即设Y(t)-l＜0时产品失效；可以利用漂移布朗运动可靠度模型对可靠度函数进行求解：
[0078]
[0079] 步骤三、退化数据预处理；
[0080] 由于需要处理的退化数据的来源不尽相同，数据间存在一定的差距，需要将它们进行一定的处理后，才可应用于评估。具体步骤如下：
[0081] 1.将所得到的退化数据转化为服从总体分布的形式；
[0082] 根据步骤二中得到的分布及新的数据形式，根据得到的新数据形式公式将退化数据进行变换，使其服从所得到的分布；
[0083] 2.确定应力折合模型；
[0084] 通常，产品的退化过程易受如
温度，振动，湿度，电应力等应力的影响。在不同的应力水平下，退化的过程是不同的，退化数据所代表的意义也因此而不同，要综合利用多源退化数据，首要解决的就是统一应力水平，将它们折合到需要评估的实际的条件中来。
[0085] 各应力都有其各自的加速模型或是折合因子，根据退化参数的具体情况及其所处的应力确定加速模型或是折合因子，作为应力折合模型。
[0086] 现列举加速模型如下：
[0087] 1)阿伦尼斯模型
[0088] 对于温度，通常使用阿伦尼斯(Arrhenius)模型，其形式如下：
[0089] ξ＝AeE/KT
[0090] 其中ξ是某寿命特征，A是一个正常数，E是与材料有关的激活能，K是波尔兹曼常数，T是绝对温度。
[0091] 2)逆幂律模型
[0092] 对于电应力，通常使用逆幂律模型，其形式如下：
[0093] ξ＝Av-c (12)
[0094] 其中ξ是某寿命特征，A是一个正常数，c是与激活能有关的正常数，v是应力。
[0096] 对于温度和
电压同时作为加速应力时，通常使用广义艾林模型，其形式如下：
[0097]
[0098] 其中ξ是某寿命特征，A、B、C、D为待定常数，k是波尔兹曼常数，c是与激活能有关的正常数，T是绝对温度，V为电压。
[0099] 4)Peck模型
[0100] 对于温度和湿度同时作为加速应力时，通常使用Peck模型，其具体形式为[0101]
[0102] 其中ξ是某寿命特征，为
相对湿度，Ea为激活能，k为波尔兹曼常数，T为绝对温度，m和A为常数。
[0103] 5)其它
[0104] 对于一些较为典型的产品以及材料，对于其所处的典型应力，都有其自身的加速模型，在很多标准及文献中都可以查到。
[0105] 同时，对于一些产品，其参数在各个应力水平间的关系(折合因子)，都以给出或作出相应得规定，在这种情况下，可以不通过加速模型，而直接应用折合因子对不同应力水平进行折合。
[0106] 因此，根据各来源退化数据以及加速退化试验中各应力水平，结合产品自身的特点，确定相应得加速模型或折合因子。
[0107] 3.将不同应力水平下的退化数据折合到同一水平下；
[0108] ①若上述步骤(2)中应力折合模型为加速模型G(s)，首先根据不同应力水平s下的退化数据，应用回归拟合的方法拟合并建立起退化过程中某一与时间有关的特征值zs与应力水平s之间的关系式：
[0109] zs＝G(s) (15)[0110] 通过公式(15)可以得到不同应力水平s下的某一与时间有关的特征值zs的值；
[0111] 其次，将应力si水平下的退化数据yi折合到所需的应力水平si下的退化数据yi，若能够建立退化数据y与某一与时间有关的特征值zs的关系：
[0112] y＝f(zs) (16)[0113] 其中，公式(16)可以通过退化模型以及某一与时间有关的特征值zs的具体意义确定；那么由公式(15)、(16)可知：
[0114]
[0115] 从而完成了不同应力水平下的退化数据间的折合；
[0116] ②若上述步骤(2)中应力折合模型为折合因子kij，通过给定的应力水平间的折合因子kij将应力水平si下的某一特征值zsi折合到应力水平sj下的某一特征值zsj：
[0117] zsi＝kij·zsj (18)[0118] 然后通过公式(16)来完成了不同应力水平间退化数据的折合；
[0119] 4.对各信息源数据的形式进行统一；
[0120] 如果步骤3中不能建立公式(16)那样的退化数据y与某一与时间有关的特征值zs的关系，则利用每个退化数据间的时间间隔Δt来建立它与某一与时间有关的特征值zs的关系，则有：
[0121] Δt＝f(zs) (19)[0122] 由公式(15)、(19)得到：
[0123]
[0124] 随着退化数据间的时间间隔Δt的改变，每个退化数据代表了不同时间内的参数退化值，需要对退化数据形式进行统一，具体方法为：
[0125] 首先，由公式(20)得到各退化数据时间间隔Δti，求出它们的最小公倍数ΔtG，得到ΔtG与Δti的倍数关系pi＝ΔtG/Δti；
[0126] 然后，将各退化数据的每pi个数据点进行合并，将pi个退化数据的退化量合并到一个点，则所有退化数据的时间间隔都相等，退化数据形式统一；
[0127] 5.相容性检验；
[0128] 将步骤4得到的各信息源数据，作为贝叶斯方法中的先验信息，将加速退化试验中经过步骤4的退化数据作为贝叶斯方法中的样本信息，将步骤二得到的分布作为总体分布；由于所选择的先验信息来自于多个不同的信息源。所以要将每个来源的先验信息与样本信息进行相容性检验，通过检验两母体是否相容，以确定先验信息是否可用。
[0129] 本发明使用秩和检验法来检验先验信息与样本信息是否相容。具体方法为：
[0130] 为第i组先验信息， 为样本信息，引入互相竞择假设：
[0131] H0：X与Yi属于同一总体；H1：X与Yi不属于同一总体。
[0132] 若Xi与Yi通过检验，认为两子样不属于同一总体；否则，两子样属于同一总体，将这一组先验信息去除。
[0133] 步骤四、先验分布、后验分布的确定及融合；
[0134] 1.根据先验信息，总体分布，样本信息的实际情况，选择先验分布的确定方法；
[0135] 根据先验信息，总体分布，样本信息的实际情况，确定先验分布。可以通过无信息先验分布法、共轭先验分布法、最大熵先验分布法、自助法和随机加权法确定先验分布法、经验贝叶斯方法等贝叶斯先验分布确定方法来确定先验分布。
[0136] 2.确定各信息源的先验分布；
[0137] 对于m个信息源的先验信息，需要根据每一组先验信息使用相同的先验分布确定方法各自确定其先验分布πi(θ)，且πi(θ)服从相同的先验分布族。
[0138] 3.稳健性检验；
[0139] 采用稳健性检验的方法来判断所确定的先验分布是否符合要求。若已知先验分布πi(θ)、总体分布f(x|θ)与样本信息x，那么边缘
密度函数的公式为：
[0140]
[0141] 通过求得的mi(x)的值来判断所确定的先验分布是否稳健，根据需要设定一个可接受值cv，若mi(x)＜cv则需要使用其它方法重新确定先验分布；
[0142] 4.通过样本信息，结合各信息源的先验分布，确定各信息源的后验分布；
[0143] 通过先验分布πi(θ)、总体分布f(x|θ)与样本信息x得到θ的后验分布πi(θ|x)，后验分布公式为：
[0144]
[0145] 5.求解并通过加权系数得到各信息源融合后的后验分布；
[0146] 对多个先验分布进行加权融合，将其融合成一个先验分布。这里使用加权融合方法中的相关系数法作为融合先验分布的方法。
[0147] 已知总体分布f(x|θ)，通过m个信息源分别得到先验分布π1(θ)，π2(θ)，...πm(θ)，那么融合后的先验分布为：
[0148]
[0149] 其中εi为加权系数：
[0150]
[0151] 采用相关系数法来求解εi；
[0152] 此时后验分布为：
[0153]
[0154] 所以融合后的后验分布为
[0155]
[0156] 其中， m(x|π)＝∑m(x|πi)；令 那么
[0157]
[0158] 那么即
[0159]
[0160] 其次，对πi(θ)取期望ui，对πi(θ|x)取期望u′i，方差σi2，对π(θ|x)取期望 通过仿真抽样，得到m+1组ui、u′i、u′的估计值，以及σ12、σ22的估计值；将数据带入方程组：
[0161]
[0162] 从而求解得到εi及ξi的值；进而公式(28)得到融合后的后验分布π(θ|x)；
[0163] 步骤五、评估产品的可靠性；
[0164] (1)通过融合后的后验分布π(θ|x)，结合贝叶斯估计方法，得到总体分布的参数θ的评估值 其中：
[0165]
[0166] (2)将 代入步骤二的可靠度函数R(t)中，求出产品在时刻t的可靠度值；如需要求时刻t的可靠度R(t)的区间估计，通过后验分布π(θ|x)得到参数θ的区间估计，从而得到产品在时刻t的可靠度R(t)的区间估计。
[0167] 实施例：
[0168] 以某型光电产品为例，通过对其步进应力加速退化试验的评估，具体介绍基于退化数据的贝叶斯可靠性综合评估方法。
[0169] 步骤一、多源退化数据的搜集；
[0170] 1.确定需要搜集的性能退化参数；
[0171] 通过分析发现光功率退化失效是此产品的主要失效原因，在各种试验中也是以此作为产品主要的性能指标，因此选择光功率作为应搜集的性能参数。
[0172] 2.搜集产品在研制阶段、生产阶段、使用阶段中与已确定的退化参数有关的退化信息；
[0173] 经过搜集分析，认为此产品在加速退化预试验中的性能参数退化数据具有较强的利用价值，各种条件与温度步进应力加速退化试验接近，且样本量较多，其退化数据可以作为步进应力加速退化试验数据的补充，能够应用到综合评估中来。因此，选择此产品的加速退化预试验中的数据作为先验信息的来源。
[0174] 将步进应力加速退化试验的数据作为样本信息的来源。
[0175] 3.对搜集到的退化信息进行筛选、整理，得到所需要的退化数据；
[0176] 选择预试验中的温度步进应力加速退化预试验中的2支产品，将其编号为1、2；选择温度恒定应力加速退化预试验中的2支产品，将其编号为3、4；将这4支产品的光功率退化数据作为先验信息。选取各支产品所处的应力如表1所示
[0177] 表1各产品应力信息
[0178]
[0179] 将温度步进应力加速退化试验的数据作为样本信息。其试验温度应力分别为60℃、80℃、100℃和110℃。选择其中4支的产品数据，将其编号为5、6、7和8。
[0180] 4.根据退化数据的来源，将这些退化数据进行划分，如数据来自于m个信息源，那么就将退化数据划分为m个组；
[0181] 由于在预试验中编号为1，2的产品所表现的一致性并不强，而编号为3、4的产品属于同一母体。因此这里将它们划分为三个信息源进行处理，信息源1为编号为1的产品，信息源2为编号为2的产品，其余的产品则归为信息源3。将这3个信息源作为先验信息。
[0182] 将编号为5、6、7和8的退化数据作为样本信息。
[0183] 步骤二、建立退化统计模型并确定可靠度函数
[0184] 1.确定产品退化参数的退化模型；
[0185] 选择漂移布朗运动对产品性能参数的退化过程进行描述。已知漂移布朗运动退化模型为：
[0186] Y(t)＝σB(t)+d(s)·t+y0
[0187] 其中Y(t)为产品的性能退化过程；B(t)为均值为0，方差为时间t的标准布朗运动，B(t)□N(0，t)；σ为扩散系数，不随应力和时间而改变，为常数；d(s)为漂移系数；y0为产品性能的初始值。
[0188] 2.建立退化统计模型并确定可靠度函数；
[0189] 由于漂移布朗运动属于一种随机过程，因此对于添加随机过程描述的性能退化参数随时间变化的确定函数，可以根据随机过程的性质对函数进行处理。
[0190] 由漂移布朗运动的性质可知，单位时间Δt的退化增量ΔY服从均值为d(s)·Δt，方差为σ2Δt的正态分布，即
[0191] ΔY□N(d(s)·Δt，σ2Δt)
[0192] 因此将正态分布N(d(s)·Δt，σ2Δt)作为总体分布，将退化增量ΔY作为样本信息。
[0193] 如果设l为性能参数的失效域值，光功率的退化是递减的，即当Y(t)-l＜0时产品失效。可以利用漂移布朗运动可靠度模型对可靠度函数进行求解：
[0194]
[0195] 步骤三、退化数据的预处理；
[0196] 1.将所得到的退化数据转化为新的数据形式；
[0197] 由于退化数据的退化增量ΔY服从均值为d(s)·Δt，方差为σ2Δt的正态分布。因此将所有的退化数据转化为退化增量ΔY的形式。
[0198] 2.确定应力折合模型；
[0199] 由于产品在步进应力试验及预试验中工作状态相同，所受的环境应力为温度，选择Arrhenius模型作为温度加速模型，即应力折合模型。此时的应力s为温度T，令[0200] d(T)＝Aexp(-B/T)。
[0201] 其中：T是绝对温度，K
[0202] B＝Ea/k，k是波尔兹曼常数8.6171×10-5eV/K，Ea是激活能，eV[0203] A是常数
[0204] 3.将不同应力水平下的退化数据折合到同一水平下；
[0205] 对先验信息进行初步的处理与折合，将不同温度下的试验数据折合到同一温度水平，为接下来的工作作准备。这里使用线性拟合的方法得到温度与退化趋势之间的关系，从而得到三个信息源在其所在温度下的退化率与其在25℃下的退化率之间的关系。
[0206] 建立线性回归方程如下
[0207] E(Y(t))＝d(T)·t+y0
[0208] E(Y(t))＝Aexp(-B/T)·t+y0
[0209] 即
[0210] E(ΔY(t))＝Aexp(-B/T)·Δt
[0211] 进而得到线性回归模型
[0212] ln(E(ΔY(t)))＝-B/T+ln(AΔt)
[0213] 4.对各信息源数据的形式进行统一；
[0214] 通过拟合可以得到加速模型的参数，这里使用检测时间折合的方法来实现各应力下退化增量向25℃下退化增量的折合。检测时间折合公式为
[0215]
[0216] 其中□t为检测时间，□tk为将应力sk折合到25℃后的检测时间。
[0217] 对于先验信息而言，原检测时间间隔□t＝6分钟，折合后信息源1和2的检测时间间隔变为120分钟；信息源3中3号产品的检测时间间隔变为240分钟，4的检测时间间隔变为300分钟。
[0218] 同理对样本信息进行折合，使用线性回归的方法，得到加速模型的参数，对检测时间进行折合。对于温度步进应力加速退化试验，其各温度下的检测时间变为：60℃的检测时间间隔变为100分钟，80℃的检测时间间隔变为300分钟，100℃的检测时间间隔变为800分钟，110℃的检测时间间隔变为1200分钟。
[0219] 由于2400约为120，240，300，100，800，1200的最小公倍数，因此选择2400分钟作为其检测时间间隔。将每2400分钟内的退化增量进行合并，使其中的多个数据成为一个数据，使所有的数据都以2400内的退化增量的形式表达，从而对各信息进行了统一折合。具体合并情况如表2所示。
[0220] 表2数据的统一折合
[0221]
[0222] 从而统一了先验数据与样本数据的数据结构，为下一步的
数据处理及评估作好了准备。
[0223] 5.相容性检验
[0224] 使用Wilcoxon-Mann-Whitney的秩和检验法来检验先验信息与样本信息是否相容。
[0225] 通过MATLAB中的秩和检验程序ranksum对各组数据进行检验，检验的结果为，信息源1与信息源3同样本信息相容，信息源2与样本信息不相容。因此将信息源2去除，选择信息源1，信息源3进入下一步的评估。
[0226] 步骤四、先验分布、后验分布的确定及融合
[0227] 1.根据先验信息，总体信息，样本信息的实际情况，选择先验分布的确定方法；
[0228] 由于共轭先验分布计算起来比较方便，同时后验分布的一些参数可得到很好的解释，本次评估选用共轭先验分布法来确定产品的先验分布。其中共轭先验分布的定义为：对于某一分布而言，假如由抽样信息算得的后验密度函数与先验分布π(θ)有相同的函数形式，则称π(θ)是θ的(自然)共轭先验分布。
[0229] 2.确定各信息源的先验分布；
[0230] 已知ΔY□N(d(s)·Δt，σ2Δt)，那么根据共轭先验分布理论，由表3可知漂移系数与时间间隔的乘积服从正态分布，扩散系数与时间间隔的乘积服从倒伽
马分布，即2 2
d(s)·□t□N(μ0，τ)，σ·□t□IG(b，a)。
[0231] 表3共轭先验分布
[0232]总体分布 参数 共轭先验分布
正态分布(方差已知) 均值 正态分布
正态分布(均值已知) 方差 倒伽玛分布
[0233] 其概率密度函数具体形式为：
[0234]
[0235]
[0236] 由先验数据Y＝(y1，y2，…ym)来确定先验分布，所以结合数据的特点选用公式[0237]
[0238]
[0239]
[0240]
[0241] 来确定先验分布的各超参数的值，从而确定了先验分布。超参数的值如表4所示[0242] 表4先验分布超参数的值
[0243]
[0244] 3.稳健性检验；
[0245] 对于扩散系数与时间乘积的先验分布σ2·□t□IG(b，a)而言，这里通过它的边缘密度函数值来判断其稳健性。
[0246] 设x＝(x1，x2，…xn)为样本信息，X□N(θ，σ2)。那么它的边缘密度函数为[0247]
[0248] 通过得到的边缘密度函数即可检验先验分布的稳健性。这里将y的均值代入进行检验。信息源1、3的检验结果为0.98，0.95。设定可接受值cv＝0.6，结果说明边缘密度函数值大于可接受值，先验分布通过了检验。2
[0249] 对于漂移系数与时间间隔乘积的先验分布d(s)·□t□N(μ0，τ)而言，其稳健性检验的具体方法如下。2
[0250] 设X＝(x1，x2，…xn)为样本信息，X□N(θ，σ)。首先计算m(x|π)，这里x为2 2
子样的均值，而其分布密度函数为N(μ0，σ/n+τ)。设中心变化量Z□N(0，1)，那么可以得到
[0251]
[0252] 通过正态分布的Z检验，就可以检验先验分布的稳健性。检验的结果为先验分布通过检验，且通过的概率分别为0.83、0.75。
[0253] 4.通过样本信息，结合各信息源的先验分布，确定各信息源的后验分布；
[0254] 求解后验分布时，由于ΔX□N(d(s)·Δt，σ2Δt)中，d(s)·Δt、σ2Δt均未知，因此首先通过正态-倒伽玛分布来确定σ2Δt后验分布的超参数，若Y＝(y1，y2，…ym)为先验信息，X＝(x1，x2，…xn)为后验信息，则：
[0255]
[0256]
[0257] 由超参数α1、β1的值可以得到σ2Δt的后验分布、期望值、方差值：
[0258]
[0259]
[0260]
[0261] 由得到的σ2Δt的评估值，结合共轭先验分布的性质可以得出d(s)·Δt后验分布的表达式为
[0262]
[0263] 其中，
[0264] 所得到的后验分布超参数如表5所示
[0265] 表5后验分布超参数的值
[0266]
[0267] 5.求解并通过加权系数得到各信息源融合后的后验分布；
[0268] 通过加权融合方法中的相关系数法分别求解两个参数相应的加权系数。2 2 2
[0269] 对于σΔt首先通过所得到的先验分布π1(σΔt)、π2(σΔt)及其后验分布2 2
π1(σΔt|x)、π2(σΔt|x)进行仿真抽样，得到其期望值的仿真评估值μ1、μ2、μ′1、
2 2
μ′2以及后验分布方差仿真评估值σ1、σ2 各3组，并通过后验分布的边缘密度函数值m(x|π1)、m(x|π2)得到：
[0270]
[0271]
[0272] 从而通过公式：
[0273] π(σ2Δt|x)＝ξ1π1(σ2Δt|x)+ξ2π2(σ2Δt|x)
[0274] 对π(σ2Δt|x)进行仿真抽样，得到3组含有未知数ε1、ε2的期望仿真评估值μ′。将所得到的3组μ′、μ1、μ2、μ′1、μ′2、σ12、σ22的值代入方程组：
[0275]2
[0276] 使用MATLAB求解方程组，从而得到σΔt的先验分布与后验分布的加权系数，信息源1为0.1685，信息源3为0.8315。同理可以得到d(s)·Δt的先验分布与后验分布的加权系数，信息源1为0.8866，信息源3为0.1134。
[0277] 因此可以得到融合后的后验分布
[0278] π(d(s)·□t|x)＝0.8866π1(d(s)·□t|x)+0.1134π2(d(s)·□t|x)[0279] π(σ2·□t|x)＝0.1685π1(σ2·□t|x)+0.8315π2(σ2·□t|x))[0280] 步骤五、评估产品的可靠性
[0281] 1.通过融合后的后验分布π(θ|x)，结合贝叶斯估计方法，得到总体分布的参数θ的评估值
[0282] 通过融合后的后验分布，由后验函数期望的性质可以得到
[0283] E(d(s)·□t|x)＝0.8866E1(d(s)·□t|x)+0.1134E2(d(s)·□t|x)[0284] E(σ2·□t|x)＝0.1685E1(σ2·□t|x)+0.8315E2(σ2·□t|x)
[0285] 而
[0286] Ei(d(s)·□t|x)＝u1i
[0287]
[0288] 那么
[0289]
[0290]
[0291] 由于2400分钟为40个小时，这里取□t＝40小时，从而得到温度为25℃时的产品退化的漂移系数及扩散系数，
[0292] 2.将 代入步骤二的可靠度函数R(t)中，即可求出产品在时刻t的可靠度值。
[0293] 若定产品的光功率退化到初始值的一半，即l＝0.5y0时，产品失效。将l的值代入可靠度函数
[0294]
[0295] 从而得到产品在25℃工作20万小时的可靠度为0.9383，10万小时的可靠度为0.9994。图2为由可靠度函数所得到的产品可靠度随着时间变化的曲线图。
