TCGA数据下载方式小结

之前对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数据）。

最后加所有样本都加入购物车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中输入下载命令：

1. gdc-client download -m gdc_manifest_20171102_125645.txt

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

• 首先我们下载后的结果是一个个文件夹，如：cb16f668-0e9c-4a0c-ae9b-eaefa0118dde；文件夹里面是count数据的压缩包，如:9b529e72-7a69-4843-828c-f6ad6d1ee99b.htseq.counts.gz；而我们需要将这个看不懂的名称改为TCGA-AR-A0TZ-01A这种形式，方便后续处理。

这种对应关系可以从Cart界面下载metadata文件中获取，然后提取出三者的对应关系。最开始是用perl的，但是感觉perl会比较繁琐，所以改用R来处理。选择的是 rjson这个包，也试了下 jsonlite包，基本功能是差不多的。

1. library(rjson)

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

3. data2 <- sapply(data, function(x){

4.        file_name <- x$file_name

5.        file_id <- x$file_id

6.        submitter_id <- x$case[[1]]$samples[[1]]$submitter_id

7.        tmp <- c(file_name, file_id, submitter_id)

8.      })

9. data2 <- t(data2)

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

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

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

• 接下来要对下载的文件进行批量解压缩，并输出到 file文件夹中，该shell脚本输入参数为上述的sampleid.txt文件。这个代码有个问题就是如果遇到样品ID相同的两个样品的，后者文件会覆盖前者文件。当然只要将每个样本的count文件输出在各自文件夹中就不会有这个错误了（但反正后续还是打算去掉重复样本，所以先这么处理着）。

1. #!/bin/bash

2. cat $1 |while read line

3. do

4.   arr=($line)

5.   filename=${arr[0]}

6.   folder=${arr[1]}

7.   submitterid=${arr[2]}

8.   gunzip -c ./${folder}/${filename} > ./file/${submitterid}.count

9. done

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

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

1. #!/bin/perl -w

2. use strict;

3. my $path = shift @ARGV;

4. opendir DIR, $path or die;

5. my @dir = readdir DIR;

6. my $header;

7. my @sample;

8. my %hash;

9. my %rm_rep;

10. foreach my $file (@dir) {

11.       if ($file =~ /^(\w+.*-\d+)[A-Z]\.count/) {

12.           my $sample_id = $1;

13.           next unless ($sample_id =~ /01$/ || $sample_id =~ /11$/);

14.           if ($rm_rep{$sample_id}) {

15.               next;

16.           }else{

17.               $rm_rep{$sample_id} = 1;

18.           }

19.           push @sample, $sample_id;

20.           $header .= "\t$sample_id";

21.           open my $fh, $path.$file or die;

22.           while (<$fh>) {

23.               chomp;

24.               next if ($_ =~ /^_/);

25.               my @array = split /\t/, $_;

26.               $hash{$array[0]} -> {$sample_id} = $array[1];

27.           }

28.           close $fh;

29.       }

30. }

31. print "$header\n";

32. map{

33.    my $gene = $_;

34.       print "$gene";

35.       foreach my $file (@sample) {

36.           print "\t".$hash{$gene} -> {$file};

37.       }

38.       print "\n";

39. }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.orgBRCA.MergernaseqilluminahiseqrnasequnceduLevel3geneexpression_data.Level3.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 转化成entrez id/gene symbol，很是方便。

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

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

生信人开发的小工具

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

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

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