基于比例危害－比例优势模型的加速寿命试验优化设计方法

摘要

本发明公开了一种基于比例危害－比例优势模型的加速寿命试验优化设计方法，包括建立基于比例危害－比例优势模型的对数似然函数；建立基于比例危害－比例优势模型的Fisher信息矩阵和渐进方差－协方差矩阵；采用基于信息的优化方法建立最优化问题目标函数和约束条件；求解最优化问题，得到加速寿命试验优化设计方案。本发明给出了包括恒定应力和步进应力加速寿命试验的优化设计方法，本发明进行加速寿命试验优化设计时基于的模型是比例危害－比例优势模型，该模型较比例危害模型和比例优势模型具有更广的适用范围和评估精度；该模型是一个非参数模型，具有无分布特性；本发明方法避免了传统优化方法中由于积分区间变化而造成的优化结果不一致的问题。

说明书

技术领域

本发明属于加速寿命试验领域，具体涉及一种基于比例危害-比例优势模型的加速寿命试验优化设计方法。

背景技术

加速寿命试验技术是使用在高应力下较短时间内获得的产品寿命数据外推产品在正常应力水平下的寿命和可靠性的技术。加速模型描述了产品寿命与应力水平之间的关系。加速模型的建立是进行外推的基础，直接影响着外推的精度。加速模型中的统计加速模型可以分为参数模型和非参数模型。参数模型需要预先确定产品的寿命分布形式，然而，如果产品寿命并不服从该寿命分布形式，则会造成较大的评估误差。非参数模型是一种无分布假设的模型，因此，更加受到研究者及工程应用人员的青睐。两个广泛使用的非参数模型是参考文献1：Cox，D.R.，Regression Models and Life-Tables，Journal of the Royal Statistical Society，Series B(Methodological)，Vol.34，No.2，pp.187-220，1972.中提出的比例危害模型和参考文献2：Brass，W.，Mortality Models and their Uses in Demography，Transactions of the Faculty of Actuaries，Vol.33，122-133，1974.中提出的比例优势模型。在医学领域，研究者基于参考文献3：Aranda-Ordaz，F.J.，On Two Families ofTransformations to Additivity for Binary Response Data，Biometrics，Vol.68，No.2，pp.357-363，1981.中的Aranda-Ordaz参数族提出了一些生存数据评估模型。参考文献4：Huang，T.，Elsayed，E.A.，and Jiang，T.，An ALT Proportional Hazard-ProportionalOdds Model，The Proceedings of the 14th ISSAT International Conference onReliability and Quality in Design，Florida，USA，Aug.7-9th，2008.将Aranda-Ordaz参数族及其衍生模型引入加速寿命试验领域，给出了加速寿命试验领域的比例危害-比例优势模型。该模型通过转移参数将比例危害模型和比例优势模型结合起来，使这两种模型成为该模型的特殊情况。比例危害-比例优势模型较比例危害模型和比例优势模型具有更广泛的适用范围和评估精度。

加速寿命试验优化设计的目的是为了减少运用加速模型进行产品寿命和可靠性外推时的误差。进行加速寿命试验优化设计通常包括以下几方面的内容：应力水平数及其量值，恒定应力试验时每个应力水平下分配的受试产品个数，或步进应力试验时的应力转换时间。目前，研究者运用了多种方法进行加速寿命试验的优化设计。对于非参数模型而言，通常使用的方法是最小化一个与可靠性相关的函数的渐进方差在一个给定区间内的积分值。例如，参考文献5：Elsayed，E.A.and Jiao，L.，Optimal Design of Proportional Hazards basedAccelerated Life Testing Plans，International J.of Materials&Product Technology，17，411-424，2002.给出的基于比例危害模型的加速寿命试验优化设计方法，该文是通过最小化危害率函数的渐进方差在一个给定区间内的积分值得到最优试验方案的；又如参考文献6：Elsayed，E.A.and Zhang，H.，Design of Optimum Reliability Test Plans underMultiple Stresses，QUALITA 2005，Quality and Dependability，Bordeaux，France，March 16-18，2005.给出的基于比例优势模型的加速寿命试验优化方法，该文是通过最小化可靠性函数的渐进方差在一个给定区间内的积分值得到最优试验方案的。然而，这两种优化方法的局限性在于，当积分区间变化时，将得到不同的优化结果。

近几年，基于信息的优化方法已应用于加速寿命试验的优化设计中，该方法可以避免上述优化方法的局限。广泛应用的基于信息的优化方法包括D-优化方法和A-优化方法。这两种优化方法的目的都是最大化信息，同时最小化与信息互为倒数的方差。这两种优化方法通过不同的途径得到优化结果，然而，究竟采用哪种优化方法取决于模型参数的特性。

D-优化方法最大化了信息矩阵的行列式值，由于信息矩阵与方差-协方差矩阵互为倒数，因此，它同时最小化了参数评估值置信区间的体积。该优化方法通常用于模型参数有较强相关性的情况。D-优化方法的目标函数是，

max det[I]

A-优化方法最小化了信息矩阵逆的迹，也就是说它最小化了方差-协方差矩阵的迹。该优化方法通常用于模型参数的相关性较弱的情况。A-优化方法的目标函数是，

min tr[I^-1]即，min tr[∑]

