KeepNotes blog

Stay hungry, Stay Foolish.

0%

TCGA数据下载方式小总结

Posted on Edited on In ,
Symbols count in article: 4.6k Reading time ≈ 4 mins.

之前对TCGA做了简单的了解，粗略了解了什么是TCGA，TCGA是做什么的等情况，接下来肯定是要学会如何下载TCGA数据，毕竟只有下载了数据才能接下去继续实战学习

官网常规下载

TCGA自2016年改版后，下载方式变得大为不同，数据都整合在GDC（Genomic Data Commons）的DATA PORTAL中，所以我们要先进网址https://portal.gdc.cancer.gov/。 我以乳腺癌（BRCA）为例，进入TCGA-BRCA界面，这次选择下载的是RNA-Seq数据，点击RNA-Seq的1092个case，然后在点击Files自行设定筛选条件（我选择下载count数据，记得还要选open的数据），最后加所有样本都加入购物车Cart（总共会有1222个文件）。为什么会比case数多呢，说明一些case不止对应1个RNA-Seq样本的count文件，有些是重复。

接下来就稍微有点不同了，如果只是比较少量的样本，直接点击Download即可下载，每个样本都是一个压缩文件。但是，TCGA数据库在数据下载有规定：让Cart文件夹大于50M时（这个依据网络情况，和下载用户数目），只能通过Data Transfer Tool工具进行下载。所以我这次要使用Data Transfer Tool工具来下载数据。

首先是要安装Data Transfer Tool，总共工具很简单，我是选择windows系统来安装的，下载后解压缩即可使用。至于怎么用，我刚开始看到是.exe文件，就以为是软件界面的，后来才发现原来是通过命令行下载的。先从刚才是Cart界面点击Download选择下载Manifest，然后在windows的cmd中输入下载命令

gdc-client download -m gdc_manifest_20171102_125645.txt

然后就是慢慢下载了，如果中间有报错断了，一般就是网络不好的原因，重新下载或者换个时间段就行了，由于下载的结果是一个个样本的压缩文件，1222个样本就是1222个压缩包，所以我要解压缩这些文件并合并成1个矩阵文件，主要可以分为3步：

  1. 将样本名称改为TCGA-C8-A12Q-01这种形式，这种ID可以从Cart界面下载metadata文件。然后在metadata文件中提取下载文件名和样本名的对应关系，最开始是用perl的，但是感觉perl会比较繁琐，所以改用R来处理。选择的是rjson这个包，也试了下jsonlite包，基本功能是差不多的

     library(rjson)

 data <- fromJSON(file = "metadata.cart.2017-10-30T08-39-11.222982.json")

 data2 <- sapply(data, function(x){
   file_name <- x$file_name
   file_id <- x$file_id
   submitter_id <- x$case[[1]]$samples[[1]]$submitter_id
   tmp <- c(file_name, file_id, submitter_id)
 })

 data2 <- t(data2)
 colnames(data2) <- c("filename", "fileid", "submitterid")

 write.table(data2, file = "sample2id.txt", sep = "\t", col.names = F, quote = F, row.names = F)

    输出的sample2id.txt文件中总共有3列，第一列为我们需要解压的文件名，第二列为我们下载下来的每个压缩包所在的文件夹的名称，第3列则是submitter id，就是我们需要的样品名

  2. 接下来要对下载的文件进行批量解压缩，并输出到file文件夹中，该shell脚本输入参数为上述的sampleid.txt文件。这个代码有个问题就是如果遇到样品名相同的两个样品的，后者文件会覆盖前者文件。当然只要将每个样本的count文件输出在各自文件夹中就不会有这个错误了（考虑到这样做的话，后面还是要去掉重复样本，所以暂时不这样做了）

     #!/bin/bash

 cat $1 |while read line
 do
   arr=($line)
   filename=${arr[0]}
   folder=${arr[1]}
   submitterid=${arr[2]}
   gunzip -c ./${folder}/${filename} > ./file/${submitterid}.count
 done

    运行上述bash extract_files.sh ./sample2id.txt脚本后，1222个TCGA乳腺癌的样本变成了1217个样本（说明有5个样本的名称是重复的，所以这些文件被覆盖了），这是在file文件夹中有1217个样本的count文件

  3. 最后将每个样本的count数进行合并，这次我用perl来整合成一个count矩阵文件，并且只保留01A和11A结尾的样本（01A表示肿瘤细胞，11表示癌旁样本）。如果以B结尾的样品名，则表明是重复，其实还有C结尾的，这次都先选择去除

     #!/bin/perl -w
 use strict;

 my $path = shift @ARGV;

 opendir DIR, $path or die;
 my @dir = readdir DIR;

 my $header;
 my @sample;
 my %hash;
 my %rm_rep;

 foreach my $file (@dir) {
     if ($file =~ /^(\w+.*-\d+)[A-Z]\.count/) {
         my $sample_id = $1;
         next unless ($sample_id =~ /01$/ || $sample_id =~ /11$/);
         if ($rm_rep{$sample_id}) {
             next;
         }else{
             $rm_rep{$sample_id} = 1;
         }
         push @sample, $sample_id;
         $header .= "\t$sample_id";

         open my $fh, $path.$file or die;
         while (<$fh>) {
             chomp;
             next if ($_ =~ /^_/);
             my @array = split /\t/, $_;
             $hash{$array[0]} -> {$sample_id} = $array[1];
         }
         close $fh;
     }
 }

 print "$header\n";
 map{
     my $gene = $_;
     print "$gene";
     foreach my $file (@sample) {
         print "\t".$hash{$gene} -> {$file};
     }
     print "\n";
 }keys %hash;

    运行perl merged_count.pl ./file/，代码看起来有点繁琐，但是还是可以用的~~最后输出一个整合后的count矩阵文件，后续分析就是根据这个count数据来了

