基于记录值数据的一种浴盆型寿命分布置信集估计
Confidence Set Estimation for a Bathtub-Shaped Lifetime Distribution Based on Recorded Data
摘要: 针对一种具浴盆型失效率曲线的寿命分布,本文在记录值数据场合研究了模型参数的置信集估计问题。首先,基于记录值样本构造枢轴量,进而分别建立模型参数的精确置信区间和置信域。其次,为补充起见,进一步利用大样本原理构造了模型参数的渐近置信区间估计。最后,通过数值模拟和算例分析比较了不同置信集结果的优良性。
Abstract: This paper investigates the confidence set estimation problem of model parameters for a lifetime distribution with a bathtub-shaped failure rate curve in the case of recorded data. Firstly, based on the recorded value samples, pivotal quantities are constructed, and then exact confidence intervals and confidence regions for the model parameters are established separately. Secondly, for the sake of supplementation, the asymptotic confidence interval estimation of the model parameters was further constructed using the large sample principle. Finally, the superiority of different confidence set results was compared through numerical simulation and case analysis.
文章引用:杨艳焱. 基于记录值数据的一种浴盆型寿命分布置信集估计[J]. 统计学与应用, 2024, 13(5): 1972-1981. https://doi.org/10.12677/sa.2024.135192

1. 引言

在社会领域和工业应用中,受限于检测环境和操作机制的多样性,真实数据往往呈现出各种特征。其中,记录值数据是一种特殊的数据表达形式,其定义最初是由Chandler [1]提出,可视为一种由观测值及其出现顺序确定样本大小的次序统计量。在许多实际应用中,记录值及其相关统计量发挥着重要作用,涉及气象、水文、体育和寿命测试等领域。例如,当足够的垂直力作用于木材梁时,它会断裂,该过程中木材的裂纹数据呈现出记录值的特征;又如,在体育赛事期间,各项运动的最好成绩不断被刷新,也以出记录值数据方式呈现。简而言之,记录值数据表示在数据观测和收集过程中,仅获得并记录比以前所有值更大(或更小)的一类数据形式。假设 { X n :n=1,2, } 是一组独立同分布的随机变量,令 Y n =max{ X 1 , X 2 ,, X n } ,如果对任意的 j>i T j > T i ,则称 Y j 是一个上记录值,否则,可获得相应的下记录值。记录值数据在实际应用领域具有广泛应用,很多学者都针对此类数据进行了深入讨论。例如,Hassan等[2]、Azhad等[3]基于上记录值研究了产品的强度应力可靠性推断问题;龙沁怡和徐丽平[4]、程丹等[5]讨论了不同分布模型场合记录值数据的统计推断方法。更详细讨论可参考Nagaraja [6]的工作。

在生产实践中,产品在实际使用时不仅存在早期失效(高故障率),也可能存在长期稳定运行期(低故障率),以及后期的再次增加故障率。针对实际产品的这一故障特征,实践中经常使用具有浴盆型失效率曲线的分布模型进行建模研究,具有这一特征的分布模型在机械、机电、电子、能源等领域具有广泛应用,深入理解和分析具有浴盆型曲线的寿命分布,可为实际数据分析提供有力支持。例如,Nassar和Dobbah [7]、Wang和Wu [8]、Chen和Gui等[9]在不同数据场合,讨论了多种具有浴盆型失效率分布的统计推断问题。文献[10]-[12]探讨了具有浴盆型失效率分布模型的贝叶斯估计方法。

假设随机变量X服从威布尔半逻辑分布(Weibull-Half-Logistic Distribution, WHL),其分布函数(Cumulative distribution function, CDF)和密度函数(probability density function, PDF)可表示为

F( x )=1 e α ( e x 1 2 ) β ,x,α,β>0 f( x )= αβ 2 ( e x 1 2 ) β1 e xα ( e x 1 2 ) β ,x>0 (1)