参考文献7：Ng，H.K.T.，Chan，P.S.，and Balakrishnan，N.，Optimal ProgressiveCensoring Plans for the Weibull Distribution，American Statistical Association andthe American Society for Quality Technometrics，Vol.46，No.4，2004.将D-优化方法和A-优化方法应用于威布尔分布的产品的加速寿命试验优化设计中；参考文献8：Ng，H.K.T.，Balakrishnan，N.，and Chan，P.S.，Optimal Sample Size Allocation forTests with Multiple Levels of Stress with Extreme Value Regression，WileyInterScience，2006.将这两种优化方法应用于处理极值分布情况的产品样本量分配问题。然而，这些基于信息的优化方法目前主要应用于基于参数模型的加速寿命试验优化设计中。

发明内容

本发明将基于信息的优化方法引入基于非参数模型的加速寿命试验优化设计中，分别针对恒定应力和步进应力这两种试验方法，给出了运用基于信息的优化方法进行基于比例危害-比例优势模型的加速寿命试验优化设计的方法。

所述的加速寿命试验优化方法包括如下步骤：

步骤一、建立基于比例危害-比例优势模型的对数似然函数；

步骤二、建立基于比例危害-比例优势模型的Fisher信息矩阵和渐进方差-协方差矩阵；

步骤三、采用基于信息的优化方法建立最优化问题目标函数和约束条件；

步骤四、求解最优化问题，得到加速寿命试验优化设计方案。

本发明的优点在于：

(1)本发明进行加速寿命试验优化设计时基于的模型是一个非参数模型，非参数模型较参数模型的优势在于非参数模型具有无分布特性；

(2)本发明进行加速寿命试验优化设计时基于的模型是比例危害-比例优势模型，该模型较比例危害模型和比例优势模型这两个广泛应用的非参数模型具有更广的适用范围和评估精度；

(3)本发明采用的基于信息的优化方法避免了传统优化方法中由于积分区间变化而造成的优化结果不一致的问题。

附图说明

图1为加速寿命试验优化设计方法流程图；

图2为多应力类型恒定应力加速寿命试验剖面；

图3为多应力类型步进应力加速寿命试验剖面。

具体实施方式

本发明提供了一种基于比例危害-比例优势模型的加速寿命试验优化设计方法，包括如下步骤(如图1)：

步骤一、建立基于比例危害-比例优势模型的对数似然函数；

步骤二、建立基于比例危害-比例优势模型的Fisher信息矩阵和渐进方差-协方差矩阵；

步骤三、采用基于信息的优化方法建立最优化问题目标函数和约束条件；

步骤四、求解最优化问题，得到加速寿命试验优化设计方案。

为了具有通用性，本发明给出了针对恒定应力和步进应力两种试验方法的多应力类型加速寿命试验的优化设计方法，其中，基于比例危害-比例优势模型的恒定应力加速寿命试验优化设计具体方法如下：

假设n个受试产品进行恒定应力加速寿命试验，有k种应力类型(例如：温度应力，电应力等)，每种应力类型有q个不同应力水平。那么，加速寿命试验将在k×q个不同应力组合下进行。定义了分配给应力水平$z_{i_1{i}_2...i_k}={(z_{i_1},z_{i_2},···,z_{i_k})}^t$的受试产品个数与参与试验的总产品个数n的比值，其中，i_j＝1，2，...q，j＝1，2，...，k。应力水平为的试验在预先确定的时间截尾。图1给出了当k＝2，q＝2时的试验剖面，每个应力类型都有介于设计应力和工作极限之间的高、低两个应力水平，因此，应力类型1和应力类型2的高、低应力水平两两组合构成了该产品加速寿命试验的试验剖面。优化目标是合理的选取应力水平和每个应力水平分配的产品比例i_j＝1，2，...q，j＝1，2，...，k，使得在正常应力水平下对产品寿命和可靠性的评估结果最精确，具体的优化方法如下：

步骤一、建立基于比例危害-比例优势模型的对数似然函数。

比例危害-比例优势模型中可靠性函数R(t；z)和概率密度函数f(t；z)的定义为，

$R(t;z)=e^{\frac{1}{-c}\ln (c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)}---\left(1\right)$

$f(t;z)={(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)}^{\frac{1}{-c}-1}{\gamma }_1{\gamma }_2{t}^{{\gamma }_2-1}{e}^{{\beta }^{t}z}---\left(2\right)$

式中：β＝(β₁，β₂，…，β_k)是模型未知参数向量，其中β向量中的元素个数等于应力类型数；z是应力水平向量；γ₁＞0、γ₂＞0是模型未知参数；c∈[0，1]是转移参数，t表示时间。

定义一个指示函数I，

其中，τ为试验截尾时间。

因此，对于某一个失效数据(t，I，z)，对数似然函数可以写为，

l(t；z)＝Ilnf(t；z)+(1-I)lnR(τ；z)                        (3)

将(1)和(2)带入(3)，可以将对数似然函数改写为，

$l(t;z)=I[(\frac{1}{-c}-1)\ln (c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)+\ln {\gamma }_1+\ln {\gamma }_2+({\gamma }_2-1)\ln t+{\beta }^{t}z]+(1-I)\left(\frac{1}{-c}\right)\ln (c{\gamma }_1{\tau }^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)---\left(4\right)$

