CN104778354B - 定期抽检型产品贮存寿命评估方法 - Google Patents

Abstract

一种定期抽检型产品贮存寿命评估方法，方法步骤如下：一、用PAVA算法修正原始试验数据中的倒挂；二、求分布函数参数的极小卡方估计，计算皮尔逊卡方统计量；三、进行皮尔逊卡方拟合优度检验；四、计算产品在给定贮存条件下的可靠寿命。本发明的优点是：它保证了寿命分布参数估计的准确性和完整性，其算法对参数的初值要求较低，算法迭代快速简单，可操作性强，较之极大似然估计对样本量的要求更低，评估结果更加稳定，体现了其在处理小样本数据时的优越性。

Description

定期抽检型产品贮存寿命评估方法

技术领域

本发明涉及一种基于保序回归理论和皮尔逊卡方拟合优度理论的定期抽检型产品贮存寿命评估方法。它针对定期抽检型产品贮存寿命试验数据，通过保序回归理论中的相邻逆序量合并算法(The Pool-Adjacent-Violators Algorithm，简称PAVA算法)进行调整，对产品总体累积失效函数的分布参数进行极小卡方估计，并运用拟合优度理论检验分布函数假设的合理性，最终评估产品的贮存可靠寿命。适用于定期抽检型产品贮存寿命试验中小样本评估等领域。

背景技术

通常情况下，军用产品在批量生产之后需要在仓库中贮存一段时间，称为产品的贮存期。在贮存期内，由于贮存环境的作用，产品可靠性指标会有所下降。为了鉴定贮存期不同阶段产品质量的变化，需要在贮存过程中对产品进行检验并分析试验数据。

通常的定期抽检型产品贮存寿命试验按如下方法进行：

记t＝(t₁,…,t_k)为检测时间点，在贮存了时间t_i之后，抽取n_i件产品进行试验，发现X_i件产品失效，从而得到数据：

(t_i,n_i,X_i),i＝1,2,…,k (1)

其中，0＜t₁＜t₂＜…＜t_k，各次抽检试验相互独立。

在实际工作中，受时间和费用等约束，试验样本量n_i不能取得太大值，因此试验结果受抽样的随机性影响就会很大，有时会出现X_i+1/n_i+1＜＜X_i/n_i的情况，这与样本失效概率随时间增加而降低的特性不相符合，工程上称之为数据的倒挂，严重的倒挂会使计算结果产生很大偏差。另一方面，针对定期抽检型产品贮存寿命试验得到的不完全数据，在应用极大似然估计等传统评估方法时，没有进行任何调整和删失，结果给出的是产品累积失效函数分布参数及可靠寿命的保守估计。在研究中发现，基于皮尔逊卡方拟合优度理论的极小卡方估计在评估此类数据时效果明显优于极大似然估计，通过仿真试验可知，极小卡方估计在评估小样本数据时比极大似然估计具有更小的变异系数。

基于此本发明提出一种基于保序回归理论和皮尔逊卡方拟合优度理论的定期抽检型产品贮存寿命评估方法，综合考虑试验数据类型和评估方法，给出产品累积失效函数的分布参数及可靠寿命的评估方法。

发明内容

(1)本发明的目的：针对定期抽检型产品在贮存寿命试验中失效数少而产生的一系列估计问题，提供一种包含倒挂数据处理、寿命估计、拟合优度检验的完整的产品贮存寿命评估方法。通过保序回归理论中PAVA算法调整试验数据，对产品累积失效函数的分布参数进行极小卡方估计，并运用拟合优度理论检验分布函数假设的合理性，最终评估产品在给定贮存条件下的可靠寿命。

(2)技术方案：

本发明需建立如下基本设置：

设置1在进行产品贮存寿命试验时，不同批次的样本需来自同一总体，产品之间可靠性指标无显著差异，累积失效函数相同。

设置2产品贮存寿命t服从指数分布、威布尔分布、I型极大值分布和II型极大值分布中的一种，各分布的累积失效函数分别为：