其中, α 是尺度参数, β 是形状参数。为简单起见,记为 WHL( α,β ) 。相应的,WHL分布的失效率函数(hazard rate function, HRF)可表示为

h( x )= αβ 2 e x ( e x 1 2 ) β1 (2)

为说明期间,图1给出WHL分布在不同参数下的密度函数和失效率函数曲线,可以发现WHL分布具有不同的浴盆型失效率曲线。

Figure 1. The probability density (left) and failure rate function (right) of WHL distribution

1. WHL分布的概率密度(左)和失效率函数(右)的图像

鉴于具有浴盆型失效率曲线的分布模型在实际研究中的广泛应用,以及记录值数据在工程实践中的广泛性,本文在记录值数据下讨论WHL模型参数的置信集估计问题。相比传统点估计结果,本文所涉及置信区间和置信域研究具有以下优势。首先,针对记录值数据样本量比较稀疏的特征,传统点估计结果在小样本场合下效果较差,置信集估计具有更好的稳健性:其次,置信集估计可提供更多参数信息,不仅反映了对参数的估计值,还提供了估计的不确定性度量;最后,置信集估计为决策者提供了一个关于参数可能值的信任范围,有助于辅助实践中作出更加周到的决策。综上考虑,本文讨论在记录值数据下,WHL寿命模型参数的置信集估计问题,利用枢轴量法和极大似然估计方法分别建立了相应参数的精确和渐近置信区间\置信域估计,并比较了两者的优良性。

2. 置信集估计

2.1. 精确置信集估计

为了建立模型参数的精确置信集估计,利用下面给出的定理构造相应的枢轴量。

定理1 假设样本 x=( x 1 , x 2 ,, x m ) 是来自 WHL( α,β ) 的记录值数据,定义枢轴量。

γ( β )= κ/2 ε/ 2( m1 ) = 1m 1 ( e x m 1 e x 1 1 ) β λ( α,β )=2α ( e x m 1 2 ) β (3)

从而 γ( β ) 服从自由度为 ( 2,2( m1 ) ) F分布, λ( α,β ) 服从自由度为 2m 的卡方分布,并且 γ( β ) λ( α,β ) 相互独立。

证明 Y i =α ( e x i 1 2 ) β ,i=1,2,,m 是来自标准指数分布的相应的记录值样本。根据Nagaraja [6]提出的理论,做如下变换