步骤二、建立基于比例危害-比例优势模型的Fisher信息矩阵及渐进方差-协方差矩阵；

计算(4)式给出的对数似然函数对每个模型参数的二阶偏导数，

$\frac{{\partial }^{2}l}{\partial {{\beta }_1}^2}=I(-\frac{1}{c}-1)\left[\frac{\left(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^{t}z}{z_1}^2\right)(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)-{\left(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^{t}z}{z}_1\right)}^2}{{(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)}^2}\right]$

$+(1-I)(-\frac{1}{c})\left[\frac{\left(c{\gamma }_1{\tau }^{{\gamma }_2}{e}^{{\beta }^{t}z}{z_1}^2\right)(c{\gamma }_1{\tau }^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)-{\left(c{\gamma }_1{\tau }^{{\gamma }_2}{e}^{{\beta }^{t}z}{z}_1\right)}^2}{{(c{\gamma }_1{\tau }^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)}^2}\right]$

.

.

.

$\frac{{\partial }^{2}l}{\partial {{\beta }_k}^2}=I(-\frac{1}{c}-1)\left[\frac{\left(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^{t}z}{z}_k^2\right)(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)-{\left(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^{t}z}{z}_k\right)}^2}{{(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)}^2}\right]$

$+(1-I)(-\frac{1}{c})\left[\frac{\left(c{\gamma }_1{\tau }^{{\gamma }_2}{e}^{{\beta }^{t}z}{z}_k^2\right)(c{\gamma }_1{\tau }^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)-{\left(c{\gamma }_1{\tau }^{{\gamma }_2}{e}^{{\beta }^{t}z}{z}_k\right)}^2}{{(c{\gamma }_1{\tau }^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)}^2}\right]$

$\frac{{\partial }^{2}l}{{\partial c}^2}=I[(-\frac{2}{c^3})\ln (c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)+\left(\frac{2}{c^2}\right)\frac{{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^{t}z}}{c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^{t}z}+1}+(\frac{1}{c}+1)\frac{{\left({\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^{t}z}\right)}^2}{{(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)}^2}]$

$+(1-I)[(-\frac{2}{c^3})\ln (c{\gamma }_1{\tau }^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)+\left(\frac{2}{c^2}\right)\frac{{\gamma }_1{\tau }^{{\gamma }_2}{e}^{{\beta }^{t}z}}{c{\gamma }_1{\tau }^{{\gamma }_2}{e}^{{\beta }^{t}z}+1}+\left(\frac{1}{c}\right)\frac{{\left({\gamma }_1{\tau }^{{\gamma }_2}{e}^{{\beta }^{t}z}\right)}^2}{{(c{\gamma }_1{\tau }^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)}^2}]$

$\frac{{\partial }^{2}l}{\partial {{\gamma }_1}^2}=I[(\frac{1}{c}+1)\frac{{\left({\mathrm{ct}}^{{\gamma }_2}{e}^{{\beta }^{t}z}\right)}^2}{{(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)}^2}-\frac{1}{{{\gamma }_1}^2}]+(1-I)\left(\frac{1}{c}\right)\frac{{\left(c{\tau }^{{\gamma }_2}{e}^{{\beta }^{t}z}\right)}^2}{{(c{\gamma }_1{\tau }^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)}^2}$

$\frac{{\partial }^{2}l}{\partial {{\gamma }_2}^2}=I[(-\frac{1}{c}-1)\frac{\left(c{\gamma }_1{t}^{{\gamma }_2}{\left(\ln t\right)}^2{e}^{{\beta }^{t}z}\right)(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)-{\left(c{\gamma }_1{t}^{{\gamma }_2}\left(\ln t\right)e^{{\beta }^{t}z}\right)}^2}{{(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)}^2}-\frac{1}{{{\gamma }_2}^2}]$

$+(1-I)(-\frac{1}{c})\left[\frac{\left(c{\gamma }_1{\tau }^{{\gamma }_2}{\left(\ln \tau \right)}^2{e}^{{\beta }^{t}z}\right)(c{\gamma }_1{\tau }^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)-{\left(c{\gamma }_1{\tau }^{{\gamma }_2}\left(\ln \tau \right)e^{{\beta }^{t}z}\right)}^2}{{(c{\gamma }_1{\tau }^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)}^2}\right]$

$\frac{{\partial }^{2}l}{\partial {\gamma }_1{\gamma }_2}=I(-\frac{1}{c}-1)\left[\frac{\left(c{t}^{{\gamma }_2}\left(\ln t\right)e^{{\beta }^{t}z}\right)(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)-\left(c{t}^{{\gamma }_2}{e}^{{\beta }^{t}z}\right)\left(c{\gamma }_1{t}^{{\gamma }_2}\left(\ln t\right)e^{{\beta }^{t}z}\right)}{{(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)}^2}\right]$

$+(1-I)(-\frac{1}{c})\left[\frac{\left(c{\tau }^{{\gamma }_2}\left(\ln \tau \right)e^{{\beta }^{t}z}\right)(c{\gamma }_1{\tau }^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)-\left(c{\tau }^{{\gamma }_2}{e}^{{\beta }^{t}z}\right)\left(c{\gamma }_1{\tau }^{{\gamma }_2}\left(\ln \tau \right)e^{{\beta }^{t}z}\right)}{{(c{\gamma }_1{\tau }^{{\gamma }_2}{e}^{{\beta }^{t}z}+1)}^2}\right]$

