Bayesian Estimation and Prediction for the Maxwell Failure Distribution Based on Type II Censored Data

We present Bayes estimators, highest posterior density (HPD) intervals, and maximum likelihood estimators (MLEs), for the Maxwell failure distribution based on Type II censored data, i.e. using the first r lifetimes from a group of n components under test. Reliability/Hazard function estimates, Bayes predictive distributions and highest posterior density prediction intervals for a future observation are also considered. Two data examples and a Monte Carlo simulation study are used to illustrate the results and to compare the performances of the different methods.

Tyagi and Bhattacharya (1989a, b) [8] [9] considered the Maxwell distribution as a lifetime model and obtained the Bayes and the minimum variance unbiased estimators of the Maxwell parameter and its reliability function.Howlader and Hossain (1998) [10] also considered the Maxwell distribution and provided Bayes estimates for the parameter θ, its associated reliability function ( ) ; R t θ , and highest posterior density intervals (HPD) for θ.
The purpose of this paper is to derive maximum likelihood estimators (or estimators that maximize the likelihood function), Bayes estimators in terms of the mean of the posterior distribution, and highest posterior density intervals for θ (which are intervals with a posterior probability of 1 α − and minimum length), the reliability function ( ) ; R t θ and the hazard function ( ) ; h t θ based on Type II censored data.In addition, we present Bayes predictive estimators based on predictive means and HPD prediction intervals for a future observation.When sampling is expensive and time consuming, the Type II censoring scheme can be used to save on the cost of the experiment and the data collection time.From a Bayesian perspective, estimation of lifetime models under censoring has been investigated by many authors (Balakrishnan (1990) [11], Bekker and Roux (2005) [12], Chaturvedi and Rani (1998) [13]).
The Maxwell probability density function (pdf) and cumulative distribution function (cdf), are respectively given by: ( ) 4 π e 0, 0, The pdf of Y is given by ( )

Estimation of θ, h(t;θ), and R(t;θ)
We assume a group of n components have lifetimes which follow a Maxwell distribution.The failure times are recorded as they occur until a fixed (known) number r of components have failed.As it is quite common in life testing situations, (e.g., destructive tests, high cost of testing an item, etc.) only the first r lifetimes in a sample of n units can be obtained.Let , , , , x , the maximum of the observed sample.This is what we refer to as Type II censoring.Instead of x we work with ( ) ( ) , , , , , , . Since the lifetimes in x follow a Maxwell distribution, Therefore, the likelihood function for θ in terms of y can be written as, where ( ) ( ) [ ] 1 and the log-likelihood is equal, except for a constant, to The maximum likelihood estimator (MLE) of θ, MLE θ , is found by numerical maximization of the loglikelihood.Relying on asymptotic normality properties of the MLE under Type II censoring as in Kong and Fei (1996) [15], an approximate ( ) confidence interval for θ can be obtained from, ( ) where ( ) By the invariance property of MLE's, we get that the MLE of the reliability function, ( ) Also, the MLE for the failure rate or hazard function, ( )

Bayes Estimation
Combining the likelihood function and the prior ( ) p θ , the posterior density of θ is given by the expression, A squared error loss is appropriate when decisions become gradually more damaging for larger errors.The Bayes estimator of θ under this loss function is the posterior expectation of θ , ( ) ( ) x y (11) In general, this expected value does not have a closed form solution.We rely on MCMC calculations to approximate ( ) For the case of the Maxwell distribution, transformed into a Gamma model, MCMC can be done with the Openbugs software as in Lunn et al. (2012) [16].
By graphing posterior realizations of ( ) , R t θ and ( ) , h t θ as a function of t, we an provide an uncertainty band for these functions in addition to the point estimates

Intervals for θ and R(t,θ)
A credible interval for θ can be obtained by taking our sample (1)  q α − respectively.( ) . This interval approximately satisfies the condition that ( ) Analogously, for the reliability function ( ) , R t θ , we transform our MCMC samples for θ and obtain a sample for ( ) based on these values provide a ( ) In a similar way, we can produce credible intervals for the hazard function ( ) It is not automatic to compute a highest posterior density (HPD) interval from a Monte Carlo sample specially if the posterior is far from symmetric or multimodal.However, we use the method in Liu, Gelman and Zheng (2013) [17], to obtain HPD intervals from a MCMC sample based on a numerical approximation to HPD intervals integrated in the SPIn R-package.We use this method to compute HPD intervals for parameters or predictions for this paper.