①指数分布：

F(t)＝1-exp(-t/θ) (2)

②威布尔分布：

F(t)＝1-exp(-(t/η)^m) (3)

③I型极大值分布：

F(t)＝exp(-exp(-(t-μ)/σ)) (4)

④II型极大值分布：

F(t)＝exp(-(t/η1)^-m1) (5)

其中，θ为指数分布的均值寿命，η和m分别为威布尔分布的尺度参数和形状参数，μ和σ分别为I型极大值分布的位置参数和尺度参数，η1和m1分别为II型极大值分布的尺度参数和形状参数。

本发明提出的方法主要针对定期抽检型产品贮存寿命试验数据，根据保序回归理论中PAVA算法调整原始失效频率，并运用极小卡方估计的方法，对产品的可靠性参数进行评估，进而对产品总体累积失效函数进行假设检验。

基于上述假设与思路，本发明提供的定期抽检型产品贮存寿命评估方法通过如下步骤实现：

步骤一：用PAVA算法修正原始试验数据中的倒挂。

根据定期抽检型产品贮存寿命试验数据结构，计算t_i时刻样本的原始失效频率f_i：

f_i＝X_i/n_i,i＝1,2,…,k (6)

其中，n_i为样本数，X_i为样本中检测到的故障数。

用PAVA算法调整样本原始失效频率中的倒挂，使其成为满足序约束的失效频率值，即求出{f_i＝X_i/n_i}关于权值的保序回归；

将样本的保序失效频率值记为计算相应的样本保序失效数

步骤二：求分布函数参数的极小卡方估计，计算皮尔逊卡方统计量。

在拟合优度检验的相关理论中，皮尔逊提出了皮尔逊卡方统计量用于检验一组独立样本的共同分布是否属于某一具有特定性质的分布族，其中变样本复合皮尔逊卡方统计量的形式为：

其中，代表参数向量，s为参数个数，为t_i时刻的理论失效频率。

具体步骤为：

I.假设总体在时刻t_i的理论失效频率服从特定分布

II.极小化皮尔逊卡方统计量即Λ表示参数集合，将所得到的参数作为真值的最佳估计，称之为的极小卡方估计；

其中，求解等价于解方程组：

其中，表示皮尔逊卡方统计量对各个参数求偏导数所得到的结果，表示理论失效分布函数对各个参数求偏导数所得到的结果。

故的极小卡方估计是下面方程的解：

III.计算皮尔逊卡方统计量。

由于参数的真值未知，故将的极小卡方估计值代入方程(8)中，可得到近似皮尔逊卡方统计量

步骤三：进行皮尔逊卡方拟合优度检验。

皮尔逊卡方统计量描述了期望频数与观察频数之间的差异。当n_i→∞时，的极限分布是自由度为k-1的卡方分布即

在实际检验过程中，由于未知，故用近似皮尔逊卡方统计量代替 的极限分布是自由度为k-s-1的卡方分布s为参数个数。因此若给定显著性水平α，则有：

查表得到的值，比较与的大小。

当成立时，小概率事件发生，拒绝原假设；

当成立时，分布假设满足要求，计算分布对应的p值。

步骤四：计算产品在给定贮存条件下的可靠寿命。

选择满足假设检验的分布中p值最大的作为产品总体累积失效函数，记为给定可靠度R，计算对应的可靠寿命t_R：

其中，F^-1(·)为总体累积失效函数的反函数。

其中，在步骤一中所述的“保序回归”，是指：

令T＝{t₁,t₂,…,t_k}为一个有限集合，f＝(f₁,f₂,…,f_k)′为定义在T上的有界函数。若T上又定义了一种半序关系“＜”，且对t_i∈T,t_j∈T,t_i＜t_j，均有：

成立，则称函数为定义在T上相对于“＜”的保序函数。

记G为保序函数的全体，若存在f^*∈G，满足：

则称为f的保序回归，其中ω＝(ω₁,ω₁,…,ω_k)′,ω_i＞0是给定权函数。

