Bayesian Subset Selection for Weibull Populations

Abstract

Selecting the best population is an important problem in reliability analysis, particularly when lifetime data are right censored and the populations are characterized by threshold behavior. This paper presents a Bayesian subset selection procedure for multiple populations modeled by three-parameter Weibull distributions. The best population is defined as the population having the largest threshold parameter, which represents the minimum lifetime before failure can occur. Posterior probabilities of being the best are estimated based on Markov Chain Monte Carlo samples, and the selected subset is defined as the smallest collection of populations whose cumulative posterior probability of containing the best population reaches a prespecified level. Posterior inference is conducted using a random-walk Metropolis algorithm with a transformed parameterization that maintains the constraints of the three-parameter Weibull model. Convergence and sampling efficiency are evaluated through graphical and numerical diagnostics. The methodology is demonstrated using a semi-simulated fatigue-life study based on right-censored Alloy T7987 data. Model adequacy is assessed with posterior predictive methods, and posterior estimates are reported for the Weibull threshold, shape, and scale parameters for each population. In the Alloy-based study, the posterior probability of being best was distributed primarily among three populations, reflecting substantial uncertainty in the population ranking. At a target posterior probability level of 0.95, these three populations constituted the selected subset. The findings demonstrate that the proposed procedure incorporates posterior uncertainty and prevents the forced selection of a single population when the lifetime data do not offer clear separation. This approach offers a flexible Bayesian framework for selection among right-censored Weibull reliability populations.

Share and Cite:

Hodaj, J. and McDonald, G. (2026) Bayesian Subset Selection for Weibull Populations. Applied Mathematics, 17, 526-571. doi: 10.4236/am.2026.178030.

1. Introduction

Ranking and selecting among multiple populations is a key statistical problem with many uses in reliability, manufacturing, engineering design, and quality improvement. Often, the main goal is not just to estimate parameters, but to find the population that performs best based on a meaningful scientific standard. Traditional ranking-and-selection methods usually pick either a single population or a subset that is likely to include the best one. Early studies focused on populations with normal distributions and laid the groundwork for indifference-zone and subset-selection methods [1] [2]. Later research expanded these ideas to more types of distributions and applications, as discussed in detail by Gupta and Panchapakesan [3], Bechhofer et al. [4], and Gibbons et al. [5].

In reliability analysis, we usually assess how populations perform by looking at lifetime measures like average lifetime, certain percentiles, survival rates at specific times, or other parameters with clear physical meaning. Reliability data often include right-censored observations, since some units are still working at the end of the study or are removed before failing. Because of this, methods designed for complete data may not work as intended, and both parameter estimation and population comparisons need to account for censored data in the likelihood [6]-[8]. As a result, ranking and subset-selection methods for censored populations have become an important topic in reliability research [9].

The Weibull distribution is commonly used in lifetime analysis because it can model different patterns of failure rates. Its shape parameter lets it represent decreasing, constant, or increasing failure rates, while the scale parameter sets the typical lifetime for the population [10] [11]. The three-parameter Weibull model also includes a threshold, or location, parameter that marks the earliest possible failure time. In this setup, failures cannot happen before the threshold. This threshold is especially useful for reliability studies when populations differ in the earliest time failures can occur.

Adding a threshold parameter brings extra statistical challenges. The smallest observed lifetime limits the possible values for this parameter, and estimating it becomes harder when data are censored or when the likelihood behaves oddly near the boundary. Maximum likelihood estimation for Weibull models with both complete and censored data has been studied in depth [12], and several methods exist for estimating Weibull parameters and reliability percentiles [13] [14]. Still, there can be a lot of uncertainty in the threshold parameter, especially if populations have similar lifetime distributions or only a few failures are observed. So, simply picking the population with the highest threshold estimate might hide important uncertainty in the ranking.

Researchers have developed several frequentist methods for selecting or ranking Weibull populations. Kingston and Patel proposed ways to choose the best among several Weibull populations and to select restricted subsets [15] [16]. Hsu [17] and Gupta and Miescke [18] studied optimal and two-stage selection methods, while Tseng and Wu [19] looked at selecting more reliable Weibull populations when data are Type-II censored. Similar methods have been created for exponential populations, including ways to select good populations and compare threshold parameters [20]-[23]. Together, these studies highlight the need for formal selection rules that consider sampling variability, instead of just using estimated rankings.

Traditional subset-selection methods are usually set up to meet certain probability requirements through repeated sampling. For example, a method might be designed so that the chosen subset is likely to include the true best population with at least a specified probability across a range of parameter values. A Bayesian approach offers a different perspective. Once data are collected, uncertainty about the population parameters is described by their joint posterior distribution. The posterior probability that each population is best can be estimated directly, making it possible to measure the evidence for each ranking. Bayesian ranking and subset-selection methods have been used in areas like traffic safety and engineering design [24] [25], and have also been studied for censored Weibull populations and Weibull-component models [26] [27]. Unlike methods that rely only on point estimates, using posterior probabilities allows for keeping several possible candidates when the data do not clearly show which population is best.

This paper presents a Bayesian subset-selection method for several populations modeled with three-parameter Weibull distributions and right-censored data. Here, the best population is the one with the largest threshold parameter. For each population, we estimate the posterior probability of having the largest threshold using Markov Chain Monte Carlo samples. Populations are ranked by these probabilities, and the selected subset is the smallest group whose combined posterior probability reaches or exceeds a chosen target. This way, the method does not force the selection of a single population when there is still a lot of uncertainty.

We use a random-walk Metropolis algorithm for posterior inference, working in a transformed parameter space. This transformation keeps the Weibull shape and scale parameters positive and makes sure the threshold stays below the smallest observed lifetime. The Metropolis algorithm, along with Hastings’ generalization, provides a flexible way to sample from posterior distributions that do not have closed-form solutions [28]-[30]. We check convergence and sampling quality with both graphical and numerical diagnostics, and we assess model fit using posterior predictive methods [31]-[33].

We demonstrate the proposed method with a semi-simulated fatigue-life study using right-censored Alloy T7987 data. The study keeps the same sample size and main features as the original Alloy data, but creates several related reliability populations with different Weibull characteristics. This setup lets us test the subset-selection method in a setting that combines real engineering data with controlled differences between populations. We present the posterior estimates of the Weibull parameters, the posterior probabilities of being best, and the selected subset. We also check how sensitive the population ranking is to using different reliability criteria based on Weibull lifetime percentiles.