Bayes Predictive Estimator and Intervals
Let z be a future observation which has already survived ( ) r x .Given the data x, the conditional joint p.d.f of z and θ is (16) and ( ) , f z θ x takes into account censoring.To produce a Monte Carlo sample of the predictive distribution ( ) , , , M z z z  , of values that correspond to the predictive distribution.The Bayes estimator of z under squared error loss function can be approximated as .
A ( ) for z can be found by making * , , , M y y y  . The HPD interval can be computed with the SPIn method.In our numerical examples, we found that the SPIn method could provide some different results depending on how we thin the MCMC sample and select the M value.At this point, it is worth mentioning that Bayesian estimation and prediction of the two-parameter Gamma distribution has been considered in Pradhan and Kundu (2011) [18].However this other paper does not explicitly address estimation and prediction under type II censoring as we propose here.

Example 1: National Radio Astronomy Observatory Data Set
The following data represent noise levels in cryogenic microwave receivers and were obtained from Darrell Hicks at the National Radio Astronomy Observatory, Socorro, NM.We arranged the observations in ascending order and dropped the last 11 data points to induce censoring, which leads to a situation with n = 86 and r = 75.After transforming the data with x  1.After a burn-in of 1000 MCMC iterations, the posterior mean estimate of θ based on 5000 additional samples results in Bayes ˆ20.34 θ = .MCMC convergence is reached after only a few iterations and was checked with time series plots of the iterations (trace plots) and autocorrelation plots.Based on our sample and as described in Section 3, a 95% credible interval for θ is (16.98,24.36).For comparison and relying on MLE asymptotic theory, we also computed an approximate 95% confidence interval for θ which resulted in the interval (16.51, 23.58).For this calculation, the expected Fisher's information was approximated with the observed information obtained from the Hessian calculation (matrix of second derivatives) resulting from the Nelder-Mead, Quasi-Newton method.
In  ; h t θ (bottom).From this figure we observe that when t approaches the value of 10, the reliability function values are close to zero.For the hazard function, the failure rates are monotonically increasing as a function of t.Based on the sampling scheme described in Section 4, Figure 2 (bottom) also shows samples of the posterior predictive for a future observation z taking into account censoring and its approximate density function.The predictive mean for a future observation is 10.40 and a 95% predictive interval for a future observation is (9.63, 13.10).These values remain essentially the same using either Jeffrey's prior or Hartigan's prior on θ.The results obtained with the SPIn R-package to compute HPD intervals from our MCMC samples are as follows.The HPD interval for θ is (16.78278,24.48906) which is not very different compared to an interval based on the 0.25 and 0.975 sample quantiles.The HPD for the predictive distribution of a future z observation is (9.36, 12.82) which also remains similar to the one reported with 0.25, 0.975 sample quantiles.In Table 1, we report 95% HPD intervals for future observations for different levels of censoring.It should be noted that the width of the interval increases as the level of censoring increases.On the other hand, the predictive interval without censoring captures all the variability that is present in the full data set.
In Table 2, we report HPD intervals for θ corresponding to this data set under 3 different priors: Hartigan's, Jeffreys' and a Gamma(1,1) prior which is an Exponential distribution of decay parameter equal to one.The lengths of the HPD intervals for θ under Hartigan prior are slightly shorter compared to the intervals obtained with Jeffreys' and Gamma priors for both censoring times considered.
In addition, we also computed the approximate 95% confidence interval for θ with the MLE asymptotic approximation to the Normal distribution as described in Section 2. For a censoring time of 64(0.64), the interval is (0.1706, 0.2790) and for a censoring time of 68 (0.68), the interval is (0.1681, 0.2697).Again, the expected Fisher's information was approximated with the observed Fisher's information at the MLE.In contrast to the NRAO data set example, the HPD intervals provide different results compared to a MLE approximation that relies in asymptotic normality.

Monte Carlo Simulation Studies
In order to assess the performance of the estimation and prediction approaches proposed in this paper, we perform a Monte Carlo simulation study based on 5000 simulated data samples of sizes n = 15, 25, 35, 50 and 100 with 10% and 30% censoring respectively.For each of our simulated data samples, we computed the Bayes estimates using a squared error loss function after a burn-in of 500 MCMC iterations.For different values of t and using the above mentioned sample sizes, we computed the average Bayes estimates, ( ) ( ) * h t , as in Section 2. The main results are presented in Table 3 and Table 4.The true value of θ for the simulation was selected equal to 20.5.The "true" reliabilities for this θ value for t = 0.5, 2.5, 5.0 and 7.5, are respectively: 0.9990, 0.8942, 0.4864 and 0.1394 and the 'true' hazard function values are respectively: 0.0060, 0.1253, 0.3691 and 0.6312.We have done some other simulations using different values of θ and the pattern of the numerical results remain similar to those reported here.
Furthermore, for a second simulation based on 5000 samples and for various values of θ, in Figure 4, we show Mean Square Error (MSE) curves for a posterior θ computed under Hartigan, Jeffreys' prior and the MLE method respectively.For each data simulation and parameter, we produce an estimate through the 3 approaches and compute a square error difference between the true parameter value and the estimated value.An average was obtained across simulation and graphed as a function of θ.The results show that as θ gets larger, the MSE values are increasing.We also note that except for the first case, where no censoring was considered, the MSE curves under Hartigan's prior are lower than those obtained with Jeffreys or MLE.Various values of n and r were considered in this study and defined the four panels shows in Figure 4.The top left panel refers to a case where there is no censoring N = 86, r = 86, the top right panel is for N = 86, r = 81.For the bottom panel, we have N = 86, r = 70 (right panel) and N = 86, r = 65 for the left panel.We notice that the differences between the curves corresponding to each prior are more marked when the censoring level increases.

