To run a Weibull Analysis, start by loading WeibullR and
ReliaPlotR
The two-parameter Weibull distribution has CDF
\[F(t) = 1 - \exp\!\left[-\left(\frac{t}{\eta}\right)^{\!\beta}\right], \quad t > 0,\]
where \(\beta > 0\) is the shape (slope on probability paper) and \(\eta > 0\) is the scale (characteristic life, the time at which 63.2% of units have failed). A shape \(\beta < 1\) indicates infant-mortality failure (decreasing hazard rate); \(\beta = 1\) reduces to the exponential distribution (constant hazard); \(\beta > 1\) indicates wearout (increasing hazard). The three-parameter extension adds a location parameter \(\gamma\) (failure-free period): \(F(t) = 1 - \exp[-((t - \gamma)/\eta)^\beta]\).
For the lognormal distribution, \(\ln(T) \sim \mathcal{N}(\mu_{\log}, \sigma_{\log}^2)\); the probability paper transformation uses the standard normal quantile \(\Phi^{-1}(F)\).
Fitting methods. Rank regression (default in
WeibullR, method.fit = "rr-xony") linearizes the CDF on
probability paper and minimizes squared deviations. Maximum likelihood
estimation (method.fit = "mle") maximizes the
log-likelihood and is preferred for heavily censored data.
Goodness-of-fit is reported as R² for rank regression and log-likelihood
for MLE (Meeker and Escobar 1998).
Confidence intervals. The Fisher matrix method
(method.conf = "fm") approximates confidence bounds via the
asymptotic normal distribution of MLE estimates. Likelihood ratio bounds
(method.conf = "lrb") invert a chi-squared test and are
generally more accurate for small samples (Meeker
and Escobar 1998).
Next, create some failure data for 5 different machines that fail at time 30, 49, 82, 90, and 96 respectively.
Then use the WeibullR package to fit a Weibull
distribution to the data, and the plotly_wblr() function to
create a probability plot.
Let’s add right censored data to the previous example. In addition to the 5 machines that failed, add suspensions for 3 machines that did not fail (right censored) at times 100, 45, and 10 respectively.
To create an interval censored model, let’s use the inspection data from Silkworth, 2020.
inspection_data <- data.frame(
left = c(0, 6.12, 19.92, 29.64, 35.4, 39.72, 45.32, 52.32),
right = c(6.12, 19.92, 29.64, 35.4, 39.72, 45.32, 52.32, 63.48),
qty = c(5, 16, 12, 18, 18, 2, 6, 17)
)Then add suspension data for units surviving until the last inspection date.
Finally, add a fit and plot the results.
obj <- wblr(suspensions, interval = inspection_data)
obj <- wblr.fit(obj, method.fit = "mle")
obj <- wblr.conf(obj, method.conf = "fm", lty = 2)
suspensions <- as.vector(suspensions$time)
plotly_wblr(obj,
susp = suspensions, fitCol = "red", confCol = "red", intCol = "blue",
main = "Parts Cracking Inspection Interval Analysis",
ylab = "Cumulative % Cracked", xlab = "Inspection Time"
)To fit a 3P Weibull, let’s create some new failure data and plot the results.
To build a contour plot, let’s rerun the first example and use the
plotly_contour() function to create a plot.
Multiple wblr models can be overlaid on the same
probability plot by passing a list of objects. All objects must use the
same distribution family. Each model is rendered in a distinct color,
and clicking a legend entry toggles all traces for that model.
failures1 <- c(30, 49, 82, 90, 96)
failures2 <- c(20, 40, 60, 80, 100)
obj1 <- wblr.conf(wblr.fit(wblr(failures1), method.fit = "mle"), method.conf = "lrb")
obj2 <- wblr.conf(wblr.fit(wblr(failures2), method.fit = "mle"), method.conf = "lrb")
plotly_wblr(list(obj1, obj2), cols = c("steelblue", "tomato"))