计算二阶偏导数的负期望值，

$E[-\frac{{\partial }^{2}l}{{{\partial \beta }_1}^2}]={\int }_0^{\tau }-\frac{{\partial }^{2}l}{{{\partial \beta }_1}^2}f(t;z)\mathrm{dt}$

.

.

.

$E[-\frac{{\partial }^{2}l}{{{\partial \beta }_k}^2}]={\int }_0^{\tau }-\frac{{\partial }^{2}l}{{{\partial \beta }_k}^2}f(t;z)\mathrm{dt}$

$E[-\frac{{\partial }^{2}l}{{\partial c}^2}]={\int }_0^{\tau }-\frac{{\partial }^{2}l}{{\partial c}^2}f(t;z)\mathrm{dt}$

$E[-\frac{{\partial }^{2}l}{{{\partial \gamma }_1}^2}]={\int }_0^{\tau }-\frac{{\partial }^{2}l}{{{\partial \gamma }_1}^2}f(t;z)\mathrm{dt}$

$E[-\frac{{\partial }^{2}l}{{{\partial \gamma }_2}^2}]={\int }_0^{\tau }-\frac{{\partial }^{2}l}{{{\partial \gamma }_2}^2}f(t;z)\mathrm{dt}$

$E[-\frac{{\partial }^{2}l}{\partial {\gamma }_1{\gamma }_2}]={\int }_0^{\tau }-\frac{{\partial }^{2}l}{\partial {\gamma }_1{\gamma }_2}f(t;z)\mathrm{dt}$

因此，Fisher信息矩阵就可以描述为，

$F=\underset{i_1=1}{\overset{q}{\Sigma }}\underset{i_2=1}{\overset{q}{\Sigma }}···\underset{i_k=1}{\overset{q}{\Sigma }}{\mathrm{np}}_{i_1{i}_2...i_k}{F}_{i_1{i}_2...i_k}---\left(5\right)$

对Fisher信息矩阵求逆，可以得到方差-协方差矩阵，

带入模型参数评估值渐进方差-协方差矩阵可以描述为，

这里所述的模型参数评估值是在加速寿命试验之前通过预试验失效数据或者相似产品历史失效数据得到的模型参数评估值。

步骤三、采用基于信息的优化方法建立最优化问题目标函数和约束条件；

利用预试验失效数据或者相似产品历史失效数据，计算得到渐进方差-协方差矩阵，从而得出该试验条件下的模型参数相关度(即协方差)。当模型参数相关度大于1时，说明模型参数具有较强的相关性，那么，采用D-优化方法建立最优化问题目标函数和约束条件。

D-优化方法最大化了信息矩阵的行列式值，这里即Fisher信息矩阵(5)。因此，目标函数及约束条件可以写为，

max det[F]

s.t.$0<p_{i_1{i}_2...i_k}<1,$i_j＝1，2，…q，j＝1，2，…k，

$\underset{i_1=1}{\overset{q}{\Sigma }}\underset{i_2=1}{\overset{q}{\Sigma }}···\underset{i_k=1}{\overset{q}{\Sigma }}{p}_{i_1{i}_2...i_k}=1,$

.

.

.

${\mathrm{np}}_{i_1{i}_2...i_k}[1-R({\tau }_{i_1{i}_2...i_k};z_{i_1{i}_2...i_k})]\ge {\mathrm{MNF}}_{i_1{i}_2...i_k}.$

其中，是应力水平下预先确定的所需最小产品失效数。

当模型参数相关度小于1时，说明模型参数的相关性较弱，那么，采用A-优化方法建立最优化问题目标函数和约束条件。

A-优化方法最小化了方差-协方差矩阵的迹，这里即渐进方差-协方差矩阵(6)。因此，目标函数及约束条件可以写为，

min  tr

s.t.$0<p_{i_1{i}_2...i_k}<1,$i_j＝1，2，…q，j＝1，2，…k，

$\underset{i_1=1}{\overset{q}{\Sigma }}\underset{i_2=1}{\overset{q}{\Sigma }}···\underset{i_k=1}{\overset{q}{\Sigma }}{p}_{i_1{i}_2...i_k}=1,$

.

.

.

${\mathrm{np}}_{i_1{i}_2...i_k}[1-R({\tau }_{i_1{i}_2...i_k};z_{i_1{i}_2...i_k})]\ge {\mathrm{MNF}}_{i_1{i}_2...i_k}.$

其中，是应力水平下预先确定的所需最小产品失效数。

步骤四、求解最优化问题，得到加速寿命试验优化设计方案；

利用模式搜索法等最优化方法求解最优化问题，可以得到满足约束的最优解，从而给出优化的恒定应力加速寿命试验方案。

基于比例危害-比例优势模型的步进应力加速寿命试验优化设计具体方法如下：

