作者:向屿
来源:生信控
森林图在Meta分析中应用广泛,在生存分析中的应用其实并无太多高深之处,不过是将多因素COX分析的统计结果进行一个可视化的展示。
以survival包中的肠癌数据集作为测试数据:
p_load(survival)
data(colon)
查看数据:
options(width=)#增加显示宽度
head(colon)
使用?colon可以查看该数据集的具体释义。本例中,我们将用到3列分类变量:
1、sex(性别,1为男性,0为女性)
2、rx(治疗方式,Obs(ervation),Lev(amisole),Lev+5-FU)
3、adhere(附着于邻近器官,1和0分别代表是和否)
使用coxph()函数构建预后模型:
#构建COX比例风险回归模型
cox_model=coxph(Surv(time,status)~sex+rx+adhere,data=colon)
summary(cox_model)
进而,使用survminer包中的ggforest()函数绘制森林图:
Thefunctionggforest()fromthesurvminerpackagecreatesaforestplotforaCoxregressionmodelfit.Hazardratioestimatesalongwithconfidenceintervalsandp-valuesareplotterforeachvariable.
p_load(survminer)
ggforest(cox_model)
显示如下:
简单方便吧,一个命令就能完成。不过,上图中值得注意的是:
1、三者都是分类变量,但rx显示了具体的类别(显示三行),而sex和adhere则没有。原因是rx是因子型,要将sex和adhere转为因子型数据:
colon=within(colon,{
sex=factor(sex,levels=c(0,1),labels=c("female","male"))
adhere=factor(adhere,levels=c(0,1),labels=c("NO","YES"))
})
2、本例中,rx的Obs作为了参考(reference,HR=1),并计算其他两类相对于Obs的HR值,是因为因子型的第一个类别就是Obs:
Levels:ObsLevLev+5FU
所以,如果想使用Lev作为参考可以更换其顺序:
colonrx=factor(colonrx,levels=c("Lev","Obs","Lev+5FU"))
3、ggforest()函数支持一些参数调整。结合上述两条,执行代码如下:
#重新执行COX
cox_model=coxph(Surv(time,status)~sex+rx+adhere,data=colon)
#带参数的
ggforest(cox_model,main="Hazardratio",cpositions=c(0.02,0.12,0.35),fontsize=0.8,refLabel="Ref.",noDigits=3)
注意,一定要在变量转换数据类型后重新执行COX,否则报错:
错误:variables‘sex’,‘adhere’werespecifiedwithdifferenttypesfromthefit此外:Warningmessage:In.get_data(model,data=data):Thedataargumentisnotprovided.Datawillbeextractedfrommodelfit.
正确执行后,显示如下:
数据挖掘·往期精彩
数据虐我千百遍、T检验、相关性分析、全自动GEO数据下载!TCGA临床信息大整合、生存分析中预测变量的类型...
往期系列推送:
想学R?从看R界传奇老司机的免费讲座开始吧
R语言概述
R入门教程之变量
专题
R语言基础篇—列表,因子和数据框
专题
R语言基础篇—矩阵
生存分析——R语言
“R”语言基础列之(一)——学渣也能搞定高大上的“R”语言啦
R语言基础系列之(二)————给我一个爱上你的理由
R软件基本操作,第三课
第四课~!安装
第五课Packages
第六课R的内置数据集
第七课送ncRNA数据库偶!
第八课R语言数据结构
第九课向量的索引
第十课向量的计算
第十一课矩阵
第十二课数组
第十三课列表
第十四课数据框
第十五课因子
第十六课缺失数据
第十七课字符串
第十八课时间日期处理
发送生信到本