在步骤一中所述的“PAVA算法”，是指：

I.如果f∈G，则f^*＝f；

II.如果存在j使得f_j＞f_j+1，则选取适当权重ω＝(ω₁,ω₁,…,ω_k)′,ω_i＞0，令：

B＝{j,j+1} (15)

f_B＝A_V(B)＝(Σ_i∈Bf_iω_i)/(Σ_i∈Bω_i) (16)

ω_B＝ω_j+ω_j+1 (17)

并且令：

III.重复上述步骤II，直到把下标集k分解为l个小块B₁,B₂,…,B_l，并且满足条件：

A_V(B₁)＜A_V(B₂)＜…＜A_V(B_l)

则：

(3)优点和功效：本发明是一种基于保序回归理论和皮尔逊卡方拟合优度理论的定期抽检型产品贮存寿命评估方法，其优点是：

①本发明针对定期抽检型产品贮存寿命试验数据，根据保序回归理论中PAVA算法调整原始失效频率，并运用极小卡方估计的方法，对产品的可靠性参数进行评估，进而对产品总体累积失效函数进行假设检验，保证了寿命分布参数估计的准确性和完整性。

②本发明的算法对参数的初值要求较低，算法迭代快速简单，可操作性强。

③本发明较之极大似然估计对样本量的要求更低，评估结果更加稳定，体现了其在处理小样本数据时的优越性。

附图说明

图1本发明所述方法流程图。

具体实施方式

见图1,下面将结合实施例对本发明做进一步详细说明。

某复杂机电产品在平均温度为21.4℃的恒湿环境下贮存，在不同时间点抽取不同批次的样本进行检验，经分析，不同批次间产品贮存失效服从同一分布函数。统计贮存时间和失效数，综合近10年的检测结果可得试验数据，如表1所示：

表1某复杂机电产品贮存试验数据

序号	贮存时间/m	检测样本数	失效样本数
				1	30	1157	0
2	42	340	0
				3	47	875	2
4	48	484	3
				5	54	1438	0
6	59	502	1
				7	66	720	0
8	68	446	4
				9	74	1536	1
10	78	949	0
				11	84	496	0
12	88	974	3
				13	90	356	0
14	91	76	1
				15	96	1536	0
16	97	556	2
				17	102	211	0
18	113	971	9

注：m表示月。

假设该复杂机电产品的贮存寿命服从指数分布、威布尔分布、I型极大值分布和II型极大值分布中的一种，根据本说明提出了定期抽检型产品贮存寿命评估方法进行产品总体累积失效函数选择以及参数评估，并预测可靠贮存寿命。

步骤一：用PAVA算法修正原始试验数据中的倒挂。

首先计算样本的原始失效频率：f_i＝X_i/n_i；

再通过PAVA算法，计算保序回归后的失效频率值f_i ^*，其中权重：ω_i＝f_i(1-f_i/n_i)。

计算结果如表2所示：

表2样本失效频率处理前后对比

贮存时间/m

检测样本数

失效数

失效频率

保序后频率

30	1157	0	0	0
					42	340	0	0	0
47	875	2	0.002286	0.002286
					48	484	3	0.006198	0.005468
54	1438	0	0	0.005468
					59	502	1	0.001992	0.005468
66	720	0	0	0.005468
					68	446	4	0.008969	0.007758
74	1536	1	0.000651	0.007758
					78	949	0	0	0.007758
84	496	0	0	0.007758
					88	974	3	0.00308	0.007758
90	356	0	0	0.007758
					91	76	1	0.013158	0.01128
96	1536	0	0	0.01128
					97	556	2	0.003597	0.01128
102	211	0	0	0.01128
					113	971	9	0.009269	0.01128

步骤二：求分布函数参数的极小卡方估计，计算皮尔逊卡方统计量。

首先，分别假设总体在时刻t_i的理论失效频率服从指数分布、威布尔分布、I型极大值分布和II型极大值分布，通过求解方程(10)，可以得到分布参数的极小卡方估计值，结果如表3所示：

