数据库numeric_TCGA数据库：生存分析

我们前面介绍了TCGA数据库中各种数据的下载处理：

TCGA数据库：RNA-Seq数据的下载与处理

TCGA数据库：SNP数据的下载整理及其可视化

TCGA数据库：miRNA数据下载与整理

TCGA数据库：临床数据下载与整理

TCGA数据库：ATAC-Seq数据的下载整理及其可视化

一文解决TCGA数据库临床数据下载与整理

也介绍了下游的差异分析：

TCGA数据库：GDCRNATools包下载数据、处理数据以及差异分析

一文就会TCGA数据库基因表达差异分析

也介绍了：基因表达谱热图绘制

本文介绍生存分析，其实，在R中，生存分析很简单，大家在网上能找到无数的文章。利用survival包就可以。就是按照下列公式就可以完成简单的生存分析。

fit ~ 分组, data=数据框)

我们这里就结合基因的表达量，来进行分析。

首先加载我们的数据。

options(stringsAsFactors = F)#加载表达数据load("F:/TCGA/HTSeq-FPKM/Rdata/data/TCGA-COAD-Exp.Rdata")#加载临床数据load("K:/TCGA/clinicalData/tidyAllCancerData/TCGA-COAD -Clindata.Rdata")

表达数据和临床数据，我之前已经上传到网盘：TCGA数据库33个Project的RNA-Seq转录组数据为你整理打包好了

处理方式可以参考文章：

【1】TCGA数据库：临床数据下载与整理；

【2】TCGA数据库：RNA-Seq数据的下载与处理

之前处理后的数据进行简单的处理，其实就是去掉正常组织的样本，再把列名变成3个字段。因为原来表达矩阵中病人的barcode长，"TCGA-AA-3662-11A-01R-1723-07"，而临床数据中的只有前3段。具体关于barcode这个在前面有介绍，也可以参考文章TCGAbiolinks包介绍，里面有详细介绍。

##提取肿瘤样本的表达矩阵exp "allGenesExp"]][[tumExp exp[,transomeData[[barcode "(.*?)-(.*?)-(.*?)-.*",colnames(tumExp) ##

同样我们也要处理一下临床数据，我们之前处理的临床数据是这样的：

我们这里也需要简单处理一下。

#提取临床数据library(dplyr)clin clin "--")clin "sample",clin[1,3]clin$follow as.numeric(clin$follow)/

其实，就是只要生存时间和生存状态这2个数据，再删除缺失值。follow除了365，单位就是年啦。

然后我们将表达矩阵与临床数据融合，因为不是每个病人的数据都是一一对应的，简单说，就是病人有表达数据，但他的临床数据就不全，我们也删除了缺失值的病人的临床数据，所以我们只需要具有临床数据又有表达数据的病人的数据。也就是取一个交集。

##融合数据interbarcode gene "LLGL2"geneExp geneExp % as.data.frame() %>% data.matrix() %>% as.data.frame() geneExp$sample mergdata by = 

这里的gene我只写了一个，所以融合后是下面这样的数据。你也可以把所有基因都融合进去，这里是案例，就演示了一个。

接下来就是添加一个分组信息。

mergdata$Group  median(mergdata[,gene]),"High","Low")

我们以表达值的中位数为分界线，高于中位值为高表达，低于或等于中位值为低表达。也可以用均值。

得到上面这样的数据后，我们就可以按照刚刚的公式进行生存分析了：

######################### 生存分析library(survival)library(survminer)fit 

绘图的话就用ggsurvplot函数。

ggsurvplot(fit,data = mergdata)

好像不是很美观，我们可以调整一下参数，比如y轴下部分很空，我们可以调整一下y轴坐标。

这样看着就好很多啦。我们在进行其他参数的调整。

ggsurvplot(fit, main = "Survival curve",           font.main = c(16, "bold", "darkblue"),           font.tickslab = c(12, "plain", "darkgreen"),           legend.title = gene,           legend.labs = c("High", "Low"),           ylim = c(0.75,1),           # Add p-value and tervals           pval = TRUE,           palette = c("#FF0000", "#2200FF"),           ggtheme = theme_bw() )

如果我们要一次批量分析很多基因的高低表达与生存的关系，写一个循环，批量绘图了，也可以参考文章GraphPad Prism绘制生存曲线和R语言批量绘制生存曲线。

---

尽管本文是介绍基因表达量的生存分析，但其他的也是一样，就看你怎么分组，比如我们前面介绍SNP的数据处理后，能否做某基因突变与野生型的生存分析呢？其实都是一样的道理，其他的也是一样。照葫芦画瓢而已，大家自己去试试。

---

一个R脚本解决某类功能基因(比如m6A甲基化)临床预后模型分析流程

肿瘤免疫细胞浸润与临床相关性分析

---

相关专辑

TCGA | 文献阅读 | R语言 | 数据库 | 理论知识

---