TCGAbiolinks数据下载TCGA数据
- 下载TCGA数据的方法有很多,但比较好用的包我认为就是TCGAbiolinks,TCGAbiolinks是一个可用于检索,下载,并准备TCGA数据用于下游分析的R包,该包发表在 著名的 核酸研究 杂志上。
- TCGAbiolinks的优点在于具备一体化的下载整合,无需再使用复杂的方法对下载的单个数据重新进行整合,新手及临床医生尤其适合,我们的目的就是分析数据,没有必要去做些非必须的事。
- UCSC xena 本来是一个很好的工具,但有一个问题是它存储的数据是log2(x 1) 转换的 RSEM nomalized data,关于这个数据有文章去使用,但我不确定能否直接用limma去做后续分析。笔者曾试图搞清楚这个数据究竟如何使用,发现众说纷纭,似乎UCSCxena本身也没有给出一个很清晰的答案,因此就表达谱数据而言,UCSCxena并不是最好的选择,当然也有说把log2转换的数据再转换回来的,总觉得并不合适。
load package
代码语言:javascript复制1if (!requireNamespace("BiocManager", quietly = TRUE))
2 install.packages("BiocManager")
3if (!requireNamespace("BiocManager", quietly = TRUE))
4 install.packages("TCGAbiolinks")
5library(TCGAbiolinks)
6library(DT)
7library(dplyr)
8library(SummarizedExperiment)
数据来源-根据TCGAbiolinks的官方说明
- 不同的数据来源 Legacy 与 harmonized
- GDC Legacy Archive:以前在CGHUB和TCGA数据门户中存储的数据的原始数据,由TCGA数据协调中心(DCC)托管,在该门户中用GRCH37(HG19)和GRCH36(HG18)作为参考基因组
- GDC harmonized database:可用数据与grch38(hg38)使用gdc生物信息学流程进行协调,该流程提供了生物标本和临床数据标准化的方法,简单讲就是对数据进行了一定标准化处理。harmonized数据库包括转录谱数据,甲基化数据,miRNA数据,但缺少芯片数据
数据检索
代码语言:javascript复制 1query.met <- GDCquery(project = "TCGA-COAD",
2 data.category = "DNA Methylation",
3 legacy = FALSE,
4 platform = c("Illumina Human Methylation 450"))
5query.exp <- GDCquery(project = "TCGA-COAD",
6 data.category = "Transcriptome Profiling",
7 data.type = "Gene Expression Quantification",
8 workflow.type = "HTSeq - FPKM-UQ")
9query.exp[1:5,1:5]
10query.met[1:5,1:5]
11datatable(getResults(query.exp, cols = c("data_type","cases")),
12 filter = 'top',
13 options = list(scrollX = TRUE, keys = TRUE, pageLength = 5),
14 rownames = FALSE)
mRNA Expression
代码语言:javascript复制 1query <- GDCquery(project = "TCGA-BRCA",# Cancer type
2 data.category = "Transcriptome Profiling",
3 data.type = "Gene Expression Quantification",
4 workflow.type = "HTSeq - Counts"# raw count
5 )
6GDCdownload(query,
7 directory = "./project",#
8 method = "api",
9 files.per.chunk = 100)#
10data <- GDCprepare(query)
11count_data=assay(data)
12count_data[1:5,1:5]
13dim(count_data)##56537 1222
14
15## clinical information
16colData(data)[1:5,1:5]
17
18## save data
19#save(count_data,file = "BRCA_count.Rdata")
clinical data
代码语言:javascript复制 1clinical <- GDCquery_clinic(project = "TCGA-BRCA", type = "clinical")
2clinical[1:5,1:5]
3dim(clinical)## 1097 68
4
5## save
6#save(clinical,file="BRCA_clinical.Rdata")
7#write.csv(clinical,file="TCGAbiolinks-BRCA-clinical.csv")
8
9## clinical-2
10clinical_2<-colData(data)
11#write.csv(clinical,file="TCGAbiolinks-BRCA-clinical.csv")
12
13## 获取所有TCGA的临床信息
14library(data.table)
15library(dplyr)
16library(regexPipes)
17clinical <- TCGAbiolinks:::getGDCprojects()$project_id %>%
18 regexPipes::grep("TCGA",value=T) %>% ## TCGA
19 sort %>%
20 plyr::alply(1,GDCquery_clinic, .progress = "text") %>%
21 rbindlist
22dim(clinical)
23clinical[1:5,1:5]
24#readr::write_csv(clinical,path = paste0("BRCA_clin_indexed.csv"))
miRNA data
代码语言:javascript复制 1query <- GDCquery(project = "TCGA-BRCA",
2 data.category = "Transcriptome Profiling",
3 data.type = "miRNA Expression Quantification",
4 workflow.type = "BCGSC miRNA Profiling")
5## 检束结果
6results<-getResults(query)
7dim(results)
8results[1:5,1:5]
9colnames(results)
10
11## download data
12GDCdownload(query,
13 method = "api",
14 files.per.chunk = 20)# 减少下载失败风险
15mir_exp<- GDCprepare(query = query,
16 summarizedExperiment=F)# set F
17
18## 数据结构
19dim(mir_exp)
20mir_exp[1:5,1:5]
21
22## save data
23#save(mir_exp,file="BRCA_miRNA_raw.Rdata")
总结
除了下载功能以外,TCGAbiolinks还包括了一些数据分析挖掘功能:
- 差异分析
- 富集分析
- 整合其它组学数据的关联分析
- 一定的可视化功能:PCA,heatmap等等
笔者认为,这些功能的好处在于能够比较方便的用简单几行代码做一些固定的分析,恰恰这也是它的缺点,这样的话包的灵活性就大打折扣了,不能够很好的整合用户的其它分析流程。考虑到这一点,我可能不会过多的使用它的其它分析功能