假设n个受试产品进行步进应力加速寿命试验，应力水平为z₁和z₂。对于每一个应力水平，有k种不同应力类型，即z₁＝(z₁₁，z₂₁，...z_k1)^t，z₂＝(z₁₂，z₂₂，...z_k2)^t。该试验在低应力水平z₁下运行到τ₁时间，转换到预先确定的高应力水平z₂运行至预先确定的τ₂时间截尾。步进应力加速寿命试验剖面如图2所示，k种不同步进应力类型分别有介于设计应力和工作极限之间的高、低两个应力水平，每个应力类型在低应力水平进行到τ₁时间时同时转换到高应力水平运行至预先确定的τ₂时间截尾。优化目标是合理的选取z₁和τ₁，使得在正常应力水平下对产品寿命和可靠性的评估结果最精确。具体的优化方法如下：

步骤一、建立基于比例危害-比例优势模型的对数似然函数；

采用累积损伤模型建立步进应力加速寿命试验中的对数似然函数。累积损伤模型假设产品剩余寿命只与其损伤量有关，而与损伤量累积方式无关。

考虑试验在低应力水平z₁下运行到τ₁时间，然后在高应力水平z₂下运行到τ₂时间截尾。对于应力水平z₂来说，它有一个等效的起始时间s。也就是说，受试产品在应力水平z₁下运行τ₁时间的损伤量相当于受试产品在应力水平z₂下运行s时间的损伤量。依据累积损伤模型，s可以由下式计算，

F(τ₁；z₁)＝F(s；z₂)

式中：F(·)为累积分布函数。

也就是，

$1-e^{-\frac{1}{c}\ln (c{\gamma }_1{{\tau }_1}^{{\gamma }_2}{e}^{{\beta }^t{z}_1}+1)}=1-e^{-\frac{1}{c}\ln (c{\gamma }_1{s}^{{\gamma }_2}{e}^{{\beta }^t{z}_2}+1)}$

因此，受试产品的累积分布函数可以表达为，

$F(t;z)=\left(\begin{array}{cc}F(t;z_1)& t\le {\tau }_1\\ F(t-{\tau }_1+s;z_2)& {\tau }_1<t\le {\tau }_2\end{array}\right)$

相应的，受试产品的概率密度函数为，

$f(t;z)=\left(\begin{array}{cc}f(t;z_1)& t\le {\tau }_1\\ f(t-{\tau }_1+s;z_2)& {\tau }_1<t\le {\tau }_2\end{array}\right)$

定义两个指示函数I₁和I₂，

其中，τ₁≤τ₂。

因此，对于某一个失效数据(t，I₁，I₂，z)，对数似然函数可以表达为，

l(t；z₁，z₂)＝I₂[I₁lnf(t；z₁)+(1-I₁)lnf(t′；z₂)]+(1-I₂)lnR(τ₂′；z₂)    (7)

式中，t′＝t-τ₁+s，τ₂′＝τ₂-τ₁+s。

将(1)和(2)带入(7)，可以将对数似然函数改写为，

$l(t;z_1,z_2)=I_1{I}_2[(-\frac{1}{c}-1)\ln (c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^t{z}_1}+1)+\ln {\gamma }_1+{\ln \gamma }_2+({\gamma }_2-1)\ln t+{\beta }^t{z}_1]$

$+(1-I_1)I_2[(-\frac{1}{c}-1)\ln (c{\gamma }_1{t}^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}+1)+{\ln \gamma }_1+{\ln \gamma }_2+({\gamma }_2-1){\ln t}^{\prime }+{\beta }^t{z}_2]---\left(8\right)$

$+(1-I_2)\left[(-\frac{1}{c})\ln (c{\gamma }_1{{\tau }_2}^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}+1)\right]$

步骤二、建立基于比例危害-比例优势模型的Fisher信息矩阵及渐进方差-协方差矩阵；

计算(8)式给出的对数似然函数对每个模型参数的二阶偏导数，

$\frac{{\partial }^{2}l}{{{\partial \beta }_1}^2}=I_1{I}_2\left[(-\frac{1}{c}-1)\frac{\left(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^t{z}_1}{z_{11}}^2\right)(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^t{z}_1}+1)-{\left(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^t{z}_1}{z}_{11}\right)}^2}{{(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^t{z}_1}+1)}^2}\right]$

$+(1-I_1)I_2\left[(-\frac{1}{c}-1)\frac{\left(c{\gamma }_1{t}^{{\prime \gamma }_2}{e}^{{\beta }^t{z}_2}{z_{12}}^2\right)(c{\gamma }_1{t}^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}+1)-{\left(c{\gamma }_1{t}^{{\prime \gamma }_2}{e}^{{\beta }^t{z}_2}{z}_{12}\right)}^2}{{(c{\gamma }_1{t}^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}+1)}^2}\right]$

$+(1-I_2)\left[(-\frac{1}{c})\frac{\left(c{\gamma }_1{\tau }_2^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}{z_{12}}^2\right)(c{\gamma }_1{\tau }_2^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}+1)-{\left(c{\gamma }_1{\tau }_2^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}{z}_{12}\right)}^2}{{(c{\gamma }_1{\tau }_2^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}+1)}^2}\right]$

.

.

.

