在之前的推文（单因素、多因素COX回归分析在R中的实现与常见问题解答 - 简书 (jianshu.com)）中我们承诺过大家，后续会详细介绍一下用来可视化单因素、多因素COX回归分析结果的森林图。今天小碗就来把坑填上。
我们会给大家介绍两种方法。方法一，不仅能用来绘制单因素COX回归分析的森林图，而且可以用来绘制多因素COX回归分析森林图，代码稍微复杂，不过不用担心，要改动的部分小碗都给大家标注出来了，到时候只要与自己的数据对应起来就好。
###方法一：绘制单因素、多因素COX回归分析森林图
##加载数据
load("方法一_示例数据.rdata")#（改）
输入数据是包含变量名（gene）、风险比、风险比95%置信区间高值与低值、p值的数据框：
##提取绘制森林图所需的数据
variable <- uniforest$gene#提取变量名（改）
HR <- sprintf("%.3f",as.numeric(uniforest$HR)) #提取各变量的HR值（改）
HRLow <- sprintf("%.3f",as.numeric(uniforest$L95CI)) #提取各变量的HR95%置信区间低值（改）
HRHigh <- sprintf("%.3f",as.numeric(uniforest$H95CI)) #提取各变量的HR95%置信区间高值（改）
HR95 <- paste0(HR,"(",HRLow,"-",HRHigh,")") #将各变量的HR、95%置信区间高低值组合在一起
pValue <- ifelse(uniforest$pvalue<0.001, "<0.001", sprintf("%.3f", as.numeric(uniforest$pvalue)))#（改）
#如果p值小于0.001就写成<0.001，否则就显示原数值
## 开始绘制森林图
pdf(file="森林图.pdf", width = 12,height = 8)
n <- nrow(uniforest)#（改）
nRow <- n+1
ylim <- c(1,nRow) #设置y轴长度
#一张森林图可以分成左右两张图，左图是文字及数值部分，右图是图形部分，我们是分开画左右图
layout(matrix(c(1,2),nc=2),width=c(2.5,2)) #设置页面排版，即1排两图，左图1，右图2
## 先画左图，即文字及数值部分
xlim = c(0,2.5) #x轴长度
par(mar=c(4,2.5,2,1)) #页面边界
plot(0,xlim=xlim,ylim=ylim,type="n",axes=F,xlab="",ylab="")
text.cex=0.8 #设置字体大小
text(0,n:1,variable,adj=0,cex=text.cex) #第一列基因名
text(1.1,n:1,pValue,adj=1,cex=text.cex)
text(1.1,n+1,'pValue',cex=text.cex,font=2,adj=1) #第二列为p值
text(2.5,n:1,HR95,adj=1,cex=text.cex)
text(2.5,n+1,'HR(95% CI)',cex=text.cex,font=2,adj=1,) #第三列是HR及其95%置信区间
## 再画右图，即图形部分
par(mar=c(4,1,2,1),mgp=c(2,0.5,0))
xlim = c(0,max(as.numeric(HRLow),as.numeric(HRHigh))) #设置x轴长度
plot(0,xlim=xlim,ylim=ylim,type="n",axes=F,ylab="",xaxs="i",xlab="Hazard ratio")
arrows(as.numeric(HRLow),n:1,as.numeric(HRHigh),n:1,angle=90,code=3,length=0.05,col="darkblue",lwd=2.5) #表示HR95%置信区间的横线
abline(v=1,col="black",lty=2,lwd=2) #在横坐标1的地方画一条竖线
boxcolor = ifelse(as.numeric(HR) > 1, 'red', 'green') #设置颜色
points(as.numeric(HR), n:1, pch = 15, col = boxcolor, cex=1.3) #HR处显示为方块(通过参数pch修改)，颜色为上面设置的颜色。
axis(1) #显示横坐标
dev.off()
然后，你就会得到一张如下图所示的精美森林图啦！
那么问题来了，森林图该怎么看呢？左半部分，即数值与文本部分不用过多解释，三列内容分别是变量、每个变量对应的p值、HR及其95%置信区间。右半部分，即图形部分，比较复杂，我们分开看：深蓝色的线段代表HR的95%置信区间，绿色或红色的方块代表的就是每个变量的HR值，从图形最下方的x轴可以读出它们的大概数值；在x=1处有一条垂直的虚线，可以看到，有些变量的深蓝色线段与虚线相交而有些不相交，只有不相交才代表变量的HR的95%置信区间不包含1，p值小于0.05，变量对患者生存有影响且差异有统计学意义。
接下来，我们来介绍一下方法二。与方法一相比，方法二简直是傻瓜式操作，并且绘制出来的图高级感满满，但只能用来绘制多因素COX回归分析森林图。输入数据是拟合的多因素COX比例风险模型，以及各样本的生存信息，为了让大家看得更明白，这里是从多因素COX分析开始演示。
###方法二：多因素COX回归分析也可以直接用ggforest()画森林图
##安装包
install.packages("survminer")
install.packages("survival")
##加载包
library(survminer)
library(survival)
##加载数据
load("方法二_示例数据.rdata")
输入数据行名是样本名，前两列分别是样本的生存状态与生存时间，后面就是各个变量：
##多因素COX回归分析----
multicox <- coxph(Surv(time = os_time,event = os_status) ~ ., data = df)#（改）
## 绘制森林图
pdf(file="多因素COX森林图.pdf", width = 12,height = 8)
ggforest(multicox, #使用coxph函数拟合的多因素COX比例风险模型（改）
data = df, #包含各样本生存数据的数据框（改）
main = "Hazard ratio",
cpositions = c(0.02,0.12, 0.25),
fontsize = 0.8, #字体大小
refLabel = "reference",
noDigits = 2)
dev.off()
之后，就会得到下图所示的森林图：
左边数值与文本部分，第一列依旧是变量名，第二列是样本数，第三列是HR及其95%置信区间，右边图形各部分的含义与方法一绘制的森林图一样，最右边的一列数值就是各变量的p值了。
好了，森林图就给大家介绍到这里，如果这期推文对你有帮助的话就留下一个免费的赞吧！你的支持就是小碗码字的动力哦！