The rest of the paper is organized as follows. Section 2 covers the three-parameter Weibull model, the likelihood for right-censored data, and the transformed parameterization. Section 3 explains the prior and posterior distributions. Section 4 defines the best population, the posterior probabilities, and the Bayesian subset-selection rule. Section 5 describes the random-walk Metropolis method and how we check for convergence. Section 6 presents the Alloy T7987-based study and its results. Finally, we discuss the advantages, limitations, and possible extensions of our approach.

2. Three-Parameter Weibull Model

2.1. Three-Parameter Weibull Distribution

For population i , let the lifetime observations follow a three-parameter Weibull distribution with threshold parameter γ i , scale parameter η i , and shape parameter β i . The probability density function is

f i ( t| γ i , η i , β i )= β i η i ( t γ i η i ) β i 1 exp[ ( t γ i η i ) β i ],t> γ i , (2.1)

where η i >0 and β i >0 . The corresponding distribution function is

F i ( t| γ i , η i , β i )=1exp[ ( t γ i η i ) β i ],t> γ i , (2.2)

and the survival function is

S i ( t| γ i , η i , β i )=exp[ ( t γ i η i ) β i ],t> γ i . (2.3)

The threshold parameter γ i determines the lower endpoint of the lifetime distribution. Under the model, failures cannot occur before γ i . The scale parameter η i controls the scale of the lifetimes above the threshold, while the shape parameter β i determines the shape of the distribution and the behavior of its hazard function [10] [11].

The hazard function is

h i ( t| γ i , η i , β i )= f i ( t| γ i , η i , β i ) S i ( t| γ i , η i , β i ) = β i η i ( t γ i η i ) β i 1 ,t> γ i . (2.4)

Thus, the hazard rate is decreasing when β i <1 , constant when β i =1 , and increasing when β i >1 . This flexibility makes the Weibull distribution useful in fatigue-life and other reliability applications.

2.2. Likelihood for Right-Censored Data

The semi-simulated Alloy reliability study includes right-censored observations. Let t ij denote the observed time for unit j in population i , and let

