「转」真实世界研究公开课高分文章复现代码_20230825
分享文章
方法学关键词
IPTW,SMD,weight-KM,multiple imputation,subgroup analysis,p for interaction
课程目录（2023年8月版本）
第一周（真实世界介绍）
1.真实世界研究的概念与高分文章简析 2.真实世界研究因果估计方法学介绍（包括高级生存分析方法） 3.整洁数据的一般要求及数据整理（tidydata） 4.学习R语言的必要性及准备工作（环境构建及包安装）
第二周（R语言基础及实用技能）
5.R语言基础（1）-一维变量 6.R语言基础（2）-二维变量及高维变量 7.函数的一般用法介绍 8.出现报错的解决方案
第三周（单因素批量分析及变量筛选）
9.tidyverse体系初探 10.组间基线表格快速制作及统计学考量 11.单因素分析批量实现及统计学考虑 12.连续变量截断值选择的方法总结（包括生存资料） 13.关键因素鉴定得一些方法（P值法，机器学习方法，效应值改变法等）
第四周（多因素分析-协变量调整）
14.有向无环图的应用及快速实现 15.多因素分析批量实现及调整协变量的敏感性分析 16.倾向性评分匹配实现及实战注意事项（PSM） 17.逆概率加权的实现及解释 18.三组的倾向性评分的实现及解释
第五周（探索性数据分析）
19.探索性数据分析及其注意事项 20.缺失值处理的策略及实战 21.未测量混杂因素的统计学测量-E值 22.亚组分析及其森林图的快速绘制（包括交互作用及趋势性检验）
第六周（生存分析）
23.生存资料的整理及常规生存分析方法（KM,COX,生存曲线及累计风险曲线绘制） 24.高级生存分析方法实战（landmark分析，RMST,加权生存分析） 25.纵向数据分析的一些方法 26.win ratio的概念及实现，结果解读
第七周（二分类预测模型得构建）
27.二分类结果变量的预测模型构建（一） 28.二分类结果变量的预测模型构建（二） 29.二分类结果变量的预测模型构建（三）
第八周（生存模型得构建）
30.生存结果变量的预测模型构建（一） 31.生存结果变量的预测模型构建（二） 32.生存结果变量的预测模型构建（三）
第9周
33.机器学习及人工智能概览（基于mlr3）
价格及报名流程
1.支付方式为淘宝店铺，一共9周的课程，整体流程如下：第一步.同学拍下淘宝并付款;第二步.学习前3周的课程，没问题后店铺发货同学确认收货后继续学习后续课程；第三步，如果不满意3周内可以申请退款。 2.课程最终定价2500元，9月底开始，总共9周课程，至少30小时。- 3.课程内容包括答疑部分，答疑内容主要是课程涉及方法学及相关文章及代码排错问题。
具体课程内容及预报名可加以下微信咨询
代码分享
说明
1.代码仅供测试上课说明使用，每项的内容会在课程中仔细说明以提供发表级图片与表格的自动化生成方案。 2.多重插补后多个数据集的IPTW及其后续计算比较复杂，本次不做演示，第三节后内容为单数据集演示。
0.加载相关包
library(tidyverse)#数据科学
library(jskm)#加权KM
library(survival)#生存分基础包
library(WeightIt)#IPTW计算权重
library(survey)#复杂抽样权重
library(cobalt)#评估SMD
library(gtsummary)#快速基线表格
library(visdat)#EDA
library(SmartEDA)#EDA
library(mice)#多重插补
library(MatchThem)#多重插补加权
library(tableone)#SMD表格
library(openxlsx)#读取xlsx格式文件
library(Publish)#亚组分析及交互作用
library(EValue)#E值计算
1.载入数据及探索性数据分析
rm(list = ls())
#数据分布情况
df <- read.xlsx("df.xlsx")
glimpse(df)
#初步可视化
#整体游览数据visdata
vis_miss(df,sort_miss = T)
vis_dat(df,sort_type = T)
vis_cor(df[,-3])
vis_guess(df)#对于数据中存在的异常值检测
#相对密度分布情况
vis_value(select(df, where(is.numeric)),viridis_option = "A")
#相对定量分析
# smartEDA
## overview of the data;
ExpData(data=df,type=1)
## structure of the data
ExpData(data=df,type=2)
names(df)
## Summary statistics by – category
ExpNumStat(df,by="GA",gp="rx",Qnt=seq(0,1,0.1),MesofShape=2,Outlier=TRUE,round=2)
## 分类比例情况
cat.table<- ExpCTable(df)
#整体报告
ExpReport(df,
Target="rx",
label=NULL,
op_file="test.html",
op_dir=getwd(),
sc=2,
sn=2,
Rc="Yes")
#转换数据
df$rx <- ifelse(df$rx=="Obs",0,1)
df$sex <- as.factor(df$sex)
# 2. mice包进行多重插补及组合权重进入多重插补----------------------------------------------------------------------
#默认参数插补数据集
df_imp <- mice(df,m = 10)
df_imp
test <- complete(df_imp,action = 2)
vis_miss(test)
#计算插补后权重权重
df_imp_weighted <- weightthem(rx ~ obstruct+adhere+age+differ+sex+surg,
datasets = df_imp,
approach = 'within',
method = 'ps',
estimand = "ATT")
#展示加权
bal.tab(df_imp_weighted,abs = T)
love.plot(df_imp_weighted,abs = T)
3. SMD基线表绘制（以下演示采用单数据完整数据集情况）
# tableone包进行分析
vars <- dput(names(df))
df <- na.omit(df)
#计算插补后权重权重
df_weighted <- weightit(rx ~ obstruct+adhere+age+differ+sex+surg,
data = df,
approach = 'within',
method = 'ps',
estimand = "ATT")
(t1<- tableone::CreateTableOne(vars = vars ,
data = df,
strata="rx",
smd = T))
ExtractSmd(t1)
#gtsummary联合cobalt
tbl_summary_ex2 <-
df %>%
tbl_summary(
by = rx,
label = list(age ~ "Patient Age"),
statistic = list(all_continuous() ~ "{mean} ({sd})"),
digits = list(age ~ c(0, 1))
) |>
add_p()
tbl_summary_ex2
tbl_summary_ex2 |>
as_flex_table() |>
flextable::save_as_docx(path = "测试.docx")
love.plot(df_weighted,abs = T)
4. IPTW-KM
#计算插补后权重权重
df_weighted <- weightit(rx ~ obstruct+adhere+age+differ+sex+surg,
data = df,
method = "glm")
# 创建一个复杂抽样设计对象，用于进行加权分析
df_svy <- svydesign(
id = ~1, # 指定一个虚拟的标识符变量，因为我们只需要一个观察一个群组
weights = ~df_weighted$weights, # 指定权重变量为
data = df, # 数据集
)
# 使用 svykm() 函数计算加权的 Kaplan-Meier 生存曲线，按性别分组
s1 <- svykm(Surv(time, status > 0) ~ rx, design = df_svy)
svycoxph(Surv(time, status) ~ rx, design = df_svy) |>
broom::tidy(exponentiate = TRUE)
coxph(Surv(time, status) ~ rx, data = df,weights = df_weighted$weights) |>
broom::tidy(exponentiate = TRUE)
# 使用 svyjskm() 函数绘制加权的 Kaplan-Meier 生存曲线
svyjskm(s1,pval = T,table = T)
5. 亚组分析及交互作用
#拟合单因素cox回归方程
fit_cox <- coxph(Surv(time,status)~rx,
data=df,
weights = df_weighted$weights)
#进行亚组分析
df$rx <- as.factor(df$rx)
df$etype <- as.factor(df$etype)
sub_cox <- subgroupAnalysis(object = fit_cox,#单因素回归方程
data = df,#数据集
treatment="rx", #治疗因素（主要研究自变量）
subgroups=c("sex","etype")#需要进行亚组分析的变量
)
sub_cox
plot(sub_cox)
6. E值计算
coxph(Surv(time, status) ~ rx,
data = df,
weights = df_weighted$weights) |>
tbl_regression(exponentiate = T)
svycoxph(Surv(time, status) ~ rx, design = df_svy) |>
tbl_regression(exponentiate = T)
evalues.HR(est = 0.73,
lo = 0.64,
hi = 0.84,
rare = F)