第一步:进入TCGA官网下载数据
-
进入官网数据下载界面后,点击Repository
-
在Files和Cases两个界面中选择自己需要的数据
-
FPKM数据格式在Workflow选项中选择
-
选择完数据后,点击Add All Files to Cart
-
点击右上角的Cart,Cart旁边的数字表示我们选择了多少个样本的数据
-
点击Download按钮然后出来两个选项,Manifest表示使用GDC Data Transfer Tool读取Manifest文件下载数据,Cart表示直接通过浏览器链接下载数据。同时还需要点击Metadata按钮,下载Metadata文件
注:GDC Data Transfer Tool下载数据方法见 https://www.jianshu.com/p/f4e92d226e6d
第二步:使用R合成表达矩阵
- 下载下来数据是很多个文件夹,每个文件夹是一个样本的数据,因此文件夹个数应该等于Cart中的个数,如果不等,代表我们下载数据时有丢失。我们将所有的数据文件夹拷贝到一个文件夹中,命名为rawdata
- 在R中处理,代码如下:
rm(list = ls())
options(stringsAsFactors = F)
#我们自己设置好工作路径,然后将rawdata文件夹拷贝到工作路径下
dir.create("all_data")
for (dirname in dir("rawdata/")){
file <- list.files(paste0(getwd(),"/rawdata/",dirname),pattern = "*.FPKM")
file.copy(paste0(getwd(),"/rawdata/",dirname,"/",file),"all_data")
}
dir.create("unpacked_FPKM")
#所有样本的单个文件都拷贝在了all_data文件夹中,但是这些文件都是压缩格式的,然后使用解压缩软件将所有压缩文件统一解压缩到unpacked_FPKM文件夹中
metadata <- jsonlite::fromJSON("metadata.cart.2020-04-24.json")
require(dplyr)
metadata_id <- metadata %>%
dplyr::select(c(file_name,associated_entities))
naid_df <- data.frame()
for (i in 1:nrow(metadata)){
naid_df[i,1] <- substr(metadata_id$file_name[i],1,nchar(metadata_id$file_name[i])-3)
naid_df[i,2] <- metadata_id$associated_entities[i][[1]]$entity_submitter_id
}
colnames(naid_df) <- c("filename","TCGA_id")
#naid_df储存了文件名和TCGA_id的对应关系
files <- dir("unpacked_FPKM")
myfread <- function(files){
data.table::fread(paste0("unpacked_FPKM/",files))[,2]
}
f <- lapply(files,myfread)
f <- do.call(cbind,f)
rownames(naid_df) <- naid_df[,1]
naid_df <- naid_df[files,]
colnames(f) <- naid_df$TCGA_id
gene_id <- data.table::fread(paste0("unpacked_FPKM/",files[1]))$V1
expr_df <- cbind(gene_id=gene_id,f)
save(expr_df,naid_df,file = "FPKM_ENSG_exprdf.Rdata")
第三步:将ensembl数据库的ENSG编号转换成gene symbol
- 在ensembl数据库中下载数据,网址:http://asia.ensembl.org/index.html
-
点击Downloads→databases→Human选项中的GTF
-
下载图示文件
- R中处理:
rm(list=ls())
options(stringsAsFactors = F)
gtf<-rtracklayer::import('Homo_sapiens.GRCh38.100.chr.gtf')
gtf_df <- as.data.frame(gtf)
save(gtf_df,file = "gtf_df.Rdata")
load("FPKM_ENSG_exprdf.Rdata")
metadata <- naid_df[,-1]
metadata<-data.frame(TCGA_id=metadata)
require(dplyr)
require(tidyr)
expr_df_nopoint <- expr_df %>%
tidyr::separate(gene_id,into = c("gene_id"),sep="\\.")
#提取蛋白编码基因
mRNA_exprSet <- gtf_df %>%
dplyr::filter(type=="gene",gene_biotype=="protein_coding") %>%
dplyr::select(c(gene_name,gene_id,gene_biotype)) %>%
dplyr::inner_join(expr_df_nopoint,by ="gene_id") %>%
tidyr::unite(gene_id,gene_name,gene_id,gene_biotype,sep = " | ")
save(mRNA_exprSet,file = "mRNA_exprSet.Rdata")
#提前LncRNA的基因
ncRNA <- c("sense_overlapping","lincRNA","3prime_overlapping_ncRNA",
"processed_transcript","sense_intronic",
"bidirectional_promoter_lncRNA","non_coding",
"antisense_RNA")
LncRNA_exprSet <- gtf_df %>%
dplyr::filter(type=="transcript",transcript_biotype %in% ncRNA) %>%
dplyr::select(c(gene_name,gene_id,transcript_biotype)) %>%
dplyr::distinct() %>%
dplyr::inner_join(expr_df_nopoint,by ="gene_id") %>%
tidyr::unite(gene_id,gene_name,gene_id,transcript_biotype,sep = " | ")
save(LncRNA_exprSet,file = "LncRNA_exprSet.Rdata")