Executive Summary
Weibull life verdict across 299 units.
The short answer
The fitted failure rate falls with age (shape parameter 0.833, 95% CI 0.695–0.999), a pattern called infant mortality in which the weakest patients fail early and survivors become more reliable. Half of all patients are expected to have failed by 317 followup days; 10% by 33 days.
The detail
Shape (beta) is 0.833 with 95% CI 0.695 to 0.999, entirely below 1, so the failure rate falls as patients age post-discharge. Characteristic life (eta) is 491.74 followup days—the point at which 63.2% are expected to have failed. B10 life (10% failure) is 33.03 days (95% CI 22.924 to 47.591); B50 life (50% failure) is 316.749 days (95% CI 240.645 to 416.921). Mean life (MTBF) is 541.808 days (95% CI 359.84 to 815.796), but 66.2% of patients are expected to have failed by then. The data covers 299 patients with 96 observed deaths and 203 censored (still living). Censoring rate is 67.9%. The fit extends to 712 followup days—2.50 times the observed maximum of 285—which is a model projection, not an observation.
What this can't tell you
The confidence intervals are wide, especially for mean life (range 359.84 to 815.796 days), reflecting the precision limits of 96 observed failures. The visible curvature in the probability plot suggests the data may contain more than one failure mode; a single Weibull fit averages them. Stratifying by clinical subgroup would clarify whether the infant-mortality pattern holds uniformly or differs by risk.
Analysis Overview
Weibull life fit on 299 units with 96 observed failures.
The short answer
Heart failure patients show a decreasing failure risk over followup time—a pattern called infant mortality, where the weakest patients fail early and survivors become progressively more reliable. The analysis tracks 299 patients over up to 285 followup days, with 96 observed deaths and 203 still living when observation stopped.
The detail
This fit uses a two-parameter Weibull model estimated by maximum likelihood on 299 units with 96 observed failures and 67.9% censored. The shape parameter beta is 0.833 (95% CI 0.695 to 0.999), entirely below 1, which defines the infant-mortality pattern: the failure rate falls from 0.00424 per followup days at day 2 to 0.00186 at day 284. The censored units—203 patients still living—are not discarded; each contributes the full time it was observed, which prevents downward bias in life estimates. The fit reaches a probability-plot R-squared of 0.9094, meaning the data is loosely straight in Weibull coordinates with visible curvature; the largest gap to the nonparametric Kaplan-Meier curve is 3.81 percentage points.
What this can't tell you
The Weibull model assumes a single homogeneous failure process. The visible curvature in the probability plot—an R-squared of 0.9094 rather than near 1—suggests the data may contain more than one failure mode, and a single Weibull fit averages them rather than resolving them separately. Stratifying by a clinical or demographic variable would clarify whether distinct subpopulations are present. The projection past 285 followup days is model assumption, not observation; no patient was followed that far, so the curve's behavior beyond that point cannot be confirmed by the data.
Data Quality
Row exclusions, how the failure flag was read, and the censoring rate.
The short answer
All 299 rows were usable; no records were dropped for missing or invalid data. The failure flag "1" was read as death, and the data follows patients to 285 followup days, with the last observed death at 241 days. Censoring rate is 67.9%—203 patients were still alive when observation ended.
The detail
299 rows loaded; 299 used. Zero rows removed. The died column was dichotomized so that "1" = failure (death) and all other values = censored (still living). Censoring rate: 67.9%, meaning 203 of 299 units contributed right-censored survival times. Each censored patient's record contributes the fact that they survived at least their recorded followup days, which is the key to preventing bias in life estimates. The data reaches 285 followup days maximum; the last observed failure (death) occurred at 241 followup days. This span—from discharge to 285 days—is the window the data can directly speak to.
What this can't tell you
The preprocessing report does not indicate whether any patients were lost to followup or whether the 203 censored patients represent a mix of still-living and administratively censored records. If loss to followup is present and unrelated to survival risk, the analysis remains valid; if loss is informative (e.g., sicker patients transferred elsewhere), it would bias the estimates. Clarifying the censoring mechanism—administrative end of study versus ongoing followup—would strengthen confidence in the survival estimates.
Weibull Parameters & Failure Mode
Shape and scale with 95% confidence intervals, plus fit quality.
| Parameter | Estimate | CI Low | CI High | Interpretation |
|---|---|---|---|---|
| Shape (beta) | 0.8333 | 0.6954 | 0.9986 | Below 1 means the failure rate falls with age, above 1 means it rises, 1 means it is flat. Here: infant mortality (failure rate falling with age). 96 observed failures carry the shape estimate; the interval below is the honest range. |
| Scale (eta) — characteristic life | 491.7 | 356.7 | 677.9 | The followup days by which 63.2% of units have failed — the natural scale of this life distribution, not a safe operating limit. |
| Observed failures carrying the fit | 96 | — | — | Censored units (203, 67.9% of the file) contribute the time they survived, but only failures pin down the shape. |
| Probability-plot R-squared | 0.9094 | — | — | How straight the data is in Weibull coordinates: loosely straight, with visible curvature. This scores the fit, and a good fit is still a model choice, not a fact. |
| Largest gap to the observed curve (percentage points) | 3.81 | — | — | The largest distance between the fitted reliability curve and the assumption-free Kaplan-Meier curve computed from the same data. |
The short answer
The shape parameter of 0.833 is the entire story: it sits below 1 (95% CI 0.695 to 0.999), meaning failure risk drops as patients age after discharge. The 96 observed deaths carry this estimate; the fit stays within 3.81 percentage points of the non-parametric Kaplan-Meier curve.
The detail
Shape (beta) = 0.8333 (95% CI 0.6954 to 0.9986); scale (eta) = 491.736 days (95% CI 356.674 to 677.941). The p-value testing shape against 1 is 0.048. The 96 observed failures pin down the shape; 203 censored units (67.9% of the file) contribute survival time but not failure information. Probability-plot R-squared is 0.9094, showing loosely straight alignment with visible curvature. The largest gap to the Kaplan-Meier curve is 3.81 percentage points—the fitted model stays close to the observed data but imposes a Weibull structure that the curvature suggests may oversimplify.
What this can't tell you
The visible curvature in the probability plot indicates the Weibull may be averaging more than one failure mode. A single straight line in Weibull coordinates would be needed to confidently rule out mixture of failure mechanisms. Consider whether clinical subgroups (e.g., by ejection fraction or comorbidity) follow separate failure patterns.
Reliability Curve R(t)
Share of units still running against followup days, fitted and projected, with the observed curve for comparison.
The short answer
The fitted reliability curve shows 53.0% of patients still alive at 285 followup days (the data's maximum), declining to 25.6% at 712 days (a model projection 2.50 times beyond the data). The smooth curve tracks the stepped Kaplan-Meier curve closely within the observed range, indicating the Weibull assumption is doing no harm there.
The detail
The fitted Weibull curve gives reliability (share still running) at any followup time: 100% at day 0, declining to 53.0% at 285 followup days. The Kaplan-Meier estimate—computed from the same data without assuming a distribution—tracks the fitted curve closely within the observed range (0 to 285 days), with the largest gap of 3.81 percentage points. Everything drawn past 285 followup days is labeled a separate series for a reason: it is the model extended to 712 days (2.50 times the longest observation) and assumes the same failure process continues unchanged. Nothing in this file can confirm that assumption.
What this can't tell you
The projection to 712 days is purely model-based; no patient was followed that long. If the underlying failure mechanism changes after 285 days—for example, if late complications emerge or if the cohort composition shifts due to selective censoring—the projected curve will be wrong. The close agreement between Weibull and Kaplan-Meier within the observed range does not guarantee the fit will extrapolate accurately. Longer followup data would test the projection directly.
Life Metrics
B10, B50, characteristic life, MTBF, and reliability at two named times.
| Metric | Estimate | CI Low | CI High | Interpretation |
|---|---|---|---|---|
| B10 life (10% of units failed) | 33.03 | 22.92 | 47.59 | The design life engineers usually quote: 10% of units are expected to have failed by 33.0 followup days. |
| B50 life (half of units failed) | 316.7 | 240.6 | 416.9 | Half of all units are expected to have failed by 317 followup days. |
| Characteristic life (63.2% failed) | 491.7 | 356.7 | 677.9 | The Weibull scale parameter, always the 63.2% point whatever the shape. |
| Mean life (MTBF) | 541.8 | 359.8 | 815.8 | The average life implied by the fit (eta times the gamma function of 1 + 1/beta). It is the mean, not a survival guarantee: 66.2% of units are expected to have failed by then. |
| Reliability at the end of the data (285) | 53.01 | — | — | The share still running at the last followup days actually covered by this data — the furthest point the file can speak to. |
| Reliability projected to 712 | 25.61 | — | — | The fitted survival share at 712 followup days. Beyond 285 followup days there is no data — this number is the fitted model extended forward, not an observation. |
The short answer
B10 life (10% failure) is 33.0 days; B50 life (50% failure) is 317 days; mean life is 542 days. At 285 followup days (the end of observed data), 53.0% of patients are still fitted to be alive. All life estimates carry 95% confidence intervals reflecting the precision of 96 observed deaths.
The detail
B10 life is 33.03 followup days (95% CI 22.924 to 47.591)—the design-life number: 10% of patients are expected to have failed by this point. B50 life is 316.749 days (95% CI 240.645 to 416.921)—half are expected to have failed by this point. Characteristic life (the Weibull scale) is 491.736 days (95% CI 356.674 to 677.941), the point at which 63.2% are expected to have failed. Mean life (MTBF) is 541.808 days (95% CI 359.84 to 815.796), but this is the average, not a safe operating point: 66.2% of patients are expected to have failed by then. Reliability at 285 followup days (the end of observed data) is 53.01%—the furthest point the file can speak to. Reliability at 712 followup days (projected) is 25.61%—beyond 285 days there is no data, so this is the fitted model extended forward, not an observation.
What this can't tell you
The confidence intervals are wide, especially for mean life (359.84 to 815.796 days), reflecting the precision limit of 96 observed failures in a cohort of 299. The B10 interval (22.9 to 47.6 days) spans a range of nearly 2× at the lower end of the life distribution, which means individual-patient predictions at short followup times are imprecise. Consider whether the sample size or followup duration is sufficient for the clinical decision being made.
Hazard Curve h(t)
The instantaneous failure rate against followup days.
The short answer
The instantaneous failure rate falls steadily from 0.00424 per followup days at day 2 to 0.00186 at day 284. This declining hazard is the infant-mortality pattern: weaker patients fail early, and survivors face lower risk as they age.
The detail
The hazard h(t) = (beta/eta) × (t/eta)^(beta−1). With beta = 0.833 and eta = 491.74, the hazard is 0.0042 per followup days at t=2 and 0.0019 at t=284, declining throughout the observed range. A falling hazard is operationally distinct from a rising one: a falling hazard makes burn-in or screening pay (remove weak units early), while a rising hazard makes age-based scheduled replacement pay. Here the reading is infant mortality: the failure rate falls with age because the weak units fail early and the survivors are progressively more reliable. The curve past 285 followup days is drawn as a separate series because no unit was observed that long; it is a model projection.
What this can't tell you
The hazard curve is derived from the fitted Weibull model; its shape depends entirely on the assumption that beta = 0.833 throughout the followup. If the true failure mechanism changes after 285 days—for example, if late-stage complications emerge with a different failure rate—the projected hazard will be wrong. The visible curvature in the probability plot suggests the data may contain more than one failure mode, which would mean the fitted hazard averages across them rather than resolving each. Longer followup or stratification by clinical subgroup would test whether the declining hazard holds uniformly.
Weibull Probability Plot
Failure data in Weibull coordinates, with the fitted line.
The short answer
The 69 observed failures plot loosely on a straight line in Weibull coordinates (R-squared 0.909), with visible curvature that suggests the data may contain more than one failure mode being averaged together.
The detail
In Weibull coordinates (ln(time) versus ln(−ln(reliability))), a single failure mechanism produces a straight line with slope beta = 0.833. The Kaplan-Meier plotting positions account for the 67.9% censoring rate, making them valid for this high-censoring regime. The 69 observed failures score R-squared = 0.9094 against the fitted line. The curvature visible in the plot—a gentle but systematic departure from straightness—indicates the Weibull is a reasonable but not perfect fit. Systematic curvature or a break into two segments would signal two distinct failure mechanisms; here the curvature is moderate, not a clear break.
What this can't tell you
The curvature alone does not prove mixture of failure modes, only that a single Weibull may be oversimplifying. A formal test for mixture (e.g., likelihood-ratio test against a two-component Weibull) would be needed to confirm. Clinical stratification by risk group or failure mechanism would clarify whether the observed curvature reflects true patient heterogeneity or random scatter.
Weibull Failure Analysis — How Components Fail Over Time
Fits a two-parameter Weibull life distribution to time-to-failure data with right-censoring (units still running when observation stopped), and reports the two numbers reliability engineering is built on: the shape parameter beta — whose value separates infant mortality from random failure from wear-out — and the characteristic life eta.
Why This Method?
Reliability data is almost never complete: at the end of a test or a warranty window some units have failed and the rest are still running. Dropping the survivors biases every life estimate downward; treating them as failures biases it too. Maximum-likelihood Weibull fitting with right-censoring uses each survivor for exactly as long as it was observed, which is what makes B10 life, MTBF and the reliability curve honest numbers rather than optimistic ones.
What This Analysis Covers
- Shape (beta) and scale (eta) with 95% confidence intervals
- The failure-mode reading of beta, hedged by its confidence interval
- B10 and B50 life, characteristic life, and mean life (MTBF)
- The reliability curve R(t) with a clearly-labelled forward projection
- The hazard curve h(t) — whether the failure rate rises or falls
- A Weibull probability plot with the fitted line and a fit-quality score
Standard Library
Platform standard-library module (LAT-1441): runs on ANY dataset via the semantic mapping {time, event}. All narrative is derived from the user's own column names and computed values.
suppressPackageStartupMessages(library(DT))
suppressPackageStartupMessages(library(htmlwidgets))
suppressPackageStartupMessages(library(arrow))
suppressPackageStartupMessages(library(knitr))
suppressPackageStartupMessages(library(rmarkdown))
suppressPackageStartupMessages(library(dplyr))
suppressPackageStartupMessages(library(tidyr))
suppressPackageStartupMessages(library(ggplot2))
suppressPackageStartupMessages(library(stringr))
suppressPackageStartupMessages(library(lubridate))
suppressPackageStartupMessages(library(broom))
suppressPackageStartupMessages(library(Matrix))
suppressPackageStartupMessages(library(cluster))
suppressPackageStartupMessages(library(data.table))
suppressPackageStartupMessages(library(survival))Core Analysis Pipeline
Step 1: Discover mapped columns (humanized names for ALL prose)
initial_rows <- nrow(df)
for (k in c("time", "event")) {
if (!k %in% names(df)) {
stop(sprintf("column_mapping must map the '%s' column.", k))
}
}
h_time <- humanize_semantic("time", col_map)
h_event <- humanize_semantic("event", col_map)Step 2: Coerce the duration numeric (95% rule); drop invalid durations
tv <- df$time
if (!is.numeric(tv)) {
conv <- suppressWarnings(as.numeric(as.character(tv)))
n_orig <- sum(!is.na(tv) & as.character(tv) != "")
if (n_orig > 0 && sum(!is.na(conv)) >= 0.95 * n_orig) {
tv <- conv
} else {
stop(sprintf(
paste0("The life/duration column(%s) could not be read as a number. ",
"Weibull analysis needs a numeric time, cycle, or distance to failure."),
h_time))
}
}
bad_time <- is.na(tv) | !is.finite(tv) | tv <= 0
dropped_time_rows <- sum(bad_time)Step 3: Binarize the failure indicator (censored = still running)
eb <- .binarize_event(df$event, h_event)
ev <- eb$ev
failure_level <- eb$level
bad_event <- is.na(ev)
dropped_event_rows <- sum(bad_event & !bad_time)
keep <- !bad_time & !bad_event
d <- data.frame(time = tv[keep], event = as.integer(ev[keep]),
stringsAsFactors = FALSE)
final_rows <- nrow(d)
rows_removed <- initial_rows - final_rowsStep 4: Guard the fit — rows, variation, and observed failures
if (final_rows < 20) {
stop(sprintf(
paste0("Only %d usable rows across %s and %s — Weibull life analysis needs ",
"at least 20 units with a positive duration and a readable failure flag."),
final_rows, h_time, h_event))
}
if (length(unique(d$time)) < 2) {
stop(sprintf(
paste0("Every usable value of %s is the same number(%s) — a life ",
"distribution cannot be fitted to a single duration."),
h_time, .fmt_sig(d$time[1])))
}
n_units <- final_rows
n_failures <- sum(d$event)
n_censored <- n_units - n_failures
censored_pct <- 100 * n_censored / n_units
if (n_failures < 3) {
stop(sprintf(
paste0("Only %d failure(s) are recorded in %s — the Weibull shape and scale ",
"cannot both be estimated from fewer than 3 observed failures. The ",
"other %d unit(s) are still running(right-censored)."),
n_failures, h_event, n_censored))
}
if (length(unique(d$time[d$event == 1])) < 2) {
stop(sprintf(
paste0("Every observed failure in %s happened at the same %s value — the ",
"shape parameter cannot be estimated without variation between failures."),
h_event, h_time))
}
few_failures <- n_failures < 20Step 5: MLE Weibull fit with right-censoring, and the parameterization
survreg fits log(T) = mu + sigma * W (W standard extreme value). The reliability convention is R(t) = exp(-(t/eta)^beta), so shape beta = 1 / sigma scale eta = exp(mu)
fit <- tryCatch(
suppressWarnings(survreg(Surv(time, event) ~ 1, data = d, dist = "weibull")),
error = function(e) NULL
)
if (is.null(fit) || !is.finite(fit$scale) || fit$scale <= 0) {
stop(sprintf(
paste0("The Weibull model could not be fitted to %s and %s — the failure ",
"times carry too little information to identify a life distribution."),
h_time, h_event))
}
mu <- unname(coef(fit)[1])
sigma <- as.numeric(fit$scale)
beta <- 1 / sigma
eta <- exp(mu)
if (!is.finite(beta) || !is.finite(eta) || beta <= 0 || eta <= 0) {
stop(sprintf(
"The fitted Weibull parameters for %s were not finite — the life data cannot support this model.",
h_time))
}95% intervals. The fit's covariance is over (mu, log sigma); log(beta) = -log(sigma), so beta and eta transform directly, and every life quantile gets a delta-method interval on the log scale.
V <- tryCatch(vcov(fit), error = function(e) NULL)
ci_available <- !is.null(V) && is.matrix(V) && nrow(V) >= 2 &&
all(is.finite(V[1:2, 1:2])) && V[1, 1] > 0 && V[2, 2] > 0
if (ci_available) {
se_mu <- sqrt(V[1, 1]); se_logsigma <- sqrt(V[2, 2])
beta_lo <- exp(log(beta) - 1.96 * se_logsigma)
beta_hi <- exp(log(beta) + 1.96 * se_logsigma)
eta_lo <- exp(mu - 1.96 * se_mu)
eta_hi <- exp(mu + 1.96 * se_mu)
beta_p <- 2 * pnorm(-abs(log(beta) / se_logsigma))
} else {
beta_lo <- NA_real_; beta_hi <- NA_real_
eta_lo <- NA_real_; eta_hi <- NA_real_
beta_p <- NA_real_
}Life quantile t_p = eta * (-log(1-p))^sigma, with a delta-method interval.
.life_q <- function(p) {
kp <- -log(1 - p)
lt <- mu + sigma * log(kp)
est <- exp(lt)
if (!ci_available) return(c(est = est, lo = NA_real_, hi = NA_real_))
g <- c(1, sigma * log(kp))
v <- as.numeric(t(g) %*% V[1:2, 1:2] %*% g)
se <- sqrt(max(v, 0))
c(est = est, lo = exp(lt - 1.96 * se), hi = exp(lt + 1.96 * se))
}Mean life = eta * Gamma(1 + sigma), same delta-method treatment.
.mean_life <- function() {
lm_ <- mu + lgamma(1 + sigma)
est <- exp(lm_)
if (!ci_available) return(c(est = est, lo = NA_real_, hi = NA_real_))
g <- c(1, digamma(1 + sigma) * sigma)
v <- as.numeric(t(g) %*% V[1:2, 1:2] %*% g)
se <- sqrt(max(v, 0))
c(est = est, lo = exp(lm_ - 1.96 * se), hi = exp(lm_ + 1.96 * se))
}
q10 <- .life_q(0.10); q50 <- .life_q(0.50); qm <- .mean_life()
b10 <- unname(q10["est"]); b10_lo <- unname(q10["lo"]); b10_hi <- unname(q10["hi"])
b50 <- unname(q50["est"]); b50_lo <- unname(q50["lo"]); b50_hi <- unname(q50["hi"])
mttf <- unname(qm["est"]); mttf_lo <- unname(qm["lo"]); mttf_hi <- unname(qm["hi"])Step 6: The failure-mode reading of beta — decided by the CI, not the point
mode_key <- if (!ci_available) {
"indeterminate"
} else if (beta_hi < 1) {
"infant_mortality"
} else if (beta_lo > 1) {
"wear_out"
} else {
"indeterminate"
}
ci_phrase <- if (ci_available) {
paste0("95% CI ", .fmt_sig(beta_lo), " to ", .fmt_sig(beta_hi))
} else {
"no confidence interval could be computed"
}
lean_word <- if (beta > 1) "above 1" else if (beta < 1) "below 1" else "at 1"
mode_label <- switch(mode_key,
infant_mortality = "infant mortality(failure rate falling with age)",
wear_out = "wear-out(failure rate rising with age)",
"not distinguishable from a constant failure rate")
mode_sentence <- switch(mode_key,
infant_mortality = paste0(
"The shape parameter is ", .fmt_sig(beta), " (", ci_phrase,
"), entirely below 1, so the fitted failure rate falls as units age — the ",
"infant-mortality pattern, in which the weak units fail early and the ",
"survivors are more reliable than the population started out. Burn-in or ",
"screening acts on this pattern; scheduled replacement does not."),
wear_out = paste0(
"The shape parameter is ", .fmt_sig(beta), " (", ci_phrase,
"), entirely above 1, so the fitted failure rate rises as units age — the ",
"wear-out pattern, in which older units are progressively more likely to ",
"fail. Age-based preventive replacement acts on this pattern; it does not ",
"help when the rate is flat."),
paste0(
"The shape parameter is ", .fmt_sig(beta), " (", ci_phrase,
"), and that interval spans 1. On this data the failure rate cannot be ",
"told apart from a constant one, so no failure mode is claimed: the point ",
"estimate sits ", lean_word, ", but the evidence does not separate infant ",
"mortality, random failure, and wear-out."))Step 7: Observed follow-up, and the horizon the projection reaches
t_max_obs <- max(d$time)
t_last_failure <- max(d$time[d$event == 1])
t_r05 <- unname(.life_q(0.95)["est"]) # time by which 95% have failed
t_horizon <- min(max(t_r05, 1.25 * t_max_obs), 2.5 * t_max_obs)
if (!is.finite(t_horizon) || t_horizon <= t_max_obs) t_horizon <- 1.25 * t_max_obs
extrapolation_factor <- t_horizon / t_max_obs
.reliability <- function(t) exp(-(t / eta)^beta)
.hazard <- function(t) (beta / eta) * (t / eta)^(beta - 1)
r_at_obs_end <- .reliability(t_max_obs)
r_at_horizon <- .reliability(t_horizon)Step 8: Curves — fitted R(t), the projected segment, and observed KM
t_start <- max(min(d$time) / 2, t_horizon / 500)
grid <- unique(c(0, seq(t_start, t_horizon, length.out = 160)))
in_data <- grid <= t_max_obs
fit_lab <- "Weibull fit(within observed data)"
proj_lab <- "Projection beyond observed data(model assumption)"
km_lab <- "Observed(Kaplan-Meier)"
rel_fit <- data.frame(
time_point = round(grid, 4),
reliability_pct = round(100 * .reliability(grid), 3),
series = ifelse(in_data, fit_lab, proj_lab),
stringsAsFactors = FALSE
)
km_fit <- survfit(Surv(time, event) ~ 1, data = d)
km_s <- summary(km_fit)
km_t <- km_s$time; km_surv <- km_s$surv
keep_km <- which(is.finite(km_t) & is.finite(km_surv))
km_t <- km_t[keep_km]; km_surv <- km_surv[keep_km]
idx <- if (length(km_t) > 140) unique(round(seq(1, length(km_t), length.out = 140))) else seq_along(km_t)
rel_km <- data.frame(
time_point = round(c(0, km_t[idx]), 4),
reliability_pct = round(100 * c(1, km_surv[idx]), 3),
series = km_lab,
stringsAsFactors = FALSE
)
reliability_df <- rbind(rel_fit, rel_km)
rownames(reliability_df) <- NULL
hazard_grid <- seq(t_start, t_horizon, length.out = 160)
hazard_df <- data.frame(
time_point = round(hazard_grid, 4),
hazard_rate = signif(.hazard(hazard_grid), 6),
series = ifelse(hazard_grid <= t_max_obs, fit_lab, proj_lab),
stringsAsFactors = FALSE
)Step 9: Weibull probability plot + fit-quality diagnostics
Plotting positions come from the Kaplan-Meier estimate, which is what makes them valid when units are censored. A Weibull population is a straight line in ln(t) vs ln(-ln R) coordinates, with slope beta.
ok_pp <- which(km_surv > 0 & km_surv < 1 & km_t > 0)
pp_r2 <- NA_real_; km_max_gap_pp <- NA_real_
probplot_df <- data.frame(log_time = numeric(0), log_log_failure = numeric(0),
series = character(0), stringsAsFactors = FALSE)
if (length(ok_pp) >= 3) {
x_all <- log(km_t[ok_pp])
y_all <- log(-log(km_surv[ok_pp]))
good <- is.finite(x_all) & is.finite(y_all)
x_all <- x_all[good]; y_all <- y_all[good]
if (length(x_all) >= 3 && stats::var(x_all) > 0 && stats::var(y_all) > 0) {
pp_r2 <- suppressWarnings(stats::cor(x_all, y_all))^2
sidx <- if (length(x_all) > 400) unique(round(seq(1, length(x_all), length.out = 400))) else seq_along(x_all)
obs_pp <- data.frame(
log_time = round(x_all[sidx], 4),
log_log_failure = round(y_all[sidx], 4),
series = "Observed failures(Kaplan-Meier positions)",
stringsAsFactors = FALSE
)
line_x <- seq(min(x_all), max(x_all), length.out = 60)
line_pp <- data.frame(
log_time = round(line_x, 4),
log_log_failure = round(beta * (line_x - log(eta)), 4),
series = "Fitted Weibull line",
stringsAsFactors = FALSE
)
probplot_df <- rbind(obs_pp, line_pp)
rownames(probplot_df) <- NULL
}How far the fitted curve sits from the observed KM curve, in points
gaps <- abs(100 * km_surv[ok_pp] - 100 * .reliability(km_t[ok_pp]))
gaps <- gaps[is.finite(gaps)]
if (length(gaps) > 0) km_max_gap_pp <- max(gaps)
}
fit_quality_word <- if (is.na(pp_r2)) {
"could not be scored"
} else if (pp_r2 >= 0.98) {
"very close to a straight line"
} else if (pp_r2 >= 0.95) {
"close to a straight line"
} else if (pp_r2 >= 0.90) {
"loosely straight, with visible curvature"
} else {
"clearly curved, which argues against a single Weibull population"
}Step 10: Tables
ident_note <- if (few_failures) {
paste0("Only ", n_failures, " failure(s) were observed, so the shape is ",
"poorly identified — read the interval, not the point estimate.")
} else {
paste0(n_failures, " observed failures carry the shape estimate; the ",
"interval below is the honest range.")
}
params_df <- data.frame(
parameter = c("Shape(beta)",
"Scale(eta) — characteristic life",
"Observed failures carrying the fit",
"Probability-plot R-squared",
"Largest gap to the observed curve(percentage points)"),
estimate = c(round(beta, 4), round(eta, 3), n_failures,
if (is.na(pp_r2)) NA_real_ else round(pp_r2, 4),
if (is.na(km_max_gap_pp)) NA_real_ else round(km_max_gap_pp, 2)),
ci_low = c(if (is.na(beta_lo)) NA_real_ else round(beta_lo, 4),
if (is.na(eta_lo)) NA_real_ else round(eta_lo, 3),
NA_real_, NA_real_, NA_real_),
ci_high = c(if (is.na(beta_hi)) NA_real_ else round(beta_hi, 4),
if (is.na(eta_hi)) NA_real_ else round(eta_hi, 3),
NA_real_, NA_real_, NA_real_),
interpretation = c(
paste0("Below 1 means the failure rate falls with age, above 1 means it rises, ",
"1 means it is flat. Here: ", mode_label, ". ", ident_note),
paste0("The ", h_time, " by which 63.2% of units have failed — the natural ",
"scale of this life distribution, not a safe operating limit."),
paste0("Censored units(", format(n_censored, big.mark = ","),
", ", .fmt_fixed(censored_pct, 1), "% of the file) contribute the ",
"time they survived, but only failures pin down the shape."),
paste0("How straight the data is in Weibull coordinates: ", fit_quality_word,
". This scores the fit, and a good fit is still a model choice, not a fact."),
paste0("The largest distance between the fitted reliability curve and the ",
"assumption-free Kaplan-Meier curve computed from the same data.")
),
stringsAsFactors = FALSE
)
proj_note <- paste0(
"Beyond ", .fmt_sig(t_max_obs), " ", h_time,
" there is no data — this number is the fitted model extended forward, not an observation.")
life_df <- data.frame(
metric = c("B10 life(10% of units failed)",
"B50 life(half of units failed)",
"Characteristic life(63.2% failed)",
"Mean life(MTBF)",
paste0("Reliability at the end of the data(", .fmt_sig(t_max_obs), ")"),
paste0("Reliability projected to ", .fmt_sig(t_horizon))),
estimate = c(round(b10, 3), round(b50, 3), round(eta, 3), round(mttf, 3),
round(100 * r_at_obs_end, 2), round(100 * r_at_horizon, 2)),
ci_low = c(if (is.na(b10_lo)) NA_real_ else round(b10_lo, 3),
if (is.na(b50_lo)) NA_real_ else round(b50_lo, 3),
if (is.na(eta_lo)) NA_real_ else round(eta_lo, 3),
if (is.na(mttf_lo)) NA_real_ else round(mttf_lo, 3),
NA_real_, NA_real_),
ci_high = c(if (is.na(b10_hi)) NA_real_ else round(b10_hi, 3),
if (is.na(b50_hi)) NA_real_ else round(b50_hi, 3),
if (is.na(eta_hi)) NA_real_ else round(eta_hi, 3),
if (is.na(mttf_hi)) NA_real_ else round(mttf_hi, 3),
NA_real_, NA_real_),
interpretation = c(
paste0("The design life engineers usually quote: 10% of units are expected to ",
"have failed by ", .fmt_sig(b10), " ", h_time, "."),
paste0("Half of all units are expected to have failed by ", .fmt_sig(b50), " ", h_time, "."),
paste0("The Weibull scale parameter, always the 63.2% point whatever the shape."),
paste0("The average life implied by the fit(eta times the gamma function of ",
"1 + 1/beta). It is the mean, not a survival guarantee: ",
.fmt_fixed(100 * (1 - .reliability(mttf)), 1),
"% of units are expected to have failed by then."),
paste0("The share still running at the last ", h_time,
" actually covered by this data — the furthest point the file can speak to."),
paste0("The fitted survival share at ", .fmt_sig(t_horizon), " ", h_time, ". ", proj_note)
),
stringsAsFactors = FALSE
)
metrics <- list(
`Units` = n_units,
`Failures` = as.integer(n_failures),
`Censored %` = round(censored_pct, 1),
`Shape(beta)` = round(beta, 3),
`Scale(eta)` = round(eta, 2),
`B10 Life` = round(b10, 2),
`B50 Life` = round(b50, 2),
`Mean Life(MTBF)` = round(mttf, 2),
`Failure Mode` = mode_label
)
cens_phrase <- if (n_censored == 0) {
paste0("every unit in the file failed, so there is no right-censoring to correct for")
} else {
paste0(format(n_censored, big.mark = ","), " of ", format(n_units, big.mark = ","),
" units(", .fmt_fixed(censored_pct, 1),
"%) were still running when observation stopped and are counted as censored")
}
json_output <- list(
answer = paste0(
"Weibull life analysis of ", format(n_units, big.mark = ","), " units with ",
n_failures, " observed failures — ", cens_phrase, ". Fitted shape beta = ",
.fmt_sig(beta),
if (ci_available) paste0(" (95% CI ", .fmt_sig(beta_lo), " to ", .fmt_sig(beta_hi), ")") else "",
", scale eta = ", .fmt_sig(eta), " ", h_time,
", which reads as ", mode_label, ". B10 life ", .fmt_sig(b10), ", B50 life ",
.fmt_sig(b50), ", mean life ", .fmt_sig(mttf), " ", h_time, ". ",
"The data itself only reaches ", .fmt_sig(t_max_obs), " ", h_time,
"; the reliability curve is projected to ", .fmt_sig(t_horizon),
" (", .fmt_fixed(extrapolation_factor, 2),
" times the observed follow-up), and everything past the end of the data is ",
"the fitted model rather than evidence.",
if (few_failures) paste0(" With only ", n_failures,
" failures the shape is poorly identified — read the interval, not the point estimate.") else ""
),
cards = lapply(
c("tldr", "overview", "preprocessing", "weibull_parameters",
"reliability_curve", "life_metrics", "hazard_curve", "probability_plot"),
function(cid) list(id = cid, metrics = metrics)
)
)
list(
initial_rows = initial_rows, final_rows = final_rows, rows_removed = rows_removed,
h_time = h_time, h_event = h_event, failure_level = failure_level,
n_units = n_units, n_failures = n_failures, n_censored = n_censored,
censored_pct = censored_pct, few_failures = few_failures,
dropped_time_rows = dropped_time_rows, dropped_event_rows = dropped_event_rows,
beta = beta, beta_lo = beta_lo, beta_hi = beta_hi, beta_p = beta_p,
eta = eta, eta_lo = eta_lo, eta_hi = eta_hi, sigma = sigma, mu = mu,
ci_available = ci_available,
mode_key = mode_key, mode_label = mode_label, mode_sentence = mode_sentence,
b10 = b10, b10_lo = b10_lo, b10_hi = b10_hi,
b50 = b50, b50_lo = b50_lo, b50_hi = b50_hi,
mttf = mttf, mttf_lo = mttf_lo, mttf_hi = mttf_hi,
t_max_obs = t_max_obs, t_last_failure = t_last_failure, t_horizon = t_horizon,
extrapolation_factor = extrapolation_factor,
r_at_obs_end = r_at_obs_end, r_at_horizon = r_at_horizon,
pp_r2 = pp_r2, km_max_gap_pp = km_max_gap_pp, fit_quality_word = fit_quality_word,
reliability_df = reliability_df, hazard_df = hazard_df, probplot_df = probplot_df,
params_df = params_df, life_df = life_df,
metrics = metrics, json_output = json_output
)
}