表3分布参数的极小卡方估计值

其次，将求得的分布参数值分别代入方程(11)，可以求得近似皮尔逊卡方统计量，结果如表4所示：

表4四种分布下的近似皮尔逊卡方统计量

步骤三：进行皮尔逊卡方拟合优度检验。

取显著性水平α＝0.05，查表得的值并与近似皮尔逊卡方统计量作比较，结果如表5所示：

表5假设检验结果

步骤四：计算产品在给定贮存条件下的可靠寿命。

由表5可知，四种分布函数均满足假设检验。其中，I型极大值分布的p值最大，故选取I型极大值分布作为该复杂机电产品贮存寿命分布函数。当可靠度R＝0.95时，由I型极大值分布函数可解得，该复杂机电产品的可靠寿命为t_0.95＝16.7年。

所以基于前面所述计算评估得产品的95％的可靠贮存寿命为16.7年。

结果表明，采用本发明方法可以精确评估产品的贮存寿命，达到预期的目的。

综上所述，本发明涉及一种基于保序回归理论和皮尔逊卡方拟合优度理论的定期抽检型产品贮存寿命评估方法。它针对定期抽检型产品贮存寿命试验，通过保序回归理论中PAVA算法调整试验数据，对产品总体累积失效函数的分布参数进行极小卡方估计，并运用皮尔逊卡方拟合优度理论检验分布函数假设的合理性，评估可靠寿命。该方法的具体步骤是：首先用PAVA算法修正原始试验数据中的倒挂，然后分别假设总体的累积失效函数服从指数分布、威布尔分布、I型极大值分布和II型极大值分布，求出分布参数的极小卡方估计并计算近似皮尔逊卡方统计量，进而开展皮尔逊卡方拟合优度检验，最后计算产品在给定贮存条件下的可靠寿命。本发明适用于定期抽检型产品贮存寿命试验中小样本评估等领域，具有较强的操作性。

Claims (1)

1.一种定期抽检型产品贮存寿命评估方法，其特征在于：设置如下：

设置1：在进行产品贮存寿命试验时，不同批次的样本需来自同一总体，产品之间可靠性指标无显著差异，累积失效函数相同；

设置2：产品贮存寿命t服从指数分布、威布尔分布、I型极大值分布和II型极大值分布中的一种，各分布的累积失效函数分别为：

①指数分布：

F(t)＝1-exp(-t/θ)…………………………………(1)

②威布尔分布：

F(t)＝1-exp(-(t/η)^m)………………………………(2)

③I型极大值分布：

F(t)＝exp(-exp(-(t-μ)/σ))………………………(3)

④II型极大值分布：

F(t)＝exp(-(t/η1)^-m1)………………………………(4)

其中，θ为指数分布的均值寿命，η和m分别为威布尔分布的尺度参数和形状参数，μ和σ分别为I型极大值分布的位置参数和尺度参数，η1和m1分别为II型极大值分布的尺度参数和形状参数；

根据保序回归理论中PAVA算法调整原始失效频率，并运用极小卡方估计的方法，对产品的可靠性参数进行评估，进而对产品总体累积失效函数进行假设检验；

定期抽检型产品贮存寿命评估方法通过如下步骤实现：

步骤一：用PAVA算法修正原始试验数据中的倒挂；

根据定期抽检型产品贮存寿命试验数据结构，计算t_i时刻样本的原始失效频率f_i：

f_i＝X_i/n_i,i＝1,2,…,k………………………………(5)

其中，n_i为样本数，X_i为样本中检测到的故障数；

用PAVA算法调整样本原始失效频率中的倒挂，使其成为满足序约束的失效频率值，即求出{f_i＝X_i/n_i}关于权值的保序回归；

将样本的保序失效频率值记为f_i ^*，计算相应的样本保序失效数

$X_i^{*}=f_i^{*}\·n_i,i=1,2,...,k...\left(6\right)$

