如何在R中按亚组计算带置信区间的年龄标准化率？

我有如下格式的数据框：

df <- data.frame ( 
  time = rep(c("2010", "2011", "2012", "2013", "2014"),4), 
  age = rep(c("40-44", "45-49", "50-54", "55-59", "60-64"),4), 
  weight = rep(c(0.38, 0.23, 0.19, 0.12, 0.08),4), 
  ethnic = rep(c(rep("M",5),rep("NM",5)),2), 
  gender = c(rep("M",10), rep("F",10)), 
  pop = round((runif(10, min = 10000, max = 99999)), digits = 0), 
  count = round((runif(10, min = 100, max = 999)), digits = 0) 
) 
df$rate = df$count / df$pop 

我希望计算直接年龄标准化发病率（发病率=count/pop）及其置信区间，且按time、gender、ethnic的所有组合亚组分别计算（注：age作为年龄分组用于标准化）。请问能否在R中实现？

我曾尝试使用epitools包中的ageadjust.direct函数：

age_adjust_test <- ageadjust.direct(count = df$count, pop = df$pop, rate = df$rate, stdpop = df$weight)

但该函数仅输出整体调整率、置信区间和粗率，请问如何按亚组获取上述结果？

当然可以实现！我们可以结合dplyr的分组功能和epitools::ageadjust.direct函数，对每个亚组单独计算标准化率及置信区间。具体步骤如下：

1. 先加载所需的包：

library(dplyr)
library(epitools)

2. 对数据按time、gender、ethnic进行分组，然后对每个组应用ageadjust.direct函数，并提取关键结果：

# 注意：务必保证每个亚组的年龄组顺序与标准权重的顺序一致
age_adjust_subgroups <- df %>%
  group_by(time, gender, ethnic) %>%
  summarize(
    # 提取直接标准化率
    adj_rate = ageadjust.direct(count = count, pop = pop, rate = rate, stdpop = weight)$rate.std,
    # 提取95%置信区间下限
    adj_rate_lower = ageadjust.direct(count = count, pop = pop, rate = rate, stdpop = weight)$lower,
    # 提取95%置信区间上限
    adj_rate_upper = ageadjust.direct(count = count, pop = pop, rate = rate, stdpop = weight)$upper,
    # 提取粗率
    crude_rate = ageadjust.direct(count = count, pop = pop, rate = rate, stdpop = weight)$rate.crude,
    .groups = "drop"  # 取消分组状态，得到普通数据框
  )

3. 查看分组后的标准化结果：

print(age_adjust_subgroups)

关键注意事项 #

• 这里的stdpop = weight使用了你数据中的weight列作为标准人口权重，必须保证所有亚组的年龄组顺序与标准权重的顺序完全匹配，否则标准化结果会出错。
• 如果你需要使用实际的标准人口数（而非权重比例），只需将stdpop替换为对应的标准人口向量即可，同样要注意年龄组顺序一致。
• ageadjust.direct返回的是一个包含多项结果的列表，我们通过$符号提取其中的标准化率、置信区间和粗率，整合到分组后的结果数据框中。

这样就能轻松得到每个time、gender、ethnic组合亚组的直接年龄标准化发病率及其置信区间啦！

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