$\frac{{\partial }^{2}l}{{{\partial \beta }_k}^2}=I_1{I}_2\left[(-\frac{1}{c}-1)\frac{\left(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^t{z}_1}{z}_{k1}^2\right)(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^t{z}_1}+1)-{\left(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^t{z}_1}{z}_{k1}\right)}^2}{{(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^t{z}_1}+1)}^2}\right]$

$+(1-I_1)I_2\left[(-\frac{1}{c}-1)\frac{\left(c{\gamma }_1{t}^{{\prime \gamma }_2}{e}^{{\beta }^t{z}_2}{z}_{k2}^2\right)(c{\gamma }_1{t}^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}+1)-{\left(c{\gamma }_1{t}^{{\prime \gamma }_2}{e}^{{\beta }^t{z}_2}{z}_{k2}\right)}^2}{{(c{\gamma }_1{t}^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}+1)}^2}\right]$

$+(1-I_2)\left[(-\frac{1}{c})\frac{\left(c{\gamma }_1{\tau }_2^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}{z}_{k2}^2\right)(c{\gamma }_1{\tau }_2^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}+1)-{\left(c{\gamma }_1{\tau }_2^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}{z}_{k2}\right)}^2}{{(c{\gamma }_1{\tau }_2^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}+1)}^2}\right]$

$\frac{{\partial }^{2}l}{{\partial c}^2}=I_1{I}_2[(-\frac{2}{c^3})\ln (c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^t{z}_1}+1)+\left(\frac{2}{c^2}\right)\frac{{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^t{z}_1}}{c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^t{z}_1}+1}+(\frac{1}{c}+1)\frac{{\left({\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^t{z}_1}\right)}^2}{{(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^t{z}_1}+1)}^2}]$

$+(1-I_1)I_2[(-\frac{2}{c^3})\ln (c{\gamma }_1{t}^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}+1)+\left(\frac{2}{c^2}\right)\frac{{\gamma }_1{t}^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}}{c{\gamma }_1{t}^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}+1}+(\frac{1}{c}+1)\frac{{\left({\gamma }_1{t}^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}\right)}^2}{{(c{\gamma }_1{t}^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}+1)}^2}]$

$+(1-I_2)[(-\frac{2}{c^3})\ln (c{\gamma }_1{\tau }_2^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}+1)+\left(\frac{2}{c^2}\right)\frac{{\gamma }_1{\tau }_2^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}}{c{\gamma }_1{\tau }_2^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}+1}+\left(\frac{1}{c}\right)\frac{{\left({\gamma }_1{\tau }_2^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}\right)}^2}{{(c{\gamma }_1{\tau }_2^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}+1)}^2}]$

$\frac{{\partial }^{2}l}{{{\partial \gamma }_1}^2}=I_1{I}_2[(\frac{1}{c}+1)\frac{{\left(c{t}^{{\gamma }_2}{e}^{{\beta }^t{z}_1}\right)}^2}{{(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^t{z}_1}+1)}^2}-\frac{1}{{{\gamma }_1}^2}]+(1-I_1)I_2[(\frac{1}{c}+1)\frac{{\left({\mathrm{ct}}^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}\right)}^2}{{(c{\gamma }_1{t}^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}+1)}^2}-\frac{1}{{{\gamma }_1}^2}]$

$+(1-I_2)\left[\left(\frac{1}{c}\right)\frac{{\left(c{\tau }_2^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}\right)}^2}{{(c{\gamma }_1{\tau }_2^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}+1)}^2}\right]$

$\frac{{\partial }^{2}l}{\partial {{\gamma }_2}^2}=I_1{I}_2[(-\frac{1}{c}-1)\frac{\left(c{\gamma }_1{t}^{{\gamma }_2}{\left(\ln t\right)}^2{e}^{{\beta }^t{z}_1}\right)(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^t{z}_1}+1)-{\left(c{\gamma }_1{t}^{{\gamma }_2}\left(\ln t\right)e^{{\beta }^t{z}_1}\right)}^2}{{(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^t{z}_1}+1)}^2}-\frac{1}{{{\gamma }_2}^2}]$

$+(1-I_1)I_2[(-\frac{1}{c}-1)\frac{\left(c{\gamma }_1{t}^{\prime {\gamma }_2}{\left({\ln t}^{\prime }\right)}^2{e}^{{\beta }^t{z}_2}\right)(c{\gamma }_1{t}^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}+1)-{\left(c{\gamma }_1{t}^{\prime {\gamma }_2}\left({\ln t}^{\prime }\right)e^{{\beta }^t{z}_2}\right)}^2}{{(c{\gamma }_1{t}^{\prime {\gamma }_2}{e}^{{\beta }^{\prime }{z}_2}+1)}^2}-\frac{1}{{{\gamma }_2}^2}]$

$+(1-I_2)\left[(-\frac{1}{c})\frac{\left(c{\gamma }_1{\tau }_2^{\prime {\gamma }_2}{\left({\ln \tau }_2^{\prime }\right)}^2{e}^{{\beta }^t{z}_2}\right)(c{\gamma }_1{\tau }_2^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}+1)-{\left(c{\gamma }_1{\tau }_2^{\prime {\gamma }_2}\left({\ln \tau }_2^{\prime }\right)e^{{\beta }^t{z}_2}\right)}^2}{{(c{\gamma }_1{\tau }_2^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}+1)}^2}\right]$

