TCGA 临床数据的提取

TCGA数据库除了有各种癌症的大样本数据外,还有很完善的临床数据,这点是其他数据库所缺少的,因此我们需要会从TCGA数据库中提取对应癌症的临床数据,然后利用这些数据来进行后续的分析(比如想知道某个基因的表达对病人的预后是否有影响)

下载方式还是跟之前的博文TCGA数据下载方式小总结一样,就是将下载选项中的RNA-Seq改为Clinical,以乳腺癌为例,总共有1097个样本的临床数据,一个样本一个文件,格式为xml

当然还有一个metadata文件,这是所有样本最近一次随访的临床数据,其格式为json

下面先以xml文件的临床数据为例,里面记录上了上百条临床信息数据,由于其是用xml来存储数据的,所以我们需要先了解我们需要提取的临床信息在xml里面叫啥名字,比如我需要发病年龄(age_at_initial_pathologic_diagnosis)、性别(gender)、TMN分期(stage_event-tnm_categories-pathologic_categories-pathologic_T, stage_event-tnm_categories-pathologic_categories-pathologic_M, stage_event-tnm_categories-pathologic_categories-pathologic_N), 生存状态(vital_status)以及生存时间(follow_ups-follow_up-days_to_death),当然还有一个是患者样本的id(bcr_patient_barcode)。临床基础有限,其他的信息暂时搞不太清楚

下面是用代码对上述信息的提取,对于xml格式,而且从文件中内容的排布来看,几乎是一条临床信息一行,因此我习惯用perl来处理这种格式

#!/usr/bin/perl -w
use strict;

my $in = shift @ARGV;

my ($patient_id,$age,$gender,$stage_pathologic_T,$stage_pathologic_M,$stage_pathologic_N,$vital_status,$days_to_death);

open my $fh, $in or die;
while (<$fh>) {
    chomp;
    $patient_id = $1 if ($_ =~ /<shared:bcr_patient_barcode.*?>(\S+?)<\/shared:bcr_patient_barcode>/);
    $age = $1 if ($_ =~ /<clin_shared:age_at_initial_pathologic_diagnosis.*?>(\S+?)<\/clin_shared:age_at_initial_pathologic_diagnosis>/);
    $gender = $1 if ($_ =~ /<shared:gender.*?>(\S+?)<\/shared:gender>/);
    $stage_pathologic_T = $1 if ($_ =~ /<shared_stage:pathologic_T.*?>(.*?)<\/shared_stage:pathologic_T>/);
    $stage_pathologic_M = $1 if ($_ =~ /<shared_stage:pathologic_M.*?>(.*?)<\/shared_stage:pathologic_M>/);
    $stage_pathologic_N = $1 if ($_ =~ /<shared_stage:pathologic_N.*?>(.*?)<\/shared_stage:pathologic_N>/);
    $vital_status =$1 if ($_ =~ /<clin_shared:vital_status.*?>(\S+?)<\/clin_shared:vital_status>/);
    if ($_ =~ /<clin_shared:days_to_death.*?>(\S+?)<\/clin_shared:days_to_death>/){
        $days_to_death = $1;
    }elsif ($_ =~ /<clin_shared:days_to_death/) {
        $days_to_death = "Not Available";
    }
}
close $fh;

print "$patient_id\t$age\t$gender\t$stage_pathologic_T\t$stage_pathologic_M\t$stage_pathologic_N\t$vital_status\t$days_to_death\n";

这里临床xml文件里面我发现其实有分两部分,前一半是旧的临床数据(或者说原始的?),后半部分是update的临床数据(我也不确定是否是这么理解),还好代码会自动以最后出现的数据为准,所以暂时是不没问题的,单个文件的结果如下,列名依次为病人样品id,年龄,性别,TMN,生存状态以及生存时间

TCGA-D8-A1Y1    80  FEMALE  T3  MX  N1a Dead    302

当然我们不可能每个样本都跑一次上述的代码,我们需要用perl自动读取目录下所有文件,然后依次处理文件,然后将结果输出到一个文件中,最后形成一个文件包含所有病人的临床信息

至于那个json的metadata文件,里面的临床信息肯定没有xml格式的完整,但是年龄,性别以及生存状态等基本信息还是有的,信息的提取可以参照TCGA数据下载方式小总结中所提到rjson这个R包来实现

Summary

对于TCGA的数据的处理还是基本摸索阶段,对于后续分析看了下资料差不多是各个组学的差异分析以及联合分析,或者还有生存曲线等方面的作图。这些东西跟常见的NGS二代测序的分析方法差不多也一致,所以TCGA可以看做一个数据库。

如果单纯的想看某个基因的TCGA中的表达情况以及生存曲线等信息,不少第三方网站都能提供信息。如果要个性化的分析,则需要会将TCGA获得的数据进行一些下游分析了。

后续的分析可能要看些文献才能了解一些”套路”了,需要学习的东西还是很多的。。。

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