专注系列化、高质量的R语言教程
本文接上篇继续介绍人群归因分数(PAF)的计算,主要介绍风险因子是多分类变量、含控制变量两种情况。
本篇目录如下:
1 多分类变量
2 含控制变量
2.1 错误计算方法
2.2 正确计算方法1
2.3 正确计算方法2(条件PAF)
2.4 正确计算方法3
3 总结
1 多分类变量
上篇举的例子中,风险因子phys
是二分变量。有读者问到风险因子是多分类变量时该如何计算PAF。
示例数据如下:
library(causalPAF)
library(tidyverse)
data <- strokedata %>%
mutate(phys = as.numeric(phys) - 1,
ahei3tert = as.numeric(ahei3tert) - 1)
ahei3tert
变量表示样本饮食习惯的AHEI得分,原数据经过上述转换后,取值范围为0、1、2,数值越大表示饮食习惯越不健康。
样本分布情况:
table <- data %>%
group_by(ahei3tert) %>%
summarise(pop = length(case),
pop_prop = pop/16623,
case_num = sum(case),
case_prop = mean(case))
table
## ahei3tert pop pop_prop case_num case_prop
## <dbl> <int> <dbl> <dbl> <dbl>
## 1 0 5114 0.308 2457 0.480
## 2 1 5824 0.350 2844 0.488
## 3 2 5685 0.342 2934 0.516
多分类情况下,PAF的计算公式仍然如下,只不过分母相加的项变多了:
## 计算相对风险度
(rr <- table$case_prop/table$case_prop[1])
## [1] 1.000000 1.016398 1.074200
## 人群暴露水平分布
(p <- table$pop_prop)
## [1] 0.3076460 0.3503579 0.3419960
## 计算PAF
1 - rr[1]/(sum(rr * p))
## [1] 0.03018201
2 含控制变量
控制变量一般刻画的是样本个体或家庭的特征,如年龄、性别、受教育程度、经济状况等,它们可能会间接影响样本的患病风险,通常不看作是风险因子,而是混杂因素(confounders)。
esex
变量表示样本性别,1表示女性,2表示男性。在本例中,将其看作控制变量。
2.1 错误计算方法
以phys
、esex
为分组变量,展示样本的分布情况:
table <- data %>%
group_by(esex, phys) %>%
summarise(pop = length(case),
pop_prop = pop/16623,
case_num = sum(case),
case_prop = mean(case))
table
## esex phys pop pop_prop case_num case_prop
## <int> <dbl> <int> <dbl> <dbl> <dbl>
## 1 1 0 759 0.0457 344 0.453
## 2 1 1 5985 0.360 2973 0.497
## 3 2 0 1538 0.0925 699 0.454
## 4 2 1 8341 0.502 4219 0.506
将上表第一行的情况看作理想暴露水平(患病人口比例最小),然后计算PAF:
rr <- table$case_prop/table$case_prop[1]
p <- table$pop_prop
1 - rr[1]/sum(rr * p)
## [1] 0.08512351
但这种计算方法是错误的。