{ Z 1 = Y 1 , Z 2 = Y 2 Y 1 , Z m = Y m Y m1 (4)

则变换后 Z 1 , Z 2 ,, Z m 是来自标准指数分布的独立同分布的随机变量。根据指数分布的可加性可得 k=2 Z 1 =2α ( e x 1 1 2 ) β

ε=2 i=2 m Z i =2α[ ( e x m 1 2 ) β ( e x 1 1 2 ) β ] (5)

k ε 相互独立,并且分别服从自由度为2和 2( m1 ) 的卡方分布。进一步,可根据抽样分布理论构造枢轴量

γ( β )= k/2 ε/ 2( m1 ) = 1m 1 ( e x m 1 e x 1 1 ) β (6)

λ( α,β )=k+ε=2α ( e x m 1 2 ) β (7)

其中 γ( β ) 服从自由度为 ( 2,2( m1 ) ) F分布, λ( α,β ) 服从自由度为 2m 的卡方分布,并且二者相互独立。

为确保可在任意水平下构造置信区间和置信域,给出如下引理。

引理1 对任意的常数 a b ,满足 0<a<b ,记 k( t )= ( e b 1 e a 1 ) t ,t>0 ,则 k( t ) 满足:

(1) k( t ) 是关于 t 递增的函数;

(2) lim t0 k( t )=1 lim t k( t )=

证明 k( t ) 关于 t 求导数可得

dk( t ) dt =t ( e b 1 e a 1 ) t1 ,t>0 (8)

由此可见, k( t ) 是关于 t 单调递增的函数,并且 lim t0 k( t )=1 lim t k( t )=

由引理1可以看出 γ( β ) 是关于 β 单调递减的。

定理2 设 T={ T 1 , T 2 ,, T m } 是来自 WHL( α,β ) 的上记录值,对于任意的 0<p<1 ,参数 β 100( 1p )% 精确置信区间(Exact confidence interval, ECI)为

( ψ( F p/2 2,2(m1) ),ψ( F 1p/2 2,2(m1) ) ) (9)

其中, F p u 1 , u 2 是自由度为 u 1 , u 2 F分布的 100p% 上侧分位数, ψ( θ ) 是等式 k( t )=θ 关于 β 的解。

证明 基于定理1可知,对于任意的 0<p<1 ,有

P{ F 1p/2 2,2(m1) <γ( β )< F p/2 2,2(m1) }=1p (10)

由引理1可知, γ( β ) 是关于 β 单调递减的。因此,上述表达式等价于

P{ ψ( F p/2 2,2(m1) )<β<ψ( F 1p/2 2,2(m1) ) }=1p (11)

所以,参数 β 100( 1p )% 精确置信区间(ECI)可以表示为

( ψ( F p/2 2,2(m1) ),ψ( F 1p/2 2,2(m1) ) ) (12)

结论成立。 □

定理3 T={ T 1 , T 2 ,, T m } 是来自 WHL( α,β ) 的上记录值样本,对于任意的 0<p<1 和给定参数 β ,参数 α 100( 1p )% ECI为

( χ 1p/2 2m B( β ) , χ p/2 2m B( β ) ) ,这里 B( β )=2 ( e x m 1 2 ) β (13)

其中, χ p u 是自由度为 u 的卡方分布的 100p% 上侧分位数。

证明 对于给定的参数 β ,基于定理1中给出的枢轴量 λ( α,β ) 的性质,利用定理2中的方法可以类似证明。此处省略了证明过程。 □

进一步,建立了参数 ( α,β ) 的精确置信域(Exact confidence region, ECR)估计。

定理4 T={ T 1 , T 2 ,, T m } 是来自 WHL( α,β ) 的上记录值样本,对于任意的 0<p<1 ,参数 ( α,β ) 100( 1p )% ECR为

{ ( α,β )|ψ( F 1 1p 2 2,2(m1) )<β<ψ( F 1+ 1p 2 2,2(m1) ), χ 1+ 1p 2 2m B( β ) <α< χ 1 1p 2 2m B( β ) } (14)

证明 由于 γ( β ) λ( α,β ) 是相互独立的,使用与定理1类似的方法,对于任意的 0<p<1 ,有

P{ F 1+ 1p 2 2,2(m1) <γ( β )< F 1 1p 2 2,2(m1) , χ 1+ 1p 2 2m <λ( α,β )< χ 1 1p 2 2m } =P{ ψ( F 1 1p 2 2,2(m1) )<β<ψ( F 1+ 1p 2 2,2(m1) ) }×P{ ( χ 1+ 1p 2 2m / B( β ) )<α<( χ 1 1p 2 2m / B( β ) ) } = 1p × 1p =1p (15)

由此便可得到参数 ( α,β ) 100( 1p )% 精确置信域估计。

2.2. 近似置信集估计

为比较前面所构造精确置信集结果的优良性,这里基于极大似然估计的渐近理论进一步给出模型参数的渐近置信区间(approximate confidence intervals, ACIs)。

对于失效时刻 T={ T 1 , T 2 ,, T m } ,设其为来自 WHL( α,β ) 的上记录值, x 1 , x 2 ,, x m 是对应的观测值。则对应的似然函数表示为

L( α,β;T )=f( x m ) i=1 m1 f( x i ) 1F( x i ) = αβ 2 ( e x m 1 2 ) β1 e x m α ( e x m 1 2 ) β i=1 m1 αβ 2 ( e x i 1 2 ) β1 e x i α ( e x i 1 2 ) β e α ( e x i 1 2 ) β = e α ( e x m 1 2 ) β i=1 m αβ 2 ( e x i 1 2 ) β1 e x i (16)

由(16)式得到相应的对数似然函数为

( α,β;T )=mlogα+mlogβ+β i=1 m log ( e x i 1 2 )α ( e x m 1 2 ) β (17)

利用(17)式分别对参数 α β 求导数,得到 β 的极大似然估计值 β ^

β ^ = [ log( e x m 1 2 ) 1 m i=1 m log ( e x i 1 2 ) ] 1 (18)

α 的极大似然估计值 α ^

α ^ =m ( e x m 1 2 ) β (19)

参数向量 ( α,β ) 的Feisher信息矩阵的表达式为

I( α,β )= [ 2 ( α,β ) αα 2 ( α,β ) αβ 2 ( α,β ) αβ 2 ( α,β ) ββ ] α= α ^ ,β= β ^ (20)

这里

2 ( α,β ) αα = m α 2 2 ( α,β ) αβ = ( e x m 1 2 ) β log( e x m 1 2 ) 2 ( α,β ) ββ = m β 2 α ( e x m 1 2 ) β [ log( e x m 1 2 ) ] 2 (21)

由渐近分布理论可知, ( ( α ^ , β ^ )( α,β ) ) T N( 0, I 1 ( α ^ , β ^ ) ) ,其中 ( α,β ) 的方差–协方差矩阵为

I 1 ( α ^ , β ^ )=[ Var( α ^ ) Cov( α ^ , β ^ ) Cov( α ^ , β ^ ) Var( β ^ ) ] (22)

因此,对任意的 0<p<1 ,参数 α,β 100( 1p )% 近似置信区间构造如下

( α ^ + v 1p/2 Var( α ^ ) , α ^ + v p/2 Var( α ^ ) ) ( β ^ + v 1p/2 Var( β ^ ) , β ^ + v p/2 Var( β ^ ) ) (23)

其中, v p 是标准正态分布的上侧 100p% 分位数。

3. 数值模拟

为说明不同估计结果的优良能,本节进行了大量的数值模拟比较。在给定显著性水平 p=0.05 、模型参数真值设定和样本量n场合,分别给出参数 α,β 不同置信区间长度(Length),以及置信域面积(Area)作为比较准则。通过进行10,000次数值模拟,相应数值结果分别由表1表2给出,其中,产生记录值样本的流程由算法1建立。

Algorithm 1. Process of generating recorded value samples

算法1. 产生记录值样本流程

步骤1 从均匀分布 U( 0,1 ) 中产生独立分布的样本 Z 1 , Z 2 ,, Z n

步骤2 通过变换 y i =ln( 1 Z i ),i=1,,n 产生来组标准指数分布的样本。

步骤3 由 w i = y 1 + y 2 ++ y i ,i=1,,n 产生来自标准指数分布的记录值样本。

步骤4 令 u i =1 e w i ,i=1,,n ,则 u i 是来自均匀分布的记录值样本。

步骤5 对于任意的CDF F( x ) ,令 r i = F 1 [ 1( 1 u i ) ],i=1,,n ,则 r i 是来自 F( x ) 的记录值样本。

Table 1. The superiority of confidence set estimation under the true value of parameters α = 1.8, β = 2.5, p = 0.05

1. 参数真值为α = 1.8,β = 2.5,p = 0.05场合置信集估计优良性

n = 3

n = 4

n = 5

n = 7

α

ACI

Low

−0.6552

−0.0439

0.2495

0.5837

Up

5.0430

4.5506

4.3139

4.0352

Length

5.6982

4.5946

4.0644

3.4515

ECI

Low

0.4524

0.6140

0.7409

0.9285

Up

5.2833

4.9389

4.6737

4.3086

Length

4.8309

4.3250

3.9328

3.3801

β

ACI

Low

−0.2964

0.1044

0.3479

0.6549

Up

4.5624

4.2124

4.0403

3.8177

Length

4.8587

4.1080

3.6924

3.1628

ECI

Low

0.1514

0.2775

0.3841

0.5479

Up

3.8493

3.8339

3.8384

3.8514

Length

3.6979

3.5564

3.4543

3.3035

Table 2. The superiority of confidence set estimation under the true value of parameters α = 2, β = 3, p = 0.05

2. 参数真值为α = 2,β = 3,p = 0.05场合置信集估计优良性

n = 3

n = 4

n = 5

n = 7

α

ACI

Low

−0.6293

−0.0969

0.2820

0.6928

Up

5.4696

5.1849

4.9975

4.8194

Length

6.0989

5.2818

4.7155

4.1266

ECI

Low

0.4991

0.6931

0.8571

1.1081

Up

5.8283

5.5759

5.4070

5.1418

Length

5.3292

4.8828

4.5499

4.0337

续表

β

ACI

Low

−0.3048

0.0837

0.3602

0.7005

Up

4.7850

4.4846

4.2861

4.0768

Length

5.0898

4.4009

3.9259

3.3763

ECI

Low

0.1557

0.2835

0.3880

0.5539

Up

3.9577

3.9176

3.8777

3.8936

Length

3.8020

3.6340

3.4897

3.3397

(α, β)

ECR

Area

26.9851

23.6258

21.4812

18.7598

表1表2中相应的数值结果中可以得到以下结论:

① 随着样本量n的增大,参数 α,β 的ACI和ECI的区间长度逐渐减小,置信域面积也有逐渐变小的趋势,反映出无论是经典似然方法还是枢轴量方法,所建立置信集估计精度均随着样本量的增加而提升。

② 在样本量n固定时,基于枢轴量得到的ECI区间长度短于基于似然方法得到的ACI区间长度,反映了利用枢轴量所得到的置信集估计具有更好的优良性。

以上数值模拟结果表明,本文所提出基于枢轴量的精确置信集估计方法相较传统基于似然的渐近置信集估计效果更好。

4. 实例分析

本节基于一组真实数据说明所提出估计方法的有效性和实用性。本节考虑洛杉矶市民中心记录的年降雨量(英寸)真实数据集,它可以在洛杉矶年鉴的网站上找到。这里使用了25年间年降雨量数据,由于这组数据数值较大,且数据之间差距较大,为了方便计算和估计,将原始数据除以8,变换后降雨量数据由表3所示。

Table 3. Annual rainfall (inches)

3. 年降雨量(英寸)

K-S距离

p

1.6025

2.2325

0.9575

0.3100

1.0100

0.9188

1.4988

2.6250

0.9200

0.20021

0.2353

1.0138

3.0438

1.5550

1.5500

3.8763

1.1363

1.4463

2.2425

0.5525

2.0525

1.1563

4.7450

1.6488

0.4013

1.6913

1.1350

在给出结果之前,首先检验WHL分布是否可用来拟合这一数据集。基于表3中的完全数据,得到拟合优度检验中的K-S距离和相应的p值,结果仍在表3中呈现。通过检验结果可知,K-S距离为0.20021,p值为0.2353,远远大于0.05,故WHL分布可以作为一个合适的模型来拟合这组真实数据。

根据表3中的完全数据,得到对应的年降雨量记录值样本如下所示:

1.6025,2.2325,2.6250,3.0438,3.8763,4.7450

相应地,表4在95%水平下利用上述记录值样本下给出参数的各种置信区间和置信域估计结果。由表4结果可以看出,参数α的ECI拥有比ACI更短的区间长度;相同地,参数β的ECI也拥有比ACI更短的区间长度,这一结果与数值模拟中的结论相一致。

Table 4. Confidence set estimation of WHL model parameters under rainfall record data

4. 降雨记录值数据下WHL模型参数的置信集估计

α

ACI

Low

0.1263

Up

1.1393

Length

1.0130

ECI

Low

0.2322

Up

1.2306

Length

0.9984

β

ACI

Low

0.1231

Up

0.9896

Length

0.8665

ECI

Low

0.1937

Up

1.5745

Length

1.3808

(α, β)

ECR

Area

1.4004

5. 结论

本文研究了一种具有浴盆型失效率寿命分布的模型参数的置信集估计问题。在记录值样本场合,通过构造枢轴量,建立了模型参数的精确置信区间和置信域估计。此外,还利用了极大似然方法构造了渐近置信集估计。数值模拟和实例分析表明,本文所提出基于枢轴量的精确置信集估计相较传统渐近置信集估计具有更好的优良性,表明了本文提出的方法是可行的。

参考文献

[1] Chandler, K.N. (1952) The Distribution and Frequency of Record Values. Journal of the Royal Statistical Society Series B: Statistical Methodology, 14, 220-228.
https://doi.org/10.1111/j.2517-6161.1952.tb00115.x
[2] Hassan, A.S., Nagy, H.F., Muhammed, H.Z. and Saad, M.S. (2020) Estimation of Multicomponent Stress-Strength Reliability Following Weibull Distribution Based on Upper Record Values. Journal of Taibah University for Science, 14, 244-253.
https://doi.org/10.1080/16583655.2020.1721751
[3] Azhad, Q.J., Arshad, M. and Khandelwal, N. (2021) Statistical Inference of Reliability in Multicomponent Stress Strength Model for Pareto Distribution Based on Upper Record Values. International Journal of Modelling and Simulation, 42, 319-334.
https://doi.org/10.1080/02286203.2021.1891496
[4] 龙沁怡, 徐丽平. 基于上记录值Lomax分布的统计推断[J]. 江西师范大学学报(自然科学版), 2023, 47(2): 216-220.
[5] 程丹, 席吉富, 罗子怡, 等. 基于下记录值的逆指数分布模型的Bayes估计[J]. 黑龙江科学, 2023, 14(12): 10-13.
[6] Nagaraja, H.N. (1988) Record Values and Related Statistics—A Review. Communications in Statistics-Theory and Methods, 17, 2223-2238.
https://doi.org/10.1080/03610928808829743
[7] Nassar, M. and Dobbah, S.A. (2020) Analysis of Reliability Characteristics of Bathtub-Shaped Distribution under Adaptive Type-I Progressive Hybrid Censoring. IEEE Access, 8, 181796-181806.
https://doi.org/10.1109/access.2020.3029023
[8] Wang, L., Wu, K., Tripathi, Y.M. and Lodhi, C. (2020) Reliability Analysis of Multicomponent Stress-Strength Reliability from a Bathtub-Shaped Distribution. Journal of Applied Statistics, 49, 122-142.
https://doi.org/10.1080/02664763.2020.1803808
[9] Chen, S. and Gui, W. (2020) Statistical Analysis of a Lifetime Distribution with a Bathtub-Shaped Failure Rate Function under Adaptive Progressive Type-II Censoring. Mathematics, 8, Article 670.
https://doi.org/10.3390/math8050670
[10] Basikhasteh, M., Lak, F. and Afshari, M. (2020) Bayesian Estimation of Stress-Strength Reliability for Two-Parameter Bathtub-Shaped Lifetime Distribution Based on Maximum Ranked Set Sampling with Unequal Samples. Journal of Statistical Computation and Simulation, 90, 2975-2990.
https://doi.org/10.1080/00949655.2020.1793155
[11] Algarni, A. and Almarashi, A.M. (2020) E‐Bayesian Estimation of Reliability Characteristics of Two‐Parameter Bathtub‐shaped Lifetime Distribution with Application. Quality and Reliability Engineering International, 37, 1635-1649.
https://doi.org/10.1002/qre.2817
[12] 徐晓岭, 邢兆飞, 王蓉华, 等. 一种具有“浴盆”失效率的寿命分布的性质与统计分析[J]. 数理统计与管理, 2017, 36(6): 1001-1015.