Fitting the Modified Weibull Distribution (MWD)

Introduction

This vignette presents the fitMWD() function, which is used to estimate the parameters of the Modified Weibull Distribution (MWD) using several estimation methods, including maximum likelihood (ML), least squares (LS), weighted least squares (WLS), and maximum product of spacings (MPS).

The MWD, introduced by Lai et al. (2003), has been widely used in reliability and survival analysis. A random variable \(X\) is said to follow an MWD if its cumulative distribution function (CDF) and probability density function (PDF) are given by \[ F(x) = 1- \exp \big( -a x^b e^{\lambda x} \big), \] and \[ f(x) = a (b + \lambda x) x^{b-1} e^{\lambda x} \exp \big( -a x^b e^{\lambda x} \big), \]
where \(x>0\), \(a>0\) is the scale parameter, \(b \ge 0\) is a shape parameter, and \(\lambda \ge 0\) is an acceleration or flexibility parameter controlling the rate of hazard growth.

Estimation methods

Let a random sample from \(X \sim MWD(a, b, \lambda)\) are observed as \(x_i, \; i=1,\dots, n\). For all estimation methods, optimization is performed using stats::optim.

The MLEs are obtained by maximizing the log-likelihood function: \[ \ell(a,b,\lambda; \mathbf{x}) = \sum_{i=1}^{n} \log \big(a (b + \lambda x_i) x_i^{b-1} e^{\lambda x_i} \big) - a \sum_{i=1}^{n} x_i^b e^{\lambda x_i}. \]

The LSE is obtained by minimizing the squared differences between the theoretical and empirical CDF values at the ordered sample points.

Let \(x_{(1)}<x_{(2)}< \dots < x_{(n)}\) denote the ordered samples. The objective function is: \[ Q(a,b,\lambda) = \sum_{i=1}^n \bigg( F(x_{(i)};a,b,\lambda) - \frac{i-0.3}{n+0.4} \bigg)^2. \]

The WLSE extends LSE by introducing weights: \[ Q^W(a,b,\lambda) = \sum_{i=1}^n w_i \bigg( F(x_{(i)};a,b,\lambda) - \frac{i-0.3}{n+0.4} \bigg)^2, \] where \(w_i = \frac{(n+1)^2(n+2)}{i(n-i+1)}, \; i=1,\dots,n.\)

The MPS estimates are obtained by maximizing: \[ M(a,b,\lambda) = \sum_{i=1}^{n+1} \frac{1}{n+1} \log \big( F(x_{(i)};a,b,\lambda) - F(x_{(i-1)};a,b,\lambda) \big), \] where \(F(x_{(0)};a,b,\lambda) = 0\) and \(F(x_{(n+1)};a,b,\lambda) = 1\).

Using fitMWD()

To illustrate the use of fitMWD(), we consider a simulated example.

# Load the package
library(SSReliabilityClaytonMWD)

# generate data from MWD(a, b, lambda)
n <- 100
a <- 0.75; b <- 1.25; lambda <- 0.60
# set seed
set.seed(123)
dat <- rMweibull(n, a, b, lambda)
# random initial points 
init <- runif(3)

MLE of the parameters:

# Fit MWD to `dat`
fit.mle <- fitMWD(data = dat, est.method = "mle", opt.method = "L-BFGS-B", starts = init,
                   lower = c(1e-05,1e-05,1e-05), upper = c(Inf,Inf,Inf), hessian = FALSE )
fit.mle$estimates
#>         a         b    lambda 
#> 0.7231634 1.2600843 0.6559157

LSE of the parameters:

fit.lse <- fitMWD(data = dat, est.method = "lse", opt.method = "L-BFGS-B", starts = init,
                   lower = c(1e-05,1e-05,1e-05), upper = c(Inf,Inf,Inf), hessian = FALSE )
fit.lse$estimates
#>         a         b    lambda 
#> 0.9299033 1.4069386 0.4020883

WLSE of the parameters:

fit.wlse <- fitMWD(data = dat, est.method = "wlse", opt.method = "L-BFGS-B", starts = init,
                    lower = c(1e-05,1e-05,1e-05), upper = c(Inf,Inf,Inf), hessian = FALSE )
fit.wlse$estimates
#>         a         b    lambda 
#> 0.9220048 1.4228131 0.4207337

MPS of the parameters:

fit.mps <- fitMWD(data = dat, est.method = "mps", opt.method = "L-BFGS-B", starts = init,
                   lower = c(1e-05,1e-05,1e-05), upper = c(Inf,Inf,Inf), hessian = FALSE )
fit.mps$estimates
#>        a        b   lambda 
#> 0.714056 1.189107 0.646072

References

Lai, C. D., Xie, M., and Murthy, D. N. P. (2003). A modified Weibull distribution.IEEE Transactions on Reliability, 52(1), 33–37.