含未婚样本的首次结婚年龄中位数估计及置信区间求解方法

嘿，这个问题很常见——你的样本里有“从未结婚”的个体，这属于右删失数据（我们不知道这些人未来会不会结婚，但至少到调查时他们还没结婚，意味着他们的首次结婚年龄肯定大于当前的年龄）。直接算普通中位数行不通，得用生存分析里的方法来搞，我给你一步步说清楚：

1. 先把数据理清楚 #

首先你得明确：那些“Never married”的个体，你得知道他们的当前年龄（这是删失的时间点）。如果没给的话，只能做保守假设（比如假设他们当前年龄等于样本里最大的结婚年龄38岁，或者更大）。

拿你的样本举例：12个个体里，10个有明确结婚年龄，2个是未结婚。假设这两个未结婚的人当前年龄是35和40，整理后的数据就是：

• 结婚年龄（事件发生）：20,20,21,22,24,26,29,30,32,38
• 删失年龄（未结婚）：35,40

2. 用Kaplan-Meier法估计生存函数，找中位数 #

这里的生存函数S(t)定义是：到年龄t时还没结婚的人的比例。我们要找的中位数，就是最小的年龄t，让S(t)=0.5——也就是刚好一半的人已经结婚的年龄。

步骤很简单：

1. 把所有年龄（结婚+删失）从小到大排序：20,20,21,22,24,26,29,30,32,35,38,40
2. 一步步算生存概率：
   • 一开始所有人都没结婚，S(0)=1.0
   • 到20岁：2个人结婚，风险集是12人 → S(20)=1*(12-2)/12≈0.833（83.3%的人还没结婚）
   • 到21岁：又1人结婚，风险集剩10人 → S(21)=0.833*(10-1)/10≈0.750
   • 到22岁：再1人结婚，风险集9人 → S(22)=0.750*(9-1)/9≈0.667
   • 到24岁：1人结婚，风险集8人 → S(24)=0.667*(8-1)/8≈0.583
   • 到26岁：1人结婚，风险集7人 → S(26)=0.583*(7-1)/7≈0.500

哎，刚好到26岁时，S(t)=0.5，所以中位数就是26岁。如果删失年龄不同，比如删失年龄是25和27，那生存概率会在26岁时降到0.5以下，中位数还是26岁（因为是第一个让S(t)≤0.5的年龄）。

3. 计算中位数的置信区间 #

有两种常用方法，适合不同场景：

方法一：对数变换法（适合大样本） #

这是统计软件默认的方法，步骤是：

1. 用Greenwood公式计算生存函数在中位数处的方差
2. 用正态近似，把生存函数的置信区间转换为年龄的置信区间。简单说，就是找两个年龄t_low和t_high，让S(t_low)=exp(-exp(ln(-ln(0.5)) + 1.96*SE))，S(t_high)=exp(-exp(ln(-ln(0.5)) - 1.96*SE))（95%置信度用1.96）。

方法二：Bootstrap百分位数法（适合小样本） #

你的样本只有12个，这个方法更靠谱：

1. 从原始样本（包含结婚年龄和删失数据）中有放回地抽12个样本，重复1000-2000次
2. 每个重复样本都用Kaplan-Meier算中位数
3. 把这些中位数排序，取第2.5%和97.5%的分位数，就是95%置信区间。

4. 用代码实操（以R为例） #

用survival包就能搞定，代码很简单：

# 加载生存分析包
library(survival)

# 整理数据：time是年龄，status=1表示已结婚，status=0表示未结婚（删失）
time <- c(21,20,35,38,24,30,40,22,20,32,26,29)
status <- c(1,1,0,1,1,1,0,1,1,1,1,1)

# 拟合Kaplan-Meier模型
km_model <- survfit(Surv(time, status) ~ 1)

# 查看结果，包含中位数和置信区间
summary(km_model)

运行后，你会看到中位数以及对应的95%置信区间，默认用的是对数变换法。如果想用Bootstrap，也可以用boot包来实现。

关键提醒 #

• 一定要有未结婚个体的当前年龄！如果没有，就做敏感性分析——比如假设他们当前年龄是38岁、40岁，看看中位数的变化范围。
• 小样本情况下，Bootstrap法的置信区间比正态近似更可信。

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