$\frac{{\partial }^{2}l}{\partial {\gamma }_1{\gamma }_2}=I_1{I}_2\left[(-\frac{1}{c}-1)\frac{\left(c{t}^{{\gamma }_2}\left(\ln t\right)e^{{\beta }^t{z}_1}\right)(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^t{z}_1}+1)-\left({\mathrm{ct}}^{{\gamma }_2}{e}^{{\beta }^t{z}_1}\right)\left(c{\gamma }_1{t}^{{\gamma }_2}\left(\ln t\right)e^{{\beta }^t{z}_1}\right)}{{(c{\gamma }_1{t}^{{\gamma }_2}{e}^{{\beta }^t{z}_1}+1)}^2}\right]$

$+(1-I_1)I_2\left[(-\frac{1}{c}-1)\frac{\left(c{t}^{\prime {\gamma }_2}\left({\ln t}^{\prime }\right)e^{{\beta }^t{z}_2}\right)(c{\gamma }_1{t}^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}+1)-\left({\mathrm{ct}}^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}\right)\left(c{\gamma }_1{t}^{\prime {\gamma }_2}\left({\ln t}^{\prime }\right)e^{{\beta }^t{z}_2}\right)}{{(c{\gamma }_1{t}^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}+1)}^2}\right]$

$+(1-I_2)\left[(-\frac{1}{c})\frac{\left(c{\tau }_2^{\prime {\gamma }_2}\left({\ln \tau }_2^{\prime }\right)e^{{\beta }^t{z}_2}\right)(c{\gamma }_1{\tau }_2^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}+1)-\left(c{\tau }_2^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}\right)\left(c{\gamma }_1{\tau }_2^{\prime {\gamma }_2}\left({\ln \tau }_2^{\prime }\right)e^{{\beta }^t{z}_2}\right)}{{(c{\gamma }_1{\tau }_2^{\prime {\gamma }_2}{e}^{{\beta }^t{z}_2}+1)}^2}\right]$

计算二阶偏导数的负期望值，

$E[-\frac{{\partial }^{2}l}{{{\partial \beta }_1}^2}]={\int }_0^{{\tau }_2}-\frac{{\partial }^{2}l}{{{\partial \beta }_1}^2}f(t;z)\mathrm{dt}$

.

.

.

$E[-\frac{{\partial }^{2}l}{{{\partial \beta }_k}^2}]={\int }_0^{{\tau }_2}-\frac{{\partial }^{2}l}{{{\partial \beta }_k}^2}f(t;z)\mathrm{dt}$

$E[-\frac{{\partial }^{2}l}{\partial c^2}]={\int }_0^{{\tau }_2}-\frac{{\partial }^{2}l}{{\partial c}^2}f(t;z)\mathrm{dt}$

$E[-\frac{{\partial }^{2}l}{{{\partial \gamma }_1}^2}]={\int }_0^{{\tau }_2}-\frac{{\partial }^{2}l}{{{\partial \gamma }_1}^2}f(t;z)\mathrm{dt}$

$E[-\frac{{\partial }^{2}l}{{{\partial \gamma }_2}^2}]={\int }_0^{{\tau }_2}-\frac{{\partial }^{2}l}{{{\partial \gamma }_2}^2}f(t;z)\mathrm{dt}$

$E[-\frac{{\partial }^{2}l}{\partial {\gamma }_1{\gamma }_2}]={\int }_0^{{\tau }_2}-\frac{{\partial }^{2}l}{\partial {\gamma }_1{\gamma }_2}f(t;z)\mathrm{dt}$

因此，Fisher信息矩阵就可以描述为，

F＝nF_i                           (9)

对Fisher信息矩阵求逆，可以得到方差-协方差矩阵，

带入模型参数评估值渐进方差-协方差矩阵可以描述为，

这里所述的模型参数评估值是在加速寿命试验之前通过预试验失效数据或者相似产品历史失效数据得到的模型参数评估值。

步骤三、采用基于信息的优化方法建立最优化问题目标函数和约束条件；

利用预试验失效数据或者相似产品历史失效数据，计算得到渐进方差-协方差矩阵，从而得出该试验条件下的模型参数相关度(即协方差)。当模型参数相关度大于1时，说明模型参数具有较强的相关性，那么，采用D-优化方法建立最优化问题目标函数和约束条件。

D-优化方法最大化了信息矩阵的行列式值，这里即Fisher信息矩阵(9)。因此，目标函数及约束条件可以写为，

max det[F]

s.t.    z₁₁≤z₁₂≤z_{1，设计应力}，

        z₂₁≤z₂₂≤z_{2，设计应力}，

     z_k1≤z_k2≤z_{k，设计应力}，

     n[1-R(τ₁；z₁)]≥MNF₁，

     0≤τ₁≤τ₂.

其中，MNF₁是应力水平z₁下预先确定的所需最小产品失效数。

当模型参数相关度小于1时，说明模型参数的相关性较弱，那么，采用A-优化方法建立最优化问题目标函数和约束条件。

A-优化方法最小化了方差-协方差矩阵的迹，这里即渐进方差-协方差矩阵(10)。因此，目标函数及约束条件可以写为，

min  tr

