【解放双手|热爱生活】从TCGA数据循环下载开始!

前段时间看到果子老师朋友圈说：如果一生一次，那就弯腰地铁请人做 | 如果一生一世，那就砸锅卖铁自己学。我梦想着要与Bioinformatics一生一世的，emo了一个多月，既然捡起来就要操练起来。解放双手，热爱生活，从循环开始。

1. 准备工作

rm(list=ls())
setwd("D:/R/BioData/TCGA")
library(TCGAbiolinks)
library(org.Hs.eg.db)
library(AnnotationDbi) 
library(SummarizedExperiment)
library(dplyr)
packageVersion("TCGAbiolinks")

projects <- getGDCprojects()
projects <- projects %>% 
  as.data.frame() %>% 
  select(project_id,tumor) %>% 
  filter(grepl(pattern="TCGA",project_id))

# sapply(projects$project_id, dir.create)
# ps:本来打算根据project id先全部建立个文件夹的，后面放到循环中了。

2. TCGA_BatchDownload定义

TCGA_BatchDownload <- function(i){
  ## 1.显示运行信息
  print(paste0("Downloading number ", i, ",project name: ", projects$project_id[i]))
  ## 2.根据project_id分别建立文件夹
  dir.create(projects$project_id[i])
  ## 3.查询信息
  query_exp = GDCquery(project = projects$project_id[i], 
                       data.category = "Transcriptome Profiling",
                       data.type = "Gene Expression Quantification",
                       workflow.type = "STAR - Counts")
  ## 4.正式下载
  GDCdownload(query_exp)
  ## 5.多个数据合并
  pre_exp = GDCprepare(query = query_exp)
  ## 6.提取表达量数据
  library(SummarizedExperiment)
  countsdata <- SummarizedExperiment::assay(pre_exp,1)
  fpkmdata <- SummarizedExperiment::assay(pre_exp,5)
  tpmdata <- SummarizedExperiment::assay(pre_exp,4)
  gene_id <- data.frame(id=rowData(pre_exp)@listData[["gene_id"]], gene_name=rowData(pre_exp)@listData[["gene_name"]],gene_type=rowData(pre_exp)@listData[["gene_type"]])
  counts <- cbind(gene_id, countsdata)
  fpkm <- cbind(gene_id, fpkmdata)
  tpm <- cbind(gene_id, tpmdata)
  
  ## 7.下载临床信息
  clinical <- GDCquery_clinic(project = projects$project_id[i], type = "clinical")

  ## 8.保存数据
  folder <- paste0("./",projects$project_id[i],"/") 
  save(counts, file = paste0(folder, projects$project_id[i], "_counts.Rdata"))
  save(fpkm, file = paste0(folder, projects$project_id[i], "_fpkm.Rdata"))
  save(tpm, file = paste0(folder, projects$project_id[i], "_tpm.Rdata"))
  write.csv(clinical, file = paste0(folder, projects$project_id[i], "_clinical.csv"),row.names = F)
}

3. 开始下载

for (i in 1:33) {
  TCGA_BatchDownload(i)
}

4. 下载结果

图1: 下载到指定文件夹下
图2: 每个项目下为4个文件