如何用R的fitdistrplus包为右删失数据直方图叠加Weibull PDF

1. 确认删失数据格式 #

fitdistcens要求删失数据为包含left和right列的数据框：右删失的观测，left设为删失阈值，right设为Inf。先模拟一份符合要求的示例数据：

set.seed(123)
# 生成原始Weibull分布数据
true_shape <- 2
true_scale <- 5
n <- 100
raw_data <- rweibull(n, shape = true_shape, scale = true_scale)
# 设定右删失阈值，大于8的观测标记为删失
cens_threshold <- 8
df <- data.frame(
  left = ifelse(raw_data <= cens_threshold, raw_data, cens_threshold),
  right = ifelse(raw_data <= cens_threshold, raw_data, Inf)
)

2. 拟合删失数据的Weibull分布 #

用fitdistcens完成分布拟合，参数会自动估计：

library(fitdistrplus)
fit_cens <- fitdistcens(df, distr = "weibull")

3. 手动绘制直方图+拟合PDF #

plotdistcens默认输出的是经验分布与拟合分布的对比图，并非你需要的直方图+PDF形式，因此需要手动绘图：

• 直方图必须使用密度刻度（freq = FALSE），确保y轴与PDF的取值范围匹配；
• 生成拟合分布的PDF序列，用lines叠加到直方图上。

代码示例：

# 绘制非删失数据的直方图（密度刻度）
non_cens_data <- raw_data[raw_data <= cens_threshold]
hist(non_cens_data, breaks = "FD", freq = FALSE,
     col = "lightblue", border = "white",
     main = "右删失数据直方图 + 拟合Weibull PDF",
     xlab = "观测值", xlim = c(0, cens_threshold + 2))

# 标记删失组区间（用空心柱展示）
cens_count <- sum(raw_data > cens_threshold)
bin_width <- diff(hist(non_cens_data)$breaks)[1]
rect(xleft = cens_threshold, xright = cens_threshold + 2,
     ybottom = 0, ytop = cens_count/(n*bin_width),
     col = NA, border = "red", lwd = 2)

# 生成拟合PDF的x轴序列
x_seq <- seq(0, cens_threshold + 2, length.out = 1000)
# 计算拟合Weibull分布的PDF值
y_pdf <- dweibull(x_seq, shape = fit_cens$estimate["shape"], 
                  scale = fit_cens$estimate["scale"])

# 叠加PDF曲线
lines(x_seq, y_pdf, col = "darkred", lwd = 2)

# 添加图例
legend("topright", legend = c("非删失数据", "删失数据", "拟合Weibull PDF"),
       col = c("lightblue", "red", "darkred"), lty = c(1,1,1), lwd = c(10,2,2))

关键注意事项 #

• 必须设置freq = FALSE：直方图y轴为密度，才能和PDF的取值范围匹配，实现正确叠加；
• fitdistcens返回的estimate中是Weibull分布的形状和尺度参数，直接传入dweibull即可计算PDF；
• 用空心柱单独标记删失区间，能更直观展示删失数据的分布特征。

内容的提问来源于stack exchange，提问作者fre1990