d ij ={ 1, if the observation is an exact failure time, 0, if the observation is right-censored. (2.5)

The likelihood contribution for observation j in population i is

[ f i ( t ij | γ i , η i , β i ) ] d ij [ S i ( t ij | γ i , η i , β i ) ] 1 d ij . (2.6)

Therefore, the likelihood for population i is

L i ( γ i , η i , β i | t i , d i )= j=1 n [ f i ( t ij | γ i , η i , β i ) ] d ij [ S i ( t ij | γ i , η i , β i ) ] 1 d ij . (2.7)

The corresponding log-likelihood is

i ( γ i , η i , β i )= j=1 n d ij log f i ( t ij | γ i , η i , β i )+ j=1 n ( 1 d ij )log S i ( t ij | γ i , η i , β i ). (2.8)

Substituting the Weibull density and survival functions gives

i ( γ i , η i , β i )= j=1 n d ij [ log( β i ) β i log( η i )+( β i 1 )log( t ij γ i ) ] j=1 n ( t ij γ i η i ) β i . (2.9)

2.3. Transformation Used for Posterior Computation

The threshold parameter must satisfy

γ i < t i( 1 ) , (2.10)

where

t i( 1 ) = min 1jn t ij . (2.11)

Transformations of constrained parameters are commonly used to improve computational behavior and simplify posterior simulation [34] [35]. To enforce the parameter constraints during posterior simulation, the following transformation is used:

θ i1 =log( β i ), θ i2 =log( η i ), θ i3 =log( δ i ), (2.12)

where

δ i = t i( 1 ) γ i ε, (2.13)

and ε is a small positive constant. Equivalently,

γ i = t i( 1 ) exp( θ i3 )ε. (2.14)

Because exp( θ i3 )>0 , this parameterization guarantees that

γ i < t i( 1 ) .

For every proposed value of θ i3 . The logarithmic transformations similarly guarantee that β i >0 and η i >0 .

The prior distributions assigned to the transformed model parameters and the resulting posterior distribution are presented in the following section.

3. Bayesian Specification

3.1. Prior Distributions

The specification of priors for the three-parameter Weibull distribution requires consideration of the positivity constraints on the shape and scale parameters, as well as the threshold constraint on the location parameter. Independent prior distributions are assigned to the Weibull parameters for each population. Since β i >0 and η i >0 , gamma prior distributions are chosen for the shape and scale parameters. The threshold parameter is handled by introducing the positive distance parameter,

δ i = t i( 1 ) γ i ε,

where ε= 10 6 is a fixed numerical offset used to keep the threshold parameter away from the boundary t i( 1 ) . This offset helps prevent numerical instability when evaluating terms involving t ij γ i , is negligible relative to the scale of the observed fatigue lifetimes, and is not an estimated model parameter. The prior is assigned directly to the transformed parameter

θ i3 =log( δ i ),

which is defined on the real line. Prior specification for Weibull reliability models and the use of formal prior distributions in Bayesian analysis are discussed in [36]-[39].

In the semi-simulated Alloy fatigue study, the prior distributions are specified as follows:

β i ~Gamma( α β i =6, θ β i =0.25 ), (3.1)

η i ~Gamma( α η i =0.01, θ η i =10000 ), (3.2)

and

θ i3 =log( δ i )~N( μ θ i3 =0, σ θ i3 2 = 2 2 ). (3.3)

In Equations (3.1) and (3.2), the second parameter represents the scale parameter of the gamma distribution.

The prior for η i is specified on a scale that aligns with the observed fatigue lifetimes, while allowing for substantial uncertainty. The gamma prior for β i respects the positivity constraint on the Weibull shape parameter while remaining sufficiently flexible for the data to inform posterior inference. The prior assigned to θ i3 is particularly important, as the threshold parameter is related to the smallest observed time through

γ i = t i( 1 ) exp( θ i3 )ε.

Larger values of δ i correspond to smaller values of γ i , whereas smaller values of δ i position the threshold nearer to the smallest observed lifetime. Therefore, the prior on θ i3 determines the proximity of the threshold to the lower bound of the observed data. The Alloy study investigates the sensitivity of posterior conclusions to the prior standard deviation of θ i3 .

3.2. Posterior Distribution on the Transformed Scale

For population i , let

θ i = ( θ i1 , θ i2 , θ i3 ) T = ( log( β i ),log( η i ),log( δ i ) ) T . (3.4)

The original model parameters can be recovered from the transformed parameters using the following relationships:

β i =exp( θ i1 ), η i =exp( θ i2 ), δ i =exp( θ i3 ), (3.5)

and

γ i = t i( 1 ) exp( θ i3 )ε. (3.6)

Applying the right-censored likelihood from Section 2.2, the posterior density for population i on the transformed scale is proportional to the product of the likelihood, the prior densities, and the Jacobian of the transformation:

π( θ i | t i , d i ) L i ( γ i , η i , β i | t i , d i )π( β i )π( η i )π( θ i3 )| ( β i , η i ) ( θ i1 , θ i2 ) |. (3.7)

Because

β i =exp( θ i1 )and η i =exp( θ i2 ),

the Jacobian term is

| ( β i , η i ) ( θ i1 , θ i2 ) |=exp( θ i1 )exp( θ i2 )= β i η i . (3.8)

No additional Jacobian term is necessary for θ i3 because the normal prior in Equation (3.3) is assigned directly to θ i3 =log( δ i ) .

Therefore, the posterior density targeted in the analysis is

π( θ i | t i , d i ) L i ( γ i , η i , β i | t i , d i )π( β i )π( η i )π( θ i3 )exp( θ i1 )exp( θ i2 ). (3.9)

Equivalently, the log-posterior density is given by

i ( θ i )=log L i ( γ i , η i , β i | t i , d i )+logπ( β i )+logπ( η i ) +logπ( θ i3 )+ θ i1 + θ i2 +C, (3.10)

where C is a constant independent of the model parameters. The posterior distribution does not have a closed form and is evaluated using the random-walk Metropolis algorithm described in Section 5.

4. Bayesian Subset Selection

The Alloy fatigue study involves a comparison of five independent three-parameter Weibull populations. The subset-selection procedure is described below for a general case of k populations, with k=5 in this application. In contrast to the frequentist probability of correct selection, which is evaluated over repeated samples with fixed parameter values, the Bayesian approach utilizes the posterior distribution obtained after observing the right-censored lifetime data. This method quantifies the posterior evidence that each population possesses the largest threshold parameter and constructs a subset whose posterior probability of containing the best population meets a prespecified criterion.

4.1. Definition of the Best Population

For population i , let

ψ i =( γ i , η i , β i ) (4.1)

denote the vector of three-parameter Weibull parameters, where γ i is the threshold parameter, η i is the scale parameter, and β i is the shape parameter.

The best population is defined as the population with the largest threshold parameter. Thus, population i is best if

γ i = max 1jk γ j . (4.2)

This definition is motivated by reliability applications, where a larger threshold parameter indicates a greater lower endpoint of the lifetime distribution and, consequently, a later earliest possible failure time.

Accordingly, this criterion is intended for engineering decisions in which protection against early failure, or equivalently a later model-implied lower lifetime bound, is the primary objective. It is not intended as a universal ordering of Weibull populations. When the engineering objective concerns median lifetime, an upper lifetime percentile, expected lifetime, reliability at a specified time, or another performance measure, the corresponding reliability functional should determine the ranking criterion. The percentile comparisons reported in this paper are therefore supplementary analyses that illustrate how the population rankings may change under alternative reliability objectives; they are not part of the formal threshold-based subset-selection rule.

Let i * denote the index of the best population:

i * =arg max 1jk γ j . (4.3)

Because the posterior distributions considered here are continuous, exact ties among the threshold parameters occur with probability zero.

4.2. Posterior Probability of Being Best

After observing the right-censored lifetime data, uncertainty about the Weibull parameters is represented by their posterior distribution. The posterior probability that population i is best is defined as

p i =Pr( γ i = max 1jk γ j |t,d ), (4.4)

Here, t and d represent the observed times and censoring indicators for all populations. The values p 1 ,, p k summarize the posterior evidence supporting each population as the best. Because one population is identified as best in each posterior draw,

i=1 k p i =1. (4.5)

In practice, these probabilities are estimated from the retained posterior draws. After burn-in and thinning, suppose M draws of the threshold parameters are available:

{ γ 1 ( m ) ,, γ k ( m ) } m=1 M .

For draw m , the population having the largest sampled threshold parameter is recorded as the best population. The posterior probability that population i is best is estimated by

p ^ i = 1 M m=1 M I ( γ i ( m ) = max 1jk γ j ( m ) ). (4.6)

Here, I( ) denotes the indicator function. Therefore, p ^ i represents the proportion of retained posterior draws in which population i has the largest threshold parameter.

4.3. Posterior Probability of Correct Selection

Let S denote a selected subset of populations. Correct selection occurs when the selected subset contains the best population:

i * S. (4.7)

Under the Bayesian model, the posterior probability of correct selection is

Pr( CS|t,d )=Pr( i * S|t,d ). (4.8)

Since p i is the posterior probability that population i is best, this probability can be written as

Pr( CS|t,d )= iS p i . (4.9)

This posterior probability differs from the frequentist probability of correct selection. The frequentist probability is a pre-data measure evaluated over repeated samples with fixed parameter values. In contrast, Equation (4.9) is conditional on the observed Alloy data and the specified prior model.

4.4. Bayesian Selected Subset

Let P * denote the desired posterior probability level. To construct the selected subset, the estimated posterior probabilities of being best are ordered from largest to smallest:

p ^ ( 1 ) p ^ ( 2 ) p ^ ( k ) , (4.10)

where p ^ ( r ) denotes the rth largest estimated posterior probability and i ( r ) denotes the population index associated with that probability.

The Bayesian selected subset is

S B ={ i ( 1 ) ,, i ( m ) }, (4.11)

where m is the smallest integer satisfying

r=1 m p ^ ( r ) P * . (4.12)

Thus, S B is the smallest subset, determined by the ordering of the estimated posterior probabilities, whose estimated posterior probability of containing the best population meets or exceeds P * .

If a single population has a dominant posterior probability of being best, the selected subset may consist solely of that population. When posterior uncertainty is higher, additional populations are included until the desired posterior probability threshold is achieved. In the Alloy study, the primary analysis sets P * =0.95 .

4.5. Percentile-Based Reliability Comparisons

Although the formal definition of the best population and the Bayesian selected subset are based on the threshold parameter γ i , the same posterior comparison can be applied to Weibull lifetime percentiles. For a three-parameter Weibull distribution, the pth percentile for population i is

t p,i = γ i + η i [ log( 1p ) ] 1/ β i ,0<p<1. (4.13)

The 10th, 50th, and 90th percentiles are examined:

t 0.10,i , t 0.50,i , t 0.90,i .

The quantity t p,i indicates the time by which a proportion p of units from population i is expected to have failed according to the Weibull model.

For each retained posterior draw, the largest value of t p,i among the populations is identified. The posterior probability that population i has the largest pth percentile is estimated by

p ^ t p ,i = 1 M m=1 M I ( t p,i ( m ) = max 1jk t p,j ( m ) ). (4.14)

This calculation is performed separately for p=0.10,0.50 , and 0.90.

These percentile-based comparisons do not define the formal Bayesian selected subset, which remains based on the threshold parameter γ i . Rather, they serve as supplementary reliability summaries to evaluate whether populations with greater posterior support for the largest threshold parameter also tend to exhibit larger lifetime percentiles.

5. Posterior Computation

Posterior inference is conducted independently for each of the five Alloy populations using a random-walk Metropolis algorithm [28] [29]. Sampling is performed on the transformed parameter scale.

θ i = ( θ i1 , θ i2 , θ i3 ) T = ( log( β i ),log( η i ),log( δ i ) ) T ,

This transformation ensures that the positivity constraints on β i , η i , and δ i , as well as the threshold constraint γ i < t i( 1 ) , are automatically satisfied. After sampling, the retained draws are mapped back to the original Weibull parameter scale using

β i =exp( θ i1 ), η i =exp( θ i2 ), δ i =exp( θ i3 ),

and

γ i = t i( 1 ) exp( θ i3 )ε.

5.1. Random-Walk Metropolis Algorithm

Let

i ( θ i )=logπ( θ i | t i , d i )

denote the log-posterior density for population i on the transformed scale. At iteration m , a candidate value is generated from

θ i * ~ N 3 ( θ i ( m1 ) , c 2 Σ ^ pilot,i ), (5.1)

where c is a tuning constant and Σ ^ pilot,i is a covariance matrix estimated from a preliminary random-walk run.

Because the proposal distribution is symmetric, the proposal-density terms cancel in the Metropolis-Hastings ratio. The candidate is accepted with probability

α RW =min{ 1, π( θ i * | t i , d i ) π( θ i ( m1 ) | t i , d i ) }. (5.2)

Equivalently, computation is performed on the log scale using

log( α RW )=min{ 0, i ( θ i * ) i ( θ i ( m1 ) ) }. (5.3)

If the proposed value is not accepted, the current value is retained for the next iteration.

5.2. Starting Values and Proposal Tuning

For each population, starting values are obtained by fitting shifted Weibull models across a grid of possible values for δ i . These values are then used to initialize

( log( β i ),log( η i ),log( δ i ) ).

A preliminary random-walk run of 8000 iterations is conducted for each chain. The empirical covariance matrix from this pilot run is used as Σ ^ pilot,i , which accounts for the posterior scales and dependencies among the three transformed parameters. In the main random-walk proposal, the tuning constant is set to c=1.2 .

Four independent chains are run for each population. Each chain contains 200,000 iterations, of which the first 20,000 are discarded as burn-in. Every third post-burn-in draw is retained, yielding 60,000 retained draws per chain. The retained draws from the four chains are pooled for posterior estimation and subset selection.

Because the five populations are modeled independently and independent priors are assigned across populations, the joint posterior distribution factors into the product of the five population-specific posterior distributions. Joint posterior draws used in the subset-selection calculations are therefore formed by combining the corresponding retained draws from the separately fitted population models.

5.3. Convergence Assessment

Convergence and sampling performance for all five populations are assessed using traceplots, autocorrelation functions, acceptance rates, rank-normalized split and folded-split R ^ diagnostics, and bulk and tail effective sample sizes. Traceplots and autocorrelation functions are examined separately within each chain for the transformed parameters

log( β i ),log( η i ),log( δ i ).

The R ^ , bulk ESS, and tail ESS diagnostics are calculated jointly across the four chains while preserving chain identity. These diagnostics are evaluated for the Weibull parameters γ i , η i , β i , and δ i , and for the selected lifetime percentiles t 0.10,i , t 0.50,i , and t 0.90,i . After convergence has been assessed, the retained draws from the four chains are combined for posterior summaries and subset-selection probabilities. The resulting numerical summaries and representative graphical displays are presented alongside the Alloy study results in Section 6.

6. Alloy T7987 Fatigue-Life Study

The proposed Bayesian subset-selection procedure is illustrated using a semi-simulated five-population study based on the Alloy T7987 fatigue-life data. The original data were reported by Meeker, Escobar, and Pascual and are included in the reliability-data collection accompanying Statistical Methods for Reliability Data, second edition [7]. The dataset contains 72 specimens, of which 67 failed during testing and five were right-censored at 300 thousand cycles. Time is measured in thousands of cycles. The original Alloy data provide a realistic censored-data baseline, while the additional populations allow the subset-selection procedure to be examined when several reliability populations have similar lifetime distributions.

6.1. Original Alloy Data and Weibull Fit

Before constructing the five-population study, a three-parameter Weibull distribution was fit to the original Alloy data using the right-censored likelihood. The maximum likelihood estimates were

γ ^ =92.9928, η ^ =93.2641, β ^ =1.3202.

The estimated shape parameter is greater than one, corresponding to an increasing fitted hazard rate. This behavior is consistent with a fatigue process in which the instantaneous risk of failure increases as the number of cycles accumulates.

Figure 1 provides two complementary graphical assessments of the fitted model. The fitted Weibull survival curve follows the Kaplan-Meier estimate closely over most of the observed range and remains positive near 300 thousand cycles, as expected from the censored specimens. On the Weibull probability scale, the transformed Kaplan-Meier points are approximately linear, with most points lying close to the fitted line or within its approximate pointwise confidence band. These results support the three-parameter Weibull distribution as a reasonable working model for the Alloy fatigue-life data.

(a) Kaplan-Meier estimate and fitted Weibull survival curve

(b) Weibull probability plot with fitted line and approximate confidence band

Figure 1. Graphical assessment of the three-parameter Weibull fit to the original Alloy T7987 fatigue-life data.

6.2. Construction of the Five Populations

A Bayesian three-parameter Weibull model was first fitted to the original Alloy T7987 fatigue-life data. The posterior mean estimates of γ , η , and β from this fit were adopted as the baseline generating values for Population P1. The remainder of the semi-simulated study was constructed relative to this Bayesian baseline. Populations P2-P5 were generated with progressively larger threshold parameters, while their scale and shape parameters were kept close to the posterior estimates obtained from the original Alloy data. Thus, the five populations retain broadly similar fatigue-life behavior, while their intended ordering is determined primarily by the threshold parameter γ .

Table 1 gives the parameter values used in the construction and the resulting censoring summaries. Each population contains 72 observations, with censoring percentages ranging from 6.94% to 12.50%. Since P5 has the largest generating threshold parameter, γ 5 =100.993 , it is the true best population under the formal threshold-based definition.

Table 1. Construction and censoring summary for the semi-simulated Alloy fatigue-life populations.

Population

γ

η

β

n

Failures

Censored

Percent censored

P1

92.993

92.7443

1.2737

72

67

5

6.94

P2

94.993

90.8894

1.2482

72

63

9

12.50

P3

96.993

92.7443

1.2737

72

67

5

6.94

P4

98.993

94.5992

1.2992

72

65

7

9.72

P5

100.993

96.4541

1.3246

72

67

5

6.94

6.3. Posterior Estimates and Threshold-Based Selection

The censored three-parameter Weibull model was fit separately to each population using the random-walk Metropolis algorithm described in Section 5. Table 2 reports the posterior means and 95% credible intervals for the principal Weibull parameters. The posterior mean of γ i is largest for P5, followed by P4 and P3. However, the credible intervals for P3-P5 overlap substantially, reflecting the similarity of the populations and the uncertainty introduced by right-censoring. The posterior intervals for η i and β i also overlap, as expected from the construction of the study.

Table 2. Random-walk posterior means and 95% credible intervals for the Weibull parameters.

Population

γ : mean and 95% CrI

η : mean and 95% CrI

β : mean and 95% CrI

P1

92.6351 [87.4915, 93.9723]

94.4843 [77.1021, 114.5755]

1.3145 [1.0572, 1.6130]

P2

99.2982 [95.3508, 100.2460]

102.5622 [81.2505, 127.6475]

1.1512 [0.9149, 1.4246]

P3

100.3579 [94.9090, 101.7160]

98.7007 [80.0642, 120.2761]

1.2832 [1.0274, 1.5845]

P4

100.3876 [96.7317, 101.2585]

102.8857 [81.3138, 128.4027]

1.1231 [0.8970, 1.3835]

P5

101.1260 [95.7891, 102.5059]

97.1267 [78.9459, 117.9522]

1.3025 [1.0450, 1.6035]

Posterior means are shown in black; 95% credible intervals are shown in brown.

Figure 2 displays the random-walk posterior densities of the five threshold parameters on a common horizontal scale. The densities shift toward larger values from P1 to the upper-threshold populations, but substantial overlap remains among P3, P4, and P5. This overlap explains why the posterior evidence is distributed across several possible best populations rather than being concentrated entirely on P5.

Figure 2. Random-walk posterior densities of the threshold parameters for the five semi-simulated Alloy populations.

The ranked posterior probabilities of being best are shown in Table 3. P5 has the largest posterior probability, 0.6702, followed by P3 with probability 0.2178 and P4 with probability 0.1062. P2 has only a small probability of being best, and the estimated probability for P1 is zero to four decimal places.

Table 3. Random-walk posterior probabilities of being best, cumulative probabilities, and inclusion in the selected subset for P * =0.95 .

Rank

Population

p ^ i

Cumulative

Included

1

P5

0.6702

0.6702

Yes

2

P3

0.2178

0.8880

Yes

3

P4

0.1062

0.9942

Yes

4

P2

0.0058

1.0000

No

5

P1

0.0000

1.0000

No

The Monte Carlo uncertainty in the estimated posterior probabilities was assessed using a binomial approximation. Let X i denote the number of retained draws within a chain for which population i has the largest threshold parameter. Then

X i ~Binomial( n, p i ), p ^ i = X i n ,

where p i is the posterior probability that population i has the largest threshold. Therefore,

SD( p ^ i )= p i ( 1 p i ) n 1 2 n .

With n=60000 retained draws per chain, the upper bound is approximately 0.00204. Thus, under the binomial approximation, the Monte Carlo uncertainty in the estimated posterior probabilities of being best is small. Serial dependence in the MCMC draws is assessed separately through the chain-aware effective sample size diagnostics.

P5 alone does not attain the target posterior probability P * =0.95 . Adding P3 raises the cumulative probability to 0.8880, which remains below the target. After P4 is added, the cumulative probability becomes 0.9942. Therefore, the Bayesian selected subset is

S B ={ P5,P3,P4 }.

The selected subset contains the true best population while acknowledging that the observed censored data do not separate P5 sufficiently from P3 and P4 to justify selecting P5 alone.

6.4. Percentile-Based Reliability Assessments

The posterior Weibull percentiles provide a complementary comparison because they depend jointly on γ i , η i , and β i . Table 4 reports the posterior mean of each percentile and the posterior probability that the corresponding population has the largest percentile.

Table 4. Random-walk posterior means and probabilities of having the largest Weibull percentile.

Population

t 0.10

t 0.50

t 0.90

Mean

P( largest )

Mean

P( largest )

Mean

P( largest )

P1

109.7321

0.0100

163.9910

0.0268

272.2456

0.0210

P2

113.8767

0.0827

173.6953

0.2339

313.3131

0.3612

P3

117.4985

0.3431

174.3818

0.2428

291.0458

0.0926

P4

114.3337

0.0971

174.4447

0.2636

318.8830

0.4603

P5

118.4326

0.4670

174.2870

0.2329

286.8809

0.0650

For the early-life percentile t 0.10 , P5 has the largest posterior mean and the largest probability of being best, although that probability is only 0.4670. The posterior means of the median t 0.50 are nearly equal for P2-P5, and the probability of having the largest median is distributed across these populations. For t 0.90 , P4 and P2 receive the largest posterior probabilities. Thus, the upper-tail ranking differs from the threshold-based ranking because upper percentiles are particularly sensitive to uncertainty in the scale and shape parameters.

The posterior probabilities that the threshold-based and percentile-based criteria identify the same population are 0.3804, 0.2095, and 0.1281 for t 0.10 , t 0.50 , and t 0.90 , respectively. The decline in agreement illustrates that the formal threshold criterion and the broader percentile comparisons answer related but distinct reliability questions. The selected subset in this paper remains defined by the largest threshold parameter.

6.5. MCMC Diagnostics

The convergence and sampling diagnostics described in Section 5.3 were evaluated for all five populations. Figure 3 shows the random-walk traceplots for P5 as a representative case. The four chains remain stable and explore common posterior regions without visible long-term drift. The log( δ 5 ) chains exhibit greater variability than the log( β 5 ) and log( η 5 ) chains because the threshold-gap parameter is directly related to the lower-endpoint constraint.

Figure 3. Random Walk traceplots for Population P5. The columns correspond to log( β 5 ) , log( η 5 ) , and log( δ 5 ) , and the rows correspond to the four chains. The horizontal axis indexes the 60,000 retained post-burn-in draws from each chain after retaining every third draw.

Table 5 provides a compact summary of the random-walk diagnostics. The rank-normalized split/folded-split R ^ , bulk effective sample size, and tail effective sample size were calculated jointly across the four chains while preserving chain identity (Vehtari et al. [40]). Bulk ESS assesses sampling efficiency in the central region of the posterior distribution, whereas tail ESS assesses sampling efficiency in the posterior tails.

The mean random-walk acceptance rates range from 0.356 to 0.367. All maximum R ^ values are at most 1.0001. Across the model parameters γ , η , β , and δ , the minimum bulk and tail ESS values are 57515.1 and 74180.7, respectively. Across the Weibull percentiles t 0.10 , t 0.50 , and t 0.90 , the corresponding minimum values are 60267.6 and 79309.0. Autocorrelation functions were also examined separately within each chain for every population and decreased with increasing lag; they are not displayed to avoid repetitive graphical material.

Table 5. Summary of chain-aware random-walk MCMC diagnostics.

Population

Mean

acceptance

Maximum

R ^

Minimum

parameter

bulk ESS

Minimum

parameter

tail ESS

Minimum

percentile

bulk ESS

Minimum

percentile

tail ESS

P1

0.364

1.0001

57515.1

74901.4

61876.6

81398.3

P2

0.367

1.0001

58130.5

74180.7

62194.6

80512.8

P3

0.357

1.0000

58281.4

77910.1

61040.0

79655.2

P4

0.367

1.0001

60755.9

77662.6

62627.0

79309.0

P5

0.356

1.0000

58868.8

77531.5

60267.6

79861.1

Taken together, the trace plots, autocorrelation functions, acceptance rates, rank-normalized R ^ values, and bulk and tail effective sample sizes indicate satisfactory convergence and sampling efficiency for the retained random-walk posterior draws.

6.6. Prior Sensitivity and Posterior Predictive Checks

Sensitivity to the prior on the threshold-gap parameter was examined by varying the standard deviation of the Normal prior for θ i3 =log( δ i ) over 1, 2, and 3. Figure 4 shows that the posterior means of γ i remain relatively stable across the three specifications. The credible intervals widen, particularly through their lower limits, as the prior becomes more diffuse. P5 retains the largest posterior mean threshold under all three prior settings, indicating that the main threshold-based conclusion is not driven by the selected prior standard deviation.

A similar sensitivity analysis was done by varying the location parameter of the Normal prior distribution for θ i3 , log( δ i ) , over log( 0.5 ) , log( 1 ) , and log( 2 ) , while the prior standard deviation was held fixed at 2. As before the posterior means of γ i remained relatively stable. The corresponding 95% credible intervals widened as the prior median of δ i increased, primarily through decreases in their lower limits. The upper credible limits remained similar, and the intervals overlapped substantially across the prior specifications. Although the posterior probabilities of being best also varied numerically, the population ordering remained P5 > P3 > P4 > P2 > P1, and the selected subset at P * =0.95 remained {P5, P3, P4}. Thus, the main subset-selection conclusion was not sensitive to the prior location values considered.

Figure 4. Prior sensitivity of the posterior threshold estimates. Each panel shows the posterior mean and 95% credible interval for γ i as the prior standard deviation of log( δ i ) varies.

Posterior predictive checks were performed for all five fitted populations. Replicated lifetime datasets were generated from the posterior predictive distribution and subjected to the same right-censoring mechanism as the observed datasets. The checks included the minimum, mean, standard deviation, numbers of failures and censored observations, median, and lower- and upper-tail quantiles.

Figure 5 shows representative checks for the minimum, mean, and standard deviation of P5. In each panel, the vertical line marks the observed statistic. The observed values lie within the main body of the corresponding posterior predictive distributions. Across all five populations and all examined statistics, the posterior predictive p-values ranged from 0.267 to 0.797, with no values near zero or one. The fitted censored Weibull models therefore reproduce the principal features of the observed datasets, including their lower endpoints, centers, spreads, tails, and censoring patterns.

Overall, the semi-simulated Alloy study demonstrates the role of Bayesian subset selection in a censored reliability setting. P5 receives the largest posterior probability of being best, but the evidence is not sufficiently concentrated to select P5 alone at P * =0.95 . The resulting subset {P5, P3, P4} preserves the target posterior probability while representing the uncertainty among the leading populations. The calculations in this section were carried out using the R-code in Appendix A.

Figure 5. Representative posterior predictive checks for Population P5. The histograms show replicated values of the minimum, mean, and standard deviation; the vertical lines mark the corresponding observed values.

7. Discussion

The semi-simulated Alloy T7987 study illustrates how Bayesian subset selection can be used when several reliability populations have similar lifetime distributions and the available data are subject to right-censoring. The posterior probability that Population P5 has the largest threshold parameter is 0.6702, which is higher than the corresponding probability for any other population but is not sufficient to select P5 alone when the target posterior content is P * =0.95 . After P3 and P4 are added according to their posterior probabilities of being best, the cumulative probability reaches 0.9942. The resulting selected subset is therefore

S B ={ P5,P3,P4 }.

This result reflects the substantial overlap among the posterior distributions of the leading threshold parameters shown in Figure 2. The inclusion of multiple populations should not be interpreted as a failure to identify a useful result. Rather, it provides a direct representation of the uncertainty supported by the observed censored data and avoids overstating the evidence in favor of a single population.

The threshold-based and percentile-based comparisons also provide different perspectives on reliability performance. The formal selection rule defines the best population as the one with the largest threshold parameter, whereas a Weibull percentile depends jointly on γ i , η i , and β i . Consequently, a population with the largest lower lifetime bound need not have the largest median or upper-tail percentile. In the Alloy study, P5 has the largest posterior mean and posterior probability of being best for the early-life percentile t 0.10 , but the ranking becomes less concentrated for t 0.50 and differs more substantially for t 0.90 , as shown in Table 4. These differences emphasize that threshold-based and percentile-based comparisons answer related but distinct reliability questions. The threshold criterion, which defines “best” throughout this paper, is appropriate when the primary objective is to compare lower lifetime bounds, while percentile comparisons may be more informative when performance over the broader lifetime distribution is of interest.

A principal advantage of the Bayesian approach is that it produces directly interpretable posterior probabilities for each population being best. These probabilities can be used to construct a selected subset with a prespecified posterior content without forcing the analysis to identify a single winner when the evidence is uncertain. The same posterior simulation also provides parameter estimates, credible intervals, reliability percentiles, ranking probabilities, and uncertainty summaries within a unified framework. In addition, the right-censored Weibull likelihood permits both exact failure times and censored observations to contribute appropriately to the analysis.

The random-walk Metropolis implementation performed well in this study. Acceptance rates were similar across populations, the Gelman-Rubin diagnostics were essentially one, and the effective sample sizes were large, as summarized in Table 5. The traceplots showed stable exploration of the posterior distributions, and the posterior predictive checks indicated that the fitted models were able to reproduce the main characteristics of the observed datasets. The prior sensitivity analyses showed that the posterior threshold summaries changed somewhat across the examined prior standard deviations and prior locations for log( δ i ) , but the main population ranking and selected subset remained unchanged. Together, these results support the numerical reliability of the posterior summaries used for the subset-selection decision.

Several limitations should nevertheless be considered. First, only P1 is directly tied to the original Alloy T7987 data; Populations P2-P5 were generated relative to the Bayesian fit of that dataset. The application therefore provides a realistic semi-simulated illustration rather than a comparison of five independently observed physical populations. Second, the analysis assumes that each population follows a three-parameter Weibull distribution. Although the graphical assessment and posterior predictive checks indicate that the three-parameter Weibull model fits the present data reasonably well, the results could be affected if this distribution is not appropriate for data from other applications. Third, the population-specific models were fitted independently, so the procedure does not borrow information across populations even when their lifetime distributions are similar. Finally, the definition of the best population is tied to the largest threshold parameter and may not align with objectives based on a median, an upper percentile, expected lifetime, or another reliability criterion.

Future work could extend the procedure in several directions. Hierarchical formulations could permit partial pooling among related populations and may improve estimation when sample sizes are small. Alternative lifetime distributions, interval censoring, progressive censoring, and competing-risk settings could also be considered. Another extension would allow the decision criterion to depend jointly on the threshold and one or more reliability percentiles. A comprehensive repeated-sampling study could evaluate the procedure’s operating characteristics across different threshold separations, censoring levels, and sample sizes. Outcomes could include the probability of correct selection, the expected selected subset size, the probability of selecting the true best population alone, and population-specific inclusion probabilities. Computational work could investigate more efficient posterior samplers for applications involving a larger number of populations or more complex models. Most importantly, applying the procedure to studies containing multiple independently observed reliability populations would provide additional evidence about its practical performance.

8. Conclusions

This paper developed a Bayesian subset-selection procedure for comparing multiple populations that follow three-parameter Weibull distributions. The best population was defined as the one with the largest threshold parameter, and posterior simulation was used to estimate the probability that each population is best. Ordering these probabilities and accumulating them until a prescribed posterior content is reached yields a selected subset that explicitly accounts for uncertainty in the population ranking.

The method was illustrated using a semi-simulated censored study based on the Alloy T7987 fatigue-life data. A Bayesian three-parameter Weibull fit to the original data supplied the baseline parameter estimates used to construct the five-population study. Under the threshold-based criterion, P5 received the largest posterior probability of being best, but the evidence was not sufficiently concentrated to select P5 alone at P * =0.95 . The selected subset {P5, P3, P4} retained the target posterior probability while representing the uncertainty among the leading populations.

The study also showed that rankings based on the Weibull threshold need not coincide with rankings based on reliability percentiles, since the latter depend on the threshold, scale, and shape parameters jointly. Overall, the proposed procedure provides a probabilistically interpretable method for comparing multiple reliability populations in the presence of right-censored lifetime data. Further development for hierarchical models, alternative lifetime distributions, unequal sample sizes, and additional censoring structures would broaden its use in practical reliability studies.

Appendix A

Conflicts of Interest

The authors declare no conflicts of interest regarding the publication of this paper.

References

[1] Bechhofer, R.E. (1954) A Single-Sample Multiple Decision Procedure for Ranking Means of Normal Populations with Known Variances. The Annals of Mathematical Statistics, 25, 16-39.[CrossRef]
[2] Gupta, S.S. (1965) On Some Multiple Decision (Selection and Ranking) Rules. Technometrics, 7, 225-245.[CrossRef]
[3] Gupta, S.S. and Panchapakesan, S. (1979) Multiple Decision Procedures. John Wiley & Sons, Inc.
[4] Bechhofer, R.E., Santner, T.J. and Goldsman, D.M. (1995) Design and Analysis of Experiments for Statistical Selection, Screening, and Multiple Comparisons. John Wiley & Sons, Inc.
[5] Gibbons, J.D., Olkin, I. and Sobel, M. (1999) Selecting and Ordering Populations: A New Statistical Methodology. Society for Industrial and Applied Mathematics.[CrossRef]
[6] Lawless, J.F. (2002) Statistical Models and Methods for Lifetime Data. Wiley.[CrossRef]
[7] Meeker, W.Q., Escobar, L.A. and Pascual, F.G. (2022) Statistical Methods for Reliability Data. 2 Edition, John Wiley & Sons, Inc.
[8] Nelson, W. (2004) Applied Life Data Analysis. John Wiley & Sons, Inc.
[9] Singh, H. and Gill, A.N. (1995) Ranking and Subset Selection Procedures for Populations with Censored Data: A Review. In: Koul, H.L. and Deshpande, J.V., Eds., Institute of Mathematical Statistics Lecture Notes-Monograph Series, Institute of Mathematical Statistics, 219-230.[CrossRef]
[10] Rinne, H. (2009) The Weibull Distribution. Chapman & Hall/CRC.
[11] Jiang, R. and Murthy, D.N.P. (2011) A Study of Weibull Shape Parameter: Properties and Significance. Reliability Engineering & System Safety, 96, 1619-1626.[CrossRef]
[12] Cohen, A.C. (1965) Maximum Likelihood Estimation in the Weibull Distribution Based on Complete and on Censored Samples. Technometrics, 7, 579-588.[CrossRef]
[13] Barbero, E., Fernández-Sáez, J. and Navarro, C. (1999) On the Estimation of Percentiles of the Weibull Distribution. Journal of Materials Science Letters, 18, 1441-1443.[CrossRef]
[14] Evans, J., Kretschmann, D. and Green, D. (2019) Procedures for Estimation of Weibull Parameters.
[15] Kingston, J.V. and Patel, J.K. (1980) Selecting the Best One of Several Meifull Populations. Communications in Statistics-Theory and Methods, 9, 383-398.[CrossRef]
[16] Kingston, J.V. and Patel, J.K. (1980) A Restricted Subset Selection Procedure for Weibull Populations. Communications in Statistics-Theory and Methods, 9, 1371-1383.[CrossRef]
[17] An Hsu, T. (1982) On Some Optimal Selection Procedures for Weibull Populations. Communications in Statistics-Theory and Methods, 11, 2657-2668.[CrossRef]
[18] Gupta, S.S. and Miescke, K.J. (1986) Optimum Two-Stage Selection Procedures for Weibull Populations. Journal of Statistical Planning and Inference, 15, 147-156.[CrossRef]
[19] Sheng-Tsaing Tseng, and Wu, H.-J. (1990) Selecting, under Type-II Censoring, Weibull Populations That Are More Reliable. IEEE Transactions on Reliability, 39, 193-198.[CrossRef]
[20] Ng, C.K. (2013) Procedures for Selecting Good Exponential Populations. Communications in Statistics-Simulation and Computation, 42, 1681-1692.[CrossRef]
[21] McDonald, G.C. and Hodaj, J. (2025) Selection Rules for Exponential Population Threshold Parameters. Applied Mathematics, 16, 1-14.[CrossRef]
[22] McDonald, G.C. and Hodaj, J. (2025) Operating Characteristics of Subset Selection Rules for Exponential Population Threshold Parameters. Applied Mathematics, 16, 441-460.[CrossRef]
[23] McDonald, G.C. and Hodaj, J. (2025) A Note Comparing Two Subset Selection Procedures for the Threshold Parameters of Two Exponential Populations. Applied Mathematics, 16, 877-883.[CrossRef]
[24] McDonald, G.C. (2016) Applications of Subset Selection Procedures and Bayesian Ranking Methods in Analysis of Traffic Fatality Data. WIREs Computational Statistics, 8, 222-237.[CrossRef]
[25] Kondo, Y., Zidek, J.V., Taylor, C.G. and van Eeden, C. (2018) Bayesian Subset Selection Methods for Finding Engineering Design Values: An Application to Lumber Strength. Sankhya A, 80, 146-172.[CrossRef]
[26] Gupta, S.S. and He, S. (1997) Selecting the Best Weibull Population Based on Type-I Censored Data: A Bayesian Approach. Technical Report TR97-03C, Department of Statistics, Purdue University.
[27] Kondo, Y. and Zidek, J.V. (2013) Bayesian Nonparametric Subset Selection Procedures with Weibull Components. Technical Report Technical Report 273, Department of Statistics, The University of British Columbia.
[28] Metropolis, N., Rosenbluth, A.W., Rosenbluth, M.N., Teller, A.H. and Teller, E. (1953) Equation of State Calculations by Fast Computing Machines. The Journal of Chemical Physics, 21, 1087-1092.[CrossRef]
[29] Hastings, W.K. (1970) Monte Carlo Sampling Methods Using Markov Chains and Their Applications. Biometrika, 57, 97-109.[CrossRef]
[30] Chib, S. and Greenberg, E. (1995) Understanding the Metropolis-Hastings Algorithm. The American Statistician, 49, 327-335.[CrossRef]
[31] Gelman, A. and Rubin, D.B. (1992) Inference from Iterative Simulation Using Multiple Sequences. Statistical Science, 7, 457-472.[CrossRef]
[32] Brooks, S.P. and Gelman, A. (1998) General Methods for Monitoring Convergence of Iterative Simulations. Journal of Computational and Graphical Statistics, 7, 434-455.[CrossRef]
[33] Gelman, A., Meng, X.L. and Stern, H. (1996) Posterior Predictive Assessment of Model Fitness via Realized Discrepancies. Statistica Sinica, 6, 733-760.
[34] Christensen, R., Johnson, W., Branscum, A. and Hanson, T.E. (2010) Bayesian Ideas and Data Analysis. CRC Press.
[35] Gelman, A., Carlin, J.B., Stern, H.S., Dunson, D.B., Vehtari, A. and Rubin, D.B. (2013) Bayesian Data Analysis. 3rd Edition, Chapman and Hall/CRC.
[36] Martz, H.F. and Waller, R.A. (1982) Bayesian Reliability Analysis. John Wiley & Sons, Inc.
[37] Hamada, M.S., Wilson, A.G., Reese, S.C. and Martz, H.F. (2008) Bayesian Reliability. Springer.
[38] Kass, R.E. and Wasserman, L. (1996) The Selection of Prior Distributions by Formal Rules. Journal of the American Statistical Association, 91, 1343-1370.[CrossRef]
[39] Sun, D. (1997) A Note on Noninformative Priors for Weibull Distributions. Journal of Statistical Planning and Inference, 61, 319-338.[CrossRef]
[40] Vehtari, A., Gelman, A., Simpson, D., Carpenter, B. and Bürkner, P. (2021) Rank-normalization, Folding, and Localization: An Improved for Assessing Convergence of MCMC (with Discussion). Bayesian Analysis, 16, 667-718. [Google Scholar] [CrossRef]

Copyright © 2026 by authors and Scientific Research Publishing Inc.

Creative Commons License

This work and the related PDF file are licensed under a Creative Commons Attribution 4.0 International License.