单因素独立预后分析，随机森林图绘制教程
[](http://mp.weixin.qq.com/s?__biz=MzAwODY5NDU0MA==&mid=2455861689&idx=2&sn=bf497c5f96fbccc071751f5e9eb087a5&chksm=8cff27dabb88aecc14142b2712e8bbd0f4ed2d6660b4814cd2031212edf682bd8c61f2fb391a&scene=21#wechat_redirect)[](http://mp.weixin.qq.com/s?__biz=MzAwODY5NDU0MA==&mid=2455861689&idx=2&sn=bf497c5f96fbccc071751f5e9eb087a5&chksm=8cff27dabb88aecc14142b2712e8bbd0f4ed2d6660b4814cd2031212edf682bd8c61f2fb391a&scene=21#wechat_redirect)**」**
### 本期教程
> 本期教程，代码和数据可在社群中获得
### 往期教程部分内容
## 写在前面
**在医学研究中，为了研究临床病理特征与风险模型的预后情况，会将各种临床病理因素纳入风险模型，如年龄（Age），肿瘤分期（stage、T、 M、N）等，进行单因素Cox独立预后分析。并基于TCGA数据集，将目标基因表达状态和临床指标分别进行的单因素和多因素Cox回归分析，从而阐述目标基因与临床病理因素之间的关系。**
**今天的教程，我们先分享独立预后单因素Cox分析，并绘制随机森林图。**
### 分析代码
1. **加载R包**
```abap
##'@加载对应R包
library(survey)
library(survival)
library(survminer)
```
1. **导入数据**
```abap
##'@导入数据
data_df <- read.csv("ARFM3M_input_data.csv", header = T, check.names = F, row.names = 1)
head(data_df)
```
**3， 单因素独立预后分析**
```
##'@单因素独立预后
rownames(data_df) <- data_df$id
head(data_df)
##
#data_df <- data_df[,-1]
data_df$Age = ifelse(data_df$Age >= 65,2,1)
##
data_df$futime <- as.numeric(data_df$futime)
data_df$fustat <- as.numeric(data_df$fustat)
#
data_df <- data_df[!data_df$futime == 0,]
##
class(data_df)
outTab <- data.frame()
for(i in colnames(data_df[,3:ncol(data_df)])){
cox <- coxph(Surv(futime, fustat) ~ data_df[,i], data = data_df)
coxSummary = summary(cox)
coxP=coxSummary$coefficients[,"Pr(>|z|)"]
outTab=rbind(outTab,
cbind(id=i,
HR=coxSummary$conf.int[,"exp(coef)"],
HR.95L=coxSummary$conf.int[,"lower .95"],
HR.95H=coxSummary$conf.int[,"upper .95"],
pvalue=coxSummary$coefficients[,"Pr(>|z|)"])
)
}
```
```abap
##'@单因素独立预后
rownames(data_df) <- data_df$id
head(data_df)
##
#data_df <- data_df[,-1]
data_df$Age = ifelse(data_df$Age >= 65,2,1)
##
data_df$futime <- as.numeric(data_df$futime)
data_df$fustat <- as.numeric(data_df$fustat)
#
data_df <- data_df[!data_df$futime == 0,]
##
class(data_df)
outTab <- data.frame()
for(i in colnames(data_df[,3:ncol(data_df)])){
cox <- coxph(Surv(futime, fustat) ~ data_df[,i], data = data_df)
coxSummary = summary(cox)
coxP=coxSummary$coefficients[,"Pr(>|z|)"]
outTab=rbind(outTab,
cbind(id=i,
HR=coxSummary$conf.int[,"exp(coef)"],
HR.95L=coxSummary$conf.int[,"lower .95"],
HR.95H=coxSummary$conf.int[,"upper .95"],
pvalue=coxSummary$coefficients[,"Pr(>|z|)"])
)
}
```
1. **绘制随机森林图**
```abap
####'@绘制图形单因素森林图
rt <- read.csv("unoCox_result.csv", header = T, row.names = 1)
gene <- rownames(rt)
hr <- sprintf("%.3f",rt$"HR")
hrLow<- sprintf("%.3f",rt$"HR.95L")
hrHigh <- sprintf("%.3f",rt$"HR.95H")
Hazard.ratio <- paste0(hr,"(",hrLow,"-",hrHigh,")")
pVal <- ifelse(rt$pvalue<0.001, "<0.001", sprintf("%.3f", rt$pvalue))
#
#输出图形
pdf(file="ATAD3A.unicox.forest.pdf", width = 7,height = 6)
n <- nrow(rt)
nRow <- n+1
ylim <- c(1,nRow)
layout(matrix(c(1,2),nc=2),width=c(3,2.5))
#绘制森林图左边的临床信息
xlim = c(0,3)
par(mar=c(4,2.5,2,1))
plot(1,xlim=xlim,ylim=ylim,type="n",axes=F,xlab="",ylab="")
text.cex=0.8
text(0,n:1,gene,adj=0,cex=text.cex)
text(1.5-0.5*0.2,n:1,pVal,adj=1,cex=text.cex);text(1.5-0.5*0.2,n+1,'pvalue',cex=text.cex,font=2,adj=1)
text(3,n:1,Hazard.ratio,adj=1,cex=text.cex);text(3,n+1,'Hazard ratio',cex=text.cex,font=2,adj=1,)
#绘制森林图
par(mar=c(4,1,2,1),mgp=c(2,0.5,0))
xlim = c(0,max(as.numeric(hrLow),as.numeric(hrHigh)))
plot(1,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)
abline(v=1,col="black",lty=2,lwd=2)
boxcolor = ifelse(as.numeric(hr) > 1, 'red', 'green')
points(as.numeric(hr), n:1, pch = 15, col = boxcolor, cex=1.3)
axis(1)
dev.off()
```
> **> 在社群中可获得本期教程全部代码和示例数据**
**若我们的教程对你有所帮助，请**点赞+收藏+转发**，这是对我们最大的支持。**
### 往期部分文章
**「1. 最全WGCNA教程（替换数据即可出全部结果与图形）」**
* (https://mp.weixin.qq.com/s?__biz=MzAwODY5NDU0MA==&mid=2455850466&idx=1&sn=6036fc2a19594f7fc38d65271f21e9cc&scene=21#wechat_redirect)
* (https://mp.weixin.qq.com/s?__biz=MzAwODY5NDU0MA==&mid=2455850727&idx=1&sn=1156c2bb0b0c9baff02838f5ffce39bf&scene=21#wechat_redirect)
* (https://mp.weixin.qq.com/s?__biz=MzAwODY5NDU0MA==&mid=2455850742&idx=1&sn=3d3eeedb3b58e536a83dc38d15725cd4&scene=21#wechat_redirect)
* (https://mp.weixin.qq.com/s?__biz=MzAwODY5NDU0MA==&mid=2455854728&idx=5&sn=bc98befb8dd0f0090bdfe69f3ffdb008&scene=21#wechat_redirect)
* (https://mp.weixin.qq.com/s?__biz=MzAwODY5NDU0MA==&mid=2455856103&idx=1&sn=774f8a084d21f757266f35c501c7155d&scene=21#wechat_redirect)
---
**「2. 精美图形绘制教程」**
* [精美图形绘制教程](https://mp.weixin.qq.com/mp/appmsgalbum?__biz=MzAwODY5NDU0MA==&action=getalbum&album_id=2614156000866385923&scene=173&from_msgid=2455848496&from_itemidx=1&count=3&nolastread=1#wechat_redirect)
**「3. 转录组分析教程」**
* **「**[转录组上游分析教程[零基础]](https://mp.weixin.qq.com/mp/appmsgalbum?__biz=MzAwODY5NDU0MA==&action=getalbum&album_id=2870608342451224581&scene=126&uin=&key=&devicetype=Windows+10+x64&version=63090719&lang=zh_CN&ascene=0)**」**
* **「**[一个转录组上游分析流程 | Hisat2-Stringtie](https://mp.weixin.qq.com/s?__biz=MzAwODY5NDU0MA==&mid=2455857417&idx=1&sn=653f3bdb0af386c22a128732ec8a01a4&scene=21#wechat_redirect)**」**
**「4. 转录组下游分析」**
* [批量做差异分析及图形绘制 | 基于DESeq2差异分析](https://mp.weixin.qq.com/s?__biz=MzAwODY5NDU0MA==&mid=2455860684&idx=1&sn=7be489c453cca737ad1092c6e4499827&scene=21#wechat_redirect)
* (https://mp.weixin.qq.com/s?__biz=MzAwODY5NDU0MA==&mid=2455859512&idx=1&sn=bac01a018f8b58afc7e2b3484f476bf4&scene=21#wechat_redirect)
* [单基因GSEA富集分析](https://mp.weixin.qq.com/s?__biz=MzAwODY5NDU0MA==&mid=2455859147&idx=1&sn=b022d80868e4b8014f64d443c82d1668&scene=21#wechat_redirect)
* [全基因集GSEA富集分析](https://mp.weixin.qq.com/s?__biz=MzAwODY5NDU0MA==&mid=2455860201&idx=1&sn=dd65c5b967123a876a6f5a38d4723ca6&scene=21#wechat_redirect)
**「小杜的生信筆記」** **，主要发表或收录生物信息学教程，以及基于R分析和可视化（包括数据分析，图形绘制等）；分享感兴趣的文献和学习资料!!**
页:
[1]