Firehose下载

第一种方式个人觉得比较可控，下载的数据的准确性也是有保证的。但是如果觉得麻烦的话，可以使用Broad Institute's GDAC Firehose，该网站不仅提供TCGA数据的下载还有其分析后的数据的展示。这里我只用过其下载的功能，所以其他的自行去网站http://gdac.broadinstitute.org/逛逛吧

进入网站后，点击上方的Dashboards下的Standard Data，然后再点击网页上的乳腺癌（BRCA）对应的open按钮，结果会出现一系列文件可供下载，为什么更好的跟第一种方法进行比较，我选择了下载RNA-Seq的count数据对应的文件，可以选择下载gdac.broadinstitute.org_BRCA.Merge_rnaseq__illuminahiseq_rnaseq__unc_edu__Level_3__gene_expression__data.Level_3.2016012800.0.0.tar.gz，其实还有个rnaseqv2对应文件，刚开始我也是下载错了，后来查了下才发现前者和后者是有区别的，官方给的解释是：rnaseqv2相比rnaseq采用不同的算法（我看了下似乎是用RSEM进行定量的？因为rnaseq是采用htseq来定量的），但是呢rnaseqv2相比rnaseq样本数多一点。。。

RNASeq Version 2 is similar to RNASeq in that it uses sequencing data to determine gene expression levels. RNASeq Version 2 uses a different set of algorithms to determine the expression levels are the results are presented in a slightly different set of files.

下载后将文件夹解压，里面有不少文件，我先只拿count文件来看看，这里就显示出Firehose下载的好处了，不用像第一种方法那样将一个个样本合并，其下载结果就是合并好的数据，而且已经将Ensembl id 转化成gene symbol，很是方便

但是！我发现一个问题，就是Firehose下载下来的文件的样本数跟第一种方法下载后并不一致，粗略看了下，大概少了几百个样本。。就算是rnaseqv2也少了几个样本。然后我将两个下载方法获得的样本名都进行去重，然后比较发现有几个样本虽然在第二个方法中出现，但是在TCGA官网上并没有RNA-Seq的数据了！！！！这时我才想起来Firehose数据并不是跟TCGA官方数据更新同步的，现在下载的Firehose的数据是2016年初的。。。

所以Firehose下载操作确实方便，但是数据并不是时时更新开放下载，如果不在意样本数目的话，这个方式确实是一个比较方便的下载TCGA数据的方法

生信人开发的小工具

我也是听别人介绍了才知道生信人公众号做了一个下载TCGA数据的小工具，初步尝试了下，感觉不错，操作蛮简单的，并且不仅提供了数据下载功能，还有合并样本、临床信息还有ID转化的功能

具体操作还是看生信人的宣传文章把https://www.shengxin.ren/article/27

这个优点我觉得在于下载数据的实时性有保证，并且它也是从官网下载的，最主要它支持断点下载。。。就算你网再不稳定。。睡一觉总能下好了，比第一种下载方式方便（因为第一种网不好容易断，导致一些文件无法被正常下载），而且最后整合的样本数目（去重后）跟第一种方法一致，这让我心里也有了底。。。。

还有一个，该工具在使用ID转化功能后，能将lncRNA区分开。。。这个就满足好多人的需求了。。但是我也准备自己来写个代码做这件事。。比较下结果。。

本文已在《生信技能树》公众号上投稿，并授权转载

本文出自于http://www.bioinfo-scrounger.com转载请注明出处

Related Posts