s.t. z₁₁≤z₁₂≤z_{1，设计应力}，

     z₂₁≤z₂₂≤z_{2，设计应力}，

            .

            .

            .

     z_k1≤z_k2≤z_{k，设计应力}，

     n[1-R(τ₁；z₁)]≥MNF₁，

     0≤τ₁≤τ₂.

其中，MNF₁是应力水平z₁下预先确定的所需最小产品失效数。

步骤四、求解最优化问题，得到加速寿命试验优化设计方案；

利用模式搜索法等最优化方法求解最优化问题，可以得到满足约束的最优解，从而给出优化的步进应力加速寿命试验方案。

实施例1：恒定应力加速寿命试验优化设计；

假设一产品进行三水平恒定应力加速寿命试验，试验应力为温度应力。可获得的受试产品数量为300个。受试产品的工作极限为200℃。对于每个应力水平，截尾时间为300小时，并且要求每个应力水平下至少有20个受试产品失效。

温度应力水平值依照阿伦尼斯模型转换为100K^-1，因此，正常应力水平25℃以及产品工作极限应力水平200℃分别转换为0.336和0.211。采用比例危害-比例优势模型对该产品加速寿命试验预试验的失效数据进行评估得到的模型参数评估值为：β＝-60，c＝0.5，α＝2，γ＝3。这些模型参数评估值将运用于加速寿命试验优化设计中。按照预试验得出的模型参数评估值，针对正常应力水平下300个受试产品运行300小时的情况，计算得到模型参数的相关度为6.40，因此，认为模型参数具有较强的相关性。

通过如下步骤来进行该产品恒定应力加速寿命试验优化设计：

(1)建立如(4)式的基于比例危害-比例优势模型的对数似然函数；

(2)建立如(5)和(6)式的基于比例危害-比例优势模型的Fisher信息矩阵和渐进方差-协方差矩阵；

(3)由于模型参数具有较强的相关性，因此，本例采用D-优化方法。目标函数是(5)式的Fisher信息矩阵的行列式值，决策变量是和i_j＝1，2，…q，j＝1，2，…k；

对于该问题，建立目标函数和约束条件如下，

max det[F]

s.t.  0＜p_i＜1，i＝1，2，3，

$\underset{i=1}{\overset{3}{\Sigma }}{p}_i=1,$

      0.211≤z₁≤z₂≤z₃≤0.336，

      np_i[1-R(τ；z_i)]≥MNF_i，i＝1，2，3，

其中，n＝300，

      MNF_i＝20，i＝1，2，3，

      τ＝300.

(4)最大化目标函数可以得到满足约束条件的决策变量最优解，从而给出优化的恒定应力加速寿命试验方案。

满足约束条件的优化结果为z₁＝0.233，z₂＝0.253，z₃＝0.265，p₁＝0.241，p₂＝0.311，p₃＝0.449。因此，该产品恒定应力加速寿命试验最优试验方案如表1所示，

表1恒定应力加速寿命试验最优试验方案

  温度  受试产品个数  104℃  135  122℃  93  156℃  72

目标函数值为3230787.265。

实施例2：步进应力加速寿命试验优化设计；

假设一产品进行两水平步进应力加速寿命试验，试验应力为温度应力。可获得的受试产品数量为300个。受试产品的工作极限为210℃。试验在低应力水平z₁进行τ₁时间后，转到高应力水平z₂(取略低于产品工作极限的应力水平值200℃)，并在τ₂＝300小时截尾。要求低应力水平z₁下至少有50个受试产品失效。

温度应力水平值依照阿伦尼斯模型转换为100K^-1，因此，正常应力水平25℃以及高应力水平200℃分别转换为0.336和0.211。采用比例危害-比例优势模型对该产品加速寿命试验预试验的失效数据进行评估得到的模型参数评估值为：β＝-60，c＝0.5，α＝2，γ＝3。这些模型参数评估值将运用于加速寿命试验优化设计中。按照预试验得出的模型参数评估值，针对正常应力水平下300个受试产品运行300小时的情况，计算得到模型参数的相关度为6.40，因此，认为模型参数具有较强的相关性。

通过如下步骤来进行该产品步进应力加速寿命试验优化设计：

(1)建立如(8)式的基于比例危害-比例优势模型的对数似然函数；

(2)建立如(9)和(10)式的基于比例危害-比例优势模型的Fisher信息矩阵和渐进方差-协方差矩阵；

(3)由于模型参数具有较强的相关性，因此，本文采用D-优化方法。目标函数是(9)式的Fisher信息矩阵的行列式值，决策变量是z₁和τ₁。

对于该问题，建立目标函数和约束条件如下，

max det[F]

s.t.  0.211≤z₁≤0.336，

      n[1-R(τ₁；z₁)]≥MNF₁，

      0≤τ₁≤τ₂，

其中，n＝300，

      MNF₁＝50，

      τ₂＝300.

(4)最大化目标函数可以得到满足约束条件的决策变量最优解，从而给出优化的步进应力加速寿命试验方案。

满足约束条件的优化结果为z₁＝0.292，τ₁＝176.5。因此，该产品步进应力加速寿命试验最优试验方案如表2所示，

表2步进应力加速寿命试验最优试验方案

  温度  应力转换时间/  截尾时间(小时)  69℃  176.5  200℃  300

目标函数值为10272030239785200000。