Conclusions
The main contribution of this paper is to obtain the parameter estimates for the Maxwell failure distribution under Type II censoring via the Gamma distribution, using Bayesian estimation under different priors and compared it with Maximum Likelihood Estimation (MLE).MLE can be thought as the maximum of the posterior distribution under an unrestricted uniform prior.We observe from Table 1 (NRAO data example) that the length of the 95% HPD prediction intervals for future observations increases as the percentage of censoring increases.The width of the HPD intervals under censoring is shorter as compared to non censored data.It can be noted from Table 2 (Burning velocity data example) that at censoring t = 64 and censoring t = 68, the HPD intervals for θ under Hartigan's prior are slightly shorter than for both Gamma and Jeffreys priors.From   4 (simulation) when censoring is set equal to 30%.The values of ( ) * h t increase when t increases for all the sample sizes studied for both 10% and 30% censoring.For the estimation of the parameter θ, with Type II censoring, the MSE calculations of Figure 4 show that Hartigan's prior has a better performance than Jefferys' prior and MLE.This is true for all of the values of the parameter θ that are studied and all values of n and r considered but it is more noticeable for cases where there is more censoring.
For estimating the parameter θ of the Maxwell distribution under Type II censoring, it appears to be clear from all our numerical results that Bayes estimation is appropriate or in some cases, superior than MLE, but with the MLE method as a good competitor.In the second data example for which the sample size n = 55, an approximate 95% confidence interval for θ provides a result that is different to HPD intervals.A limitation of our approach is the use of a square data transformation for the Gamma distribution which can pose challenges under the Type II censoring for large numerical data values.The priors considered for comparisons are obtained under the non-censored case, therefore, a similar study can be attempted computing Jeffreys' and Hartigan's priors numerically and where the Type II censoring is incorporated into the likelihood function.An Openbugs model that implements the approach described in this manuscript with the Gamma and truncated Gamma distributions, is available by request from the second author.
time of the ith component to fail.Since the remaining ( ) n r − components have not yet failed and thus, they have lifetimes greater than [ ] r

Figure 2 (
top), we show a histogram of the posterior samples of θ along with the density representation of its posterior distribution.

Figure 3
presents the posterior mean (in solid black line) and 200 posterior samples for the reliability function ( ) ; R t θ (top) and the hazard function ( )

Figure 2 .
Figure 2. NRAO data set: Histogram of posterior samples and posterior density of θ (top).Predictive samples and distribution for a future observation (bottom).

Figure 4 .
Figure 4. Simulation study: MSE comparisons for Hartigan, Jeffreys and MLE with various values of N and r and based on 5000 data simulations.by New Mexico Tech, Socorro, New Mexico, USA and research funds from the College of Arts and Sciences, The University of New Mexico,USA.G. Huerta performed all the computations and simulations applied in this paper.He also proposed the model from a Gamma distribution perspective.The Openbugs model for this paper is available from G. Huerta.
In Section 4, we describe how Bayes predictive estimators and highest posterior density prediction intervals for a future observation are produced.Finally, Sections 5 and 6 provide results for two data examples and a simulation study.The conclusions are presented in Section 7.
; h t θ .In Section 3, we discuss how HPD intervals for θ and ( ) , R t θ are obtained.
The MLE of θ obtained numerically with a Nelder-Mead, Quasi-Newton algorithm (in 100-th units) is MLE θ = 20.17.A graph of the MLEs, ( ) MLE ; R t θ and ( ) MLE ; h t θ is shown in Figure

Table 1 .
HPD-intervals for a future observation.

Table 3
(simulation study), it can be observed that the values of the Bayes estimates *R t are larger for all samples

Table 3 .
Estimates of * R t and ( ) * h t for 10% censoring.

Table 4 .
Estimates of * R t and ( ) * h t for 30% censoring.* R t decreases for all sample sizes.The same pattern is observed for Table