在临床研究统计分析中,临床基线表对于展示研究数据的信息和结构尤为重要。
本文以R语言中自带临床数据集lung为例,展示基本的COX回归与临床基线表的制作。
首先是数据的生存分析
#获取R自带数据集lung
lung <- lung
#查看数据的基本结构
head(lung)
# 如下所示:
# inst time status age sex ph.ecog ph.karno pat.karno
# 1 3 306 2 74 1 1 90 100
# 2 3 455 2 68 1 0 90 90
# 3 3 1010 1 56 1 0 90 90
# 4 5 210 2 57 1 1 90 60
# 5 1 883 2 60 1 0 100 90
# 6 12 1022 1 74 1 1 50 80
# meal.cal wt.loss
# 1 1175 NA
# 2 1225 15
# 3 NA 15
# 4 1150 11
# 5 NA 0
# 6 513 0
#----------------------生存分析(Kaplan-Meier曲线)----------------------#
# 基本流程:
# 1)Surv()提取生存时间
# 2)如果不分组,就用survfit()拟合模型
# 3)如果有分组,就用survdiff()分析差异
# 1)Surv()提取生存时间
Lusurv <- Surv(time = lung$time,event = lung$status)
# 2)不分组,就用survfit()拟合模型并作图
Lufit <- survfit(Lusurv~1)
plot(Lufit)
# 3)有分组,就用survdiff()分析分组差异再作图
Lufit <- survfit(Lusurv~lung$sex)
plot(Lufit,conf.int = "none",
col = c("red","blue"),
lwd = 2,xlab = "Time(month)",ylab = "Survival Probability",
mark.time = T)
abline(h=0.5,lty=3) #中位生存时间线
legend("bottomleft",c("Male","Female"),
col = c("red","blue"),lwd = 2)
# 查看不同性别生存时间是否差异显著
survdiff(Lusurv~lung$sex)
#如下:
# Call:
# survdiff(formula = Lusurv ~ lung$sex)
#
# N Observed Expected (O-E)^2/E (O-E)^2/V
# lung$sex=1 138 112 91.6 4.55 10.3
# lung$sex=2 90 53 73.4 5.68 10.3
#
# Chisq= 10.3 on 1 degrees of freedom, p= 0.001
#可以看出p= 0.001 ,不同性别生存时间有显著差异
然后是单因素COX回归
#-----------------单因素COX回归-----------------#
# 提取生存时间
Lusurv <- Surv(time = lung$time,event = lung$status)
# 将生存时间加入原先的lung数据(dataframe)
lung$Lusurv <- with(lung,Lusurv)
#单因素cox回归(性别)
cox_sex <- coxph(Lusurv~sex,data = lung)
# 我们需要的是风险比,95%置信区间和pvalue
cox_sexSum <- summary(cox_sex)
cox_sexSum$coefficients # exp(coef) 即风险比
# coef exp(coef) se(coef) z Pr(>|z|)
# sex -0.531 0.588 0.167 -3.18 0.00149
cox_sexSum$conf.int # lower .95 upper .95 即置信区间
# exp(coef) exp(-coef) lower .95 upper .95
# sex 0.588 1.7 0.424 0.816
# paste0d的collapse参数,可以把这些字符串拼成一个长字符串,而不是放在一个向量中
HR <- round(cox_sexSum$coefficients[,2],2)
PValue <- round(cox_sexSum$coefficients[,5],3) #保留3位
CI <- paste0(round(cox_sexSum$conf.int[,3:4],2),collapse = '-')
Unicox <- data.frame("Characteristics" = "sex",
"Hazard Ratio" = HR,
"CI95" = CI,
"P Value" = PValue)
#----------自定义制作表格函数Unicoxtable()---------#
Unicoxtable <- function(x){
FML <- as.formula(paste0("Lusurv~",x))
cox_sex <- coxph(FML,data = lung)
cox_sexSum <- summary(cox_sex)
HR <- round(cox_sexSum$coefficients[,2],2)
PValue <- round(cox_sexSum$coefficients[,5],3) #保留3位
CI <- paste0(round(cox_sexSum$conf.int[,3:4],2),collapse = '-')
Unicox <- data.frame("Characteristics" = x,
"Hazard Ratio" = HR,
"CI95" = CI,
"P Value" = PValue)
return(Unicox)
}
#将变量名输入即可,例如age
Unicoxtable("age")
#结果如下:
# Characteristics Hazard.Ratio CI95 P.Value
# 1 age 1.02 1-1.04 0.042
#---------多个变量同时处理--------#
library(plyr)
library(xlsx)
Varname <- colnames(lung)[-c(2:3,11)] #提取部分变量名
Univar <- lapply(Varname, Unicoxtable)
Univar <- ldply(Univar,data.frame)
# 将pvalue小于0.001的特殊显示
Univar$P.Value[Univar$P.Value == 0] <- "<0.001"
#最后输出文件
write.xlsx(Univar,"Unicoxtable.xlsx",
col.names = T,row.names = F,
showNA = F)
接下来是多因素COX回归
#----------------------多因素COX回归---------------------#
# 筛选pvalue小于0.05的变量用于后续多因素分析
Univar$Characteristics[Univar$P.Value<0.05]
fml <- as.formula(paste0("Lusurv~",paste0(Univar$Characteristics[Univar$P.Value<0.05],collapse = "+")))
MultiCox <- coxph(fml,data = lung)
MultiSum <- summary(MultiCox)
MultiName <- as.character(Univar$Characteristics[Univar$P.Value<0.05])
MHR <- round(MultiSum$coefficients[,2],2)
MPValue <- round(MultiSum$coefficients[,5],3)
MCILow <- round(MultiSum$conf.int[,3],2)
MCIUp <- round(MultiSum$conf.int[,4],2)
MCI <- paste0(MCILow,"-",MCIUp)
MulCox <- data.frame("Characteristics" = MultiName,
"Hazard Ratio" = MHR,
"CI95" = MCI,
"P Value" = MPValue)
MulCox$P.Value[MulCox$P.Value == 0] <- "<0.001"
#输出结果
write.xlsx(MulCox,"MulCoxtable.xlsx",
col.names = T,row.names = F,
showNA = F)
#---------------------将单因素和多因素的数据框合并---------------------#
# 两个数据库行不一样,按照变量名匹配
outPut <- merge.data.frame(Univar,MulCox,
by="Characteristics",
all = T,sort = T)
#输出结果
write.xlsx(outPut,"COXtable.xlsx",
col.names = T,row.names = F,
showNA = F)
至此,就完成了生存分析、单/多因素COX回归以及基线表的制作
附诗一首
晚风吻尽荷花叶
任我醉倒在池边
夜梦美人来相伴
乍醒犹思残酒甜