本文介绍了如何使用R语言的Nordpred和BAPC包，基于GBD数据库对未来疾病负担进行预测。作者通过复现相关研究，详细展示了数据准备、模型建立与运算的过程，包括APC模型的构建、数据整理、模型评价等，并分享了在模型应用中遇到的问题和疑惑。
最近看到一篇文献，利用GBD数据库中的数据对疾病负担的发展趋势进行预测。觉得很有趣，就抽空用R语言复盘了一下里面关键的方法。
参考文献：
Changing trends in the disease burden of esophageal cancer in China from 1990 to 2017 and its predicted level in 25 years
DOI: 10.1002/cam4.3775
1. 相关R包的安装
主要用到的R包有3个，一个一个讲：
1.1 Nordpred package
用于做 Nordpred age-period-cohort (APC) analysis
注意，这个包不在r-cran上，需要利用remotes包从github上下载
install.packages("remotes")
library(remotes)
remotes::install_github("haraldwf/nordpred")
library(nordpred)
如果下载比较慢，可以下载我存在百度云盘上的包
链接：https://pan.baidu.com/s/1DmooS88GtWVRflcelJ71sQ
提取码：jm4q
接着本地安装一下
install.packages("F:/R_packages/haraldwf-nordpred-ef83b2c.tar.gz", repos = NULL, type = "source")
1.2 INLA package
为了利用Integrated nested Laplace approximation (INLA) approach
install.packages("INLA", repos = "https://inla.r-inla-download.org/R/stable", dependencies = TRUE)
这个直接下比较麻烦，我试了好几次都半路掉线失败了，用我分享的包吧
#安装INLA前需要预装一个包 foreach
install.packages('foreach')
install.packages("F:/R_packages/INLA_21.02.23.zip", repos = NULL, type = "source")
1.3 BAPC package
为了利用Bayesian age-period-cohort models，可以产生 age-specific and age-standardized projected rates
#安装BAPC前需要预装一些包
install.packages("cmprsk")
install.packages("fanplot")
install.packages("Epi")
#接着下载安装BAPC包
install.packages("BAPC", repos = "http://R-Forge.R-project.org")
library(BAPC)
提醒：常见报错有rlang包版本低，重新卸了再装即可：
remove.packages('rlang')
install.packages('rlang')
同样，如果BAPC包下载不顺利，可以用我云盘上存的包进行本地安装：
install.packages("F:/R_packages/BAPC_0.0.34.tar.gz", repos = NULL, type = "source")
2. 数据准备
（1） GBD数据下载：http://ghdx.healthdata.org/gbdresults-tool
这里的案例用的是cancer medicine上发的中国食管癌文章，选用GBD 2017的数据；目前更新到GBD 2019的版本；在官网上按需要的条件一选，然后下载csv即可；
（2） 人口预测数据：https://population.un.org/wpp/Download/Standard/CSV/
（3） 人口结构标化，WHO2000-2025标准：https://seer.cancer.gov/stdpopulations/world.who.html
3. 模型讲解（可略过）
以下是我个人的理解，因为我只是个临床大夫，对模型的理解可能不够准确，所以这一段可以略过。
nordpred包和BAPC包的模型基础都是 **Age-Period-Cohort (APC)**模型，这个模型的理论基础是，发病率或死亡率和年龄结构和人口规模有关联。挖掘这种关联的基础还是广义线性模型（GLM）。nordpred包中有两种GLM：
- power model即指数模型，当时把这个模型应用于癌症预测的时候（Moller et al.），认为指数为5的时候最好用，nordpred包里也采用了指数为5的固定指数模型；
- Poisson log-link model，这个也是比较经典的了，log变换后满足泊松分布；
在上述基础上，APC模型可以在传统概率论框架下搭建，也可以在贝叶斯框架下进行，后者即为Bayesian age-period-cohort models（BAPC）。
有文献提到了一种在贝叶斯框架下，以随机指数模型作为链接函数的APC模型，讲得很仔细，具体的方程和评价方法都有，建议参考：
Verena Jürgens et al. A Bayesian generalized age–period–cohort power model for cancer projections；
DOI: 10.1002/sim.6248
INLA：integrated nested Laplace approximation (INLA) approach
贝叶斯推断中有一种比较有名的方法：马尔可夫链，但是这种方法效率不高，比较慢。而ILNA方法作为马尔科夫链的低计算密度的替代，用于隐高斯模型中进行近似贝叶斯推断。BAPC包中提供了在INLA框架下搭建的APC模型，官网描述：Package implementing Bayesian age-period-cohort models with the focus on projections. BAPC uses integrated nested Laplace approximations (INLA) for full Bayesian inference. BAPC generates age-specific and age-standardized projected rates. When interest lies in the predictive distribution, Poisson noise is automatically added.
4. 模型建立与运算
4.1 nordpred构建APC模型
nordpred的包详解见：https://rdrr.io/github/haraldwf/nordpred/man/nordpred.html
其中还自带一个示例，可以粘贴代码直接运行
4.1.1 数据整理
(1) 食管癌的数据，各个年龄段，1990-2017年，中国（这里以男性的发病率为例）
注意：
- nordpred默认是18个年龄组，不符合的话会报错。由于食管癌15岁以下的没数据或都是为0，所以整成一组了。
- nordpred包的模型最多预测5个值，如果你输入1990-2017年每年的数据，就可以预测2018-2022的值；所以，为了预测2042年的发病情况，这里把1990-2017年的数据按5年一组（1990-1992是3年一组），求了平均数；
esoph <- fread('F:/nordpred/esoph.csv')
unique(esoph$age_name)
ages <- c("15 to 19", "20 to 24", "25 to 29",
"30 to 34", "35 to 39", "40 to 44", "45 to 49", "50 to 54", "55 to 59",
"60 to 64", "65 to 69", "70 to 74", "75 to 79", "80 to 84", "85 to 89",
"90 to 94", "95 plus")
### for incidence for male
esoph_ages <- esoph[age_name %in% ages &
sex_name == 'Male' &
metric_name == 'Number' &
measure_name == 'Incidence',
.(age_id, age_name, year, val)]
#esoph_ages <- esoph_ages[order(age_id),]
esoph_ages_n <- dcast(data = esoph_ages, age_id + age_name ~ year)
rownames(esoph_ages_n) <- c("15-19", "20-24", "25-29", "30-34",
"35-39", "40-44", "45-49", "50-54",
"55-59", "60-64", "65-69", "70-74",
"75-79", "80-84", "85-89", "90-94","95+")
esoph_ages_n["0-14",] <- 0
esoph_ages_n <- esoph_ages_n[order(esoph_ages_n$age_id),]
#这里需要按5年一组整理发病数， 我取了平均；因为有很多重复工作，所以写了个函数来完成，函数代码见文末
#参数分别是要整理的表格名，开始年份，结束年份，作为观察数据终点的年份
esoph_ages_g <- function_year5(esoph_ages_n , 1990, 2017, 2017)
整理后效果：
(2) 中国1990-2042的人口结构整理
注意：
- 1990-2042年的数据按5年一组（1990-1992是3年一组），求了平均数；
- 源数据的单位是每1000人，所以提取的时候*1000
population <- fread('F:/nordpred/WPP2019_PopulationByAgeSex_Medium.csv')
china_population <- population[Location == 'China',]
unique(china_population$AgeGrp)
china_population_1990_2017 <- china_population[Time %in% 1990:2042,
.(Age_id = 1:21, AgeGrp, Time, PopMale = PopMale*1000)]
china_population_1990_2017_n <- dcast(data = china_population_1990_2017, Age_id + AgeGrp ~ Time)
#calculation of 95+
china_population_1990_2017_n[22, 3:55] <- china_population_1990_2017_n[20, 3:48] + china_population_1990_2017_n[21, 3:55]
china_population_1990_2017_n[1, 3:55] <- china_population_1990_2017_n[1, 3:48] + china_population_1990_2017_n[2, 3:55] + china_population_1990_2017_n[3, 3:55]
china_population_1990_2017_n <- china_population_1990_2017_n[-c(2,3,20,21),]
rownames(china_population_1990_2017_n) <- c("0-14", "15-19", "20-24", "25-29", "30-34",
"35-39", "40-44", "45-49", "50-54",
"55-59", "60-64", "65-69", "70-74",
"75-79", "80-84", "85-89", "90-94","95+")
china_population_1990_2017_g <- function_year5(china_population_1990_2017_n, 1990, 2042, 2017)
整理后的效果：
(3) 世界标准化人口
## data for age-standardization
age_stand <- read_file('F:/nordpred/AS.txt')
age_stand <- gsub(pattern = '\n', replacement = '\t', x = age_stand)
age_stand <- gsub(
最低0.47元/天 解锁文章
扫一扫
946
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
