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 XX is said to follow an MWD if its cumulative distribution function (CDF) and probability density function (PDF) are given by F(x)=1exp(axbeλx), F(x) = 1- \exp \big( -a x^b e^{\lambda x} \big), and f(x)=a(b+λx)xb1eλxexp(axbeλx), f(x) = a (b + \lambda x) x^{b-1} e^{\lambda x} \exp \big( -a x^b e^{\lambda x} \big),
where x>0x>0, a>0a>0 is the scale parameter, b0b \ge 0 is a shape parameter, and λ0\lambda \ge 0 is an acceleration or flexibility parameter controlling the rate of hazard growth.

Estimation methods

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

  • Maximum Likelihood Estimation (MLE)

The MLEs are obtained by maximizing the log-likelihood function: (a,b,λ;𝐱)=i=1nlog(a(b+λxi)xib1eλxi)ai=1nxibeλxi. \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}.

  • Least Squares Estimation (LSE)

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)<<x(n)x_{(1)}<x_{(2)}< \dots < x_{(n)} denote the ordered samples. The objective function is: Q(a,b,λ)=i=1n(F(x(i);a,b,λ)i0.3n+0.4)2. Q(a,b,\lambda) = \sum_{i=1}^n \bigg( F(x_{(i)};a,b,\lambda) - \frac{i-0.3}{n+0.4} \bigg)^2.

  • Weighted Least Squares Estimation (WLSE)

The WLSE extends LSE by introducing weights: QW(a,b,λ)=i=1nwi(F(x(i);a,b,λ)i0.3n+0.4)2, 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 wi=(n+1)2(n+2)i(ni+1),i=1,,n.w_i = \frac{(n+1)^2(n+2)}{i(n-i+1)}, \; i=1,\dots,n.

  • Maximum Product of Spacings (MPS)

The MPS estimates are obtained by maximizing: M(a,b,λ)=i=1n+11n+1log(F(x(i);a,b,λ)F(x(i1);a,b,λ)), 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,λ)=0F(x_{(0)};a,b,\lambda) = 0 and F(x(n+1);a,b,λ)=1F(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.