步骤二：求分布函数参数的极小卡方估计，计算皮尔逊卡方统计量；

在拟合优度检验的相关理论中，皮尔逊提出了皮尔逊卡方统计量用于检验一组独立样本的共同分布是否属于某一具有特定性质的分布族，

其中，变样本复合皮尔逊卡方统计量的形式为：

${\mathrm{\χ}}^2\left(\stackrel{\→}{\mathrm{\λ}}\right)=\underset{i=1}{\overset{k}{\Σ}}\frac{{(X_i-n_i{p}_i(\stackrel{\→}{\mathrm{\λ}}\left)\right)}^2}{n_i{p}_i\left(\stackrel{\→}{\mathrm{\λ}}\right)}...\left(7\right)$

其中，代表参数向量，s为参数个数，为t_i时刻的理论失效频率；

具体步骤为：

I.假设总体在时刻t_i的理论失效频率服从特定分布

II.极小化皮尔逊卡方统计量即Λ表示参数集合，将所得到的参数作为真值的最佳估计，称之为的极小卡方估计；

其中，求解等价于解方程组：

$\begin{array}{c}\frac{\∂{\mathrm{\χ}}^2\left(\stackrel{\→}{\mathrm{\λ}}\right)}{\∂{\mathrm{\λ}}_j}=\underset{i=1}{\overset{k}{\Σ}}\{-\frac{2n_i(X_i-n_i{p}_i(\stackrel{\→}{\mathrm{\λ}}))}{n_i{p}_i\left(\stackrel{\→}{\mathrm{\λ}}\right)}-\frac{(X_i-n_i{p}_i{\left(\stackrel{\→}{\mathrm{\λ}}\right)}^2}{n_i{p}_i^2\left(\stackrel{\→}{\mathrm{\λ}}\right)}\}\·\frac{\∂p_i\left(\stackrel{\→}{\mathrm{\λ}}\right)}{\∂{\mathrm{\λ}}_j}\\ =\underset{i=1}{\overset{k}{\Σ}}\left\{\frac{-X_i^2+n_i^2{p}_i^2\left(\stackrel{\→}{\mathrm{\λ}}\right)}{n_i{p}_i^2\left(\stackrel{\→}{\mathrm{\λ}}\right)}\right\}\·\frac{\∂p_i\left(\stackrel{\→}{\mathrm{\λ}}\right)}{\∂{\mathrm{\λ}}_j}\\ =\underset{i=1}{\overset{k}{\Σ}}(1-{\left(\frac{f_i}{p_i\left(\stackrel{\→}{\mathrm{\λ}}\right)}\right)}^2)\·\frac{n_i\·\∂p_i\left(\stackrel{\→}{\mathrm{\λ}}\right)}{\∂{\mathrm{\λ}}_j}=0,j=1,2,\mathrm{...},s\end{array}...\left(8\right)$

其中，表示皮尔逊卡方统计量对各个参数求偏导数所得到的结果，表示理论失效分布函数对各个参数求偏导数所得到的结果；

故的极小卡方估计是下面方程的解：

$\underset{i=1}{\overset{k}{\Σ}}\{1-{\[\frac{f_i}{p_i\left(\stackrel{\→}{\mathrm{\λ}}\right)}\]}^2\}\·\frac{n_i\·\∂p_i\left(\stackrel{\→}{\mathrm{\λ}}\right)}{\∂{\mathrm{\λ}}_j}=0,j=1,2,...,s...\left(9\right)$

III.计算皮尔逊卡方统计量；

由于参数的真值未知，故将的极小卡方估计值代入方程(7)中，得到近似皮尔逊卡方统计量

${\mathrm{\χ}}^2\left({\stackrel{\→}{\mathrm{\λ}}}^{\′}\right)=\underset{i=1}{\overset{k}{\Σ}}\frac{{(X_i^{*}-n_i{p}_i({\stackrel{\→}{\mathrm{\λ}}}^{\′}\left)\right)}^2}{n_i{p}_i\left({\stackrel{\→}{\mathrm{\λ}}}^{\′}\right)}=\underset{i=1}{\overset{k}{\Σ}}\frac{n_i{(f_i^{*}-p_i({\stackrel{\→}{\mathrm{\λ}}}^{\′}\left)\right)}^2}{p_i\left({\stackrel{\→}{\mathrm{\λ}}}^{\′}\right)}...\left(10\right)$

步骤三：进行皮尔逊卡方拟合优度检验；

皮尔逊卡方统计量描述了期望频数与观察频数之间的差异；当n_i→∞时，的极限分布是自由度为k-1的卡方分布即

在实际检验过程中，由于未知，故用近似皮尔逊卡方统计量代替的极限分布是自由度为k-s-1的卡方分布s为参数个数；因此若给定显著性水平α，则有：

$P\left({\mathrm{\χ}}^2\right({\stackrel{\→}{\mathrm{\λ}}}^{\′})\≥{\mathrm{\χ}}_{k-s-1}^2(1-\mathrm{\α}\left)\right|H_0)\≤\mathrm{\α}...\left(11\right)$

查表得到的值，比较与的大小；

当成立时，小概率事件发生，拒绝原假设；

当成立时，分布假设满足要求，计算分布对应的p值；

步骤四：计算产品在给定贮存条件下的可靠寿命；

选择满足假设检验的分布中p值最大的作为产品总体累积失效函数，记为给定可靠度R，计算对应的可靠寿命t_R：

$t_R=F^{-1}(R,\stackrel{\→}{\mathrm{\λ}})$

其中，F^-1(·)为总体累积失效函数的反函数；

其中，在步骤一中所述的“保序回归”，是指：

令T＝{t₁,t₂,…,t_k}为一个有限集合，f＝(f₁,f₂,…,f_k)′为定义在T上的有界函数；若T上又定义了一种半序关系“＜”，且对t_i∈T,t_j∈T,t_i＜t_j，均有：

$f_i^{*}=f^{*}\left(t_i\right)\≤f_j^{*}=f^{*}\left(t_j\right)...\left(12\right)$

成立，则称函数为定义在T上相对于“＜”的保序函数；

记G为保序函数的全体，若存在f^*∈G，满足：

$\underset{i=1}{\overset{k}{\Σ}}{(f_i-f_i^{*})}^2{\mathrm{\ω}}_i=\underset{\∀g\∈G}{min}\underset{i=1}{\overset{k}{\Σ}}{(f_i-g_i)}^2{\mathrm{\ω}}_i...\left(13\right)$

则称为f的保序回归，其中ω＝(ω₁,ω₁,…,ω_k)′,ω_i＞0是给定权函数；

在步骤一中所述的“PAVA算法”，是指：

I.如果f∈G，则f^*＝f；

II.如果存在j使得f_j＞f_j+1，则选取适当权重ω＝(ω₁,ω₁,…,ω_k)′,ω_i＞0，令：

B＝{j,j+1} (14)

f_B＝A_V(B)＝(∑_i∈Bf_iω_i)/(∑_i∈Bω_i)………………………………(15)

ω_B＝ω_j+ω_j+1 (16)

并且令：

$\tilde{f}={(f_1,...,f_{j-1},f_B,f_{j+2},...,f_k)}^{\′}...\left(17\right)$

$\widetilde{\mathrm{\ω}}={({\mathrm{\ω}}_1,...,{\mathrm{\ω}}_{j-1},{\mathrm{\ω}}_B,{\mathrm{\ω}}_{j+2},...,{\mathrm{\ω}}_k)}^{\′}...\left(18\right)$

III.重复上述步骤II，直到把下标集k分解为l个小块B₁,B₂,…,B_l，并且满足条件：

A_V(B₁)＜A_V(B₂)＜…＜A_V(B_l)

则：

f_i ^*＝A_V(B_t),i∈B_t,i＝1,2,…,l………………………………(19)。