for循环批量生存分析出图
尝试一下怎么批量生成图片,不过有点小bug还未解决。
library(survival)
library(survminer)
head(ALL) #我的数据
#OS OS_Status gene1 gene2 gene3 gene4 gene5 gene6 gene7 gene8 gene9 gene10
# 49.97 0 Gain Gain Diploid Diploid Gain Loss Loss Loss Gain Diploid
# 36.27 0 Diploid Diploid Diploid Diploid Diploid Diploid Diploid Diploid Diploid Diploid
# 19.68 1 Loss Diploid Gain Gain Loss Loss Loss Loss Diploid Diploid
# 31.50 0 Diploid Diploid Diploid Diploid Loss Diploid Diploid Loss Loss Diploid
# 53.45 1 Diploid Diploid Diploid Diploid Gain Diploid Diploid Diploid Loss Diploid
# 68.50 0 Diploid Diploid Diploid Diploid Diploid Loss Diploid Diploid Loss Gain
survfit估计生存率
for (i in 3:12){
surv=as.formula(paste('Surv(OS, OS_Status==1) ~ ',colnames(ALL)[i]))
KM <- survfit(surv,
data = ALL,
type = 'kaplan-meier',
conf.type = 'log')
#两两比较生存率是否有差异
p=pairwise_survdiff(surv, data = ALL,rho = 0)
#输出p<0.05的值
a=p$p.value
if (is.na(a[1])){
a[1]=10
} #这里是将NA值去掉
else if(is.na(a[2])){
a[2]=10
}
else if (a[1]<0.05 | a[2]<0.05){
print(paste(colnames(ALL)[i],p$p.value,sep = ' '))
}
}
用survminer包画生存曲线
#将每个基因提取出来作为一个List
D=list()
for (i in 3:12){
j=i-2
D[[j]]=data.frame('OS'=ALL$OS,
'OS_Status'=ALL$OS_Status,
'Type'=ALL[,i])
}
#绘图
for (i in 1:length(D)){
x=D[[i]]
KM_survfit <- survfit(Surv(OS, OS_Status==1) ~ Type
data = x,
type = 'kaplan-meier',
conf.type = 'log')
)
j=i+2 #这个2是因为我数据里从第三列开始才是基因
ggsave(filename = paste0(colnames(ALL)[j-1],'_OS','.pdf'),width = 5,height = 5) #就是这个bug,colnames我要设置成j-1而不是j最后输出结果才正确
p=ggsurvplot(KM_survfit,
data = x,
pval = T, #还有就是pairwise的p值好像不能自动生成
title=paste('Overall Survival -',colnames(ALL)[j],sep = ' ')) #这里colnames为j又没有问题
print(p)
dev.off()
}