21

肿瘤全外显子测序数据分析流程大放送

这个一个肿瘤外显子项目的文章发表并且公布的公共数据,我这里给出全套分析流程代码。只需要你肯实践,就可以运行成功。

PS:有些后起之秀自己运营公众号或者博客喜欢批评我们这些老人,一味的堆砌代码不给解释,恶意揣测我们是因为不懂代码的原理。我表示很无语,我写了3千多篇教程,如果一篇篇都重复提到基础知识,我真的做不到。比如下面的流程,包括软件的用法,软件安装,注释数据库的下载,我博客都说过好几次了,直播我的基因组系列也详细解读过,我告诉你去哪里学,你却不珍惜,不当回事,呵呵。

Continue reading

03

癌症基因的somatic mutation calling 流程的评价体系

癌症基因的somatic mutation calling 流程的评价体系

文章是:A comprehensive assessment of somatic mutation detection in cancer using whole-genome sequencing

WGS已经逐步走入临床,ICGC目前支持了74个国际项目,刻画了两万五千个癌症患者的基因组特性,希望能因此探究癌症的生物学机制。但对这些数据的分析缺乏严格论证的标准,不同的分析者有着自己独特的分析流程。

Continue reading

03

TCGA CNV全攻略

TCGA CNV全攻略

明白什么是CNV

对正常人来说,基因组应该是二倍体的,所以凡是测到非2倍体的地方都是CNV。但是CNV本身就是人群遗传物质多样性的体现,所以对癌症样本来说,是需要过滤掉正常人体内的germline的CNV,得到somatic的CNV。

Continue reading

02

TCGA计划的4个找somatic mutation的软件使用体验

TCGA计划的4个找somatic mutation的软件使用体验

体细胞突变(somatic mutation)是指患者某些组织或者器官后天性地发生了体细胞变异,虽然它不会遗传给后代个体,却可以通过细胞分裂,遗传给子代细胞。体细胞突变对肿瘤的发生发展有关键性的作用,并且它也是制定肿瘤癌症靶向治疗措施的关键所在。NGS使体细胞变异的检测更加全面,成本更低,在检测多种体细胞变异上具有很大的优势,但在使用过程中还存在着挑战:如样品降解、覆盖度不足、遗传异质性和组织污染(杂质)等问题。 为应对以上挑战,降低错误率,科学家采取了不同的算法和统计模型用于检测体细胞突变。目前最受欢迎的有Varscan、SomaticSniper、 Strelka 和MuTect2

Continue reading

24

一个标准的TCGA大文章应该做哪些数据?

很多人总是问我如何挖掘TCGA的数据,发文章!
可是他却连TCGA的数据是怎么来的都不知道,TCGA发了几十篇CNS大文章(自己测序的)了,每篇文章都有几百个左右的癌症样本的6种数据,这几年凑成了一万多个样本,都放在GDC里面可以任意下载。同时也出来了十几篇TCGA的数据挖掘大文章(主要包括亚型,driver mutation,假基因等新型研究领域)
那么一篇标准的一个标准的TCGA大文章应该自己测哪些数据?
其实稍微仔细浏览几篇文章就明白了,套路也是存在的,https://tcga-data.nci.nih.gov/docs/publications/
我们就以2013年发表在新英格兰杂志上面的Genomic and Epigenomic Landscapes of Adult De Novo Acute Myeloid Leukemia 为例子吧!

Continue reading

十二 28

TCGA表达数据的多项应用之4–求指定基因在指定癌症里面的表达量相关性矩阵,与所有的基因比较。

这个不出图,会给出TCGA里面涉及到的所有基因跟你指定的基因的表达量相关系数和P值,分别你一次性的看清楚你感兴趣的基因跟体内其它基因在该癌症种类的相关性,当然,相关非因果,请谨慎应用! Continue reading

十二 25

TCGA表达数据的多项应用之1–下载数据并且导入mysql

这个TCGA表达数据的多项应用系列帖子是应群里朋友的要求来写的,你们也可以继续提需求,我会接着写下去,其实从TCGA数据库里面下载到了数据之后,后面的所有分析都跟TCGA没有半毛钱关系了,大家要有这个想法,别三两句就问TCGA数据怎么分析,http://www.bio-info-trainee.com/?s=TCGA&submit=Search 本系列最后会形成一个shiny版本的交互式表达数据查询,处理,绘图,统计的网页APP。
我这里偷懒一下了,直接下载GEO里面的TCGA的表达数据,而不是去TCGA的官网里面下载:
它处理了目前(大概是2015年6月)TCGA收集的所有癌症样本的mRNA表达数据,并且统一处理成了count和RPKM两种表达量形式。 GEO地址:http://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE62944

Continue reading

19

用GISTIC多个segment文件来找SCNA变异

这个软件在TCGA计划里面被频繁使用者,用这个软件的目的很简单,就是你研究了很多癌症样本,通过芯片得到了每个样本的拷贝数变化信息,芯片结果一般是segment结果,可以解释为CNV区域,需要用GISTIC把样本综合起来分析,寻找somatic的CNV,并且注释基因信息。

有两个难点,一是在linux下面安装matlab工作环境,二是如何制作输入文件。

Continue reading

18

2016-TCGA数据挖掘系列文章之癌症男女有别

这是TCGA数据挖掘系列文章之一,是安德森癌症研究中心的Han Liang主导的,纯粹的生物信息学数据分析文章。
文章题目是:comprehensive characterization of molecular differences in cancer between male and female patients.
研究意义:癌症病人的性别对肿瘤发生,扩散的意义不言而喻。不仅仅是因为很多癌症本来就是有性别特异性,比如卵巢癌之于女性、前列腺癌之于男性。即使对于其它并非性别特异性的癌症种类,男女病人在肿瘤发生,扩散,以及治疗阶段的反应也大不一样。但是以前对这样分子机理研究的很有限,一般集中在某些性别相关的分子pattern,比如非小细胞肺癌女性患者的EGFR突变,但那些研究要么就局限于单一的基因,要么局限于单一的数据类型,或者研究单一的癌症。严重缺乏一个全面的,系统的分析癌症患者的性别差异。而且TCGA数据库的出现让这一个研究变成了可能,这也就是本文章的出现的原因。
数据挖掘的对象:
如表所示,涉及到13种癌症,TCGA的六种数据()都用上了,因为是2016年,所以数据量也比较全面了。

Continue reading

16

TCGA数据挖掘系列文章之-pseudogene假基因探究

这是TCGA数据挖掘系列文章之一,是安德森癌症研究中心的Han Liang主导的,纯粹的生物信息学数据分析文章。
TCGA数据库的数据量现在已经非常可观了,一万多的肿瘤样本数据,关于假基因的这篇文章是2014年发的,所以他们只研究了2,808个样本数据,也只涉及到7个癌症种类。
假基因是原来的能翻译成蛋白的基因经过各种突变导致丧失功能的基因。
比如
PTEN-->PTENP1
KRAS-->KRASP1
NANOG-->NANOGP1
很好理解,一般来说看到结尾是P1,等字眼的都是假基因,现在共有一万多假基因,我一般以http://www.genenames.org/cgi-bin/statistics (人类基因命名委员会)为标准参考。
文章主要做了6件事情。
一是重新定义及规范了假基因该研究什么就是把Yale Pseudogene database的假基因资源和GENCODE Pseudogene Resource的假基因资源结合起来,然后定义了一些过滤手段,具体流程如下。
1
二是下载了TCGA的那2,808个样本的RNA-seq的level2数据,也就是bam文件,重新提取关于假基因的表达数据。如果只是自己下载表达数据的话,关于假基因的定量并不准确,而且只有五百多个假基因。
当然,一般人没有条件下载RNA-seq的level2数据,所以想学习这个流程的话,直接下载表达矩阵吧。
Cancer type Number of nontumour samples Number of tumour samples Sequencing strategy Number of mappable reads Number of detectable pseudogenes
Breast invasive carcinoma 105 837 Paired-end 161 M 747
Kidney renal clear cell carcinoma 67 448 Paired-end 166 M 712
Lung squamous cell carcinoma 17 220 Paired-end 171 M 813
Ovarian serous cystadenocarcinoma 0 412 Paired-end 170 M 670
Glioblastoma multiforme 0 154 Paired-end 106 M 875
Colorectal carcinoma 0 228 Single-end 22 M 168
Uterine corpus endometrioid carcinoma 4 316 Single-end 26 M 181
第三件事是把假基因与其配对的野生型基因的表达数据做了相关性分析,一般来说,它们的相关性由下面三个原因决定。
(i) the sequence similarity between the pseudogene/gene pair;
(ii) the molecular mechanisms through which the pseudogene functions;
(iii) the detection sensitivity given the setting of RNA-seq experiments.
结论是不怎么相关,暗示着假基因虽然不编码蛋白产物,但仍然行使着某种功能。
第四件事是如果RNA-seq有正常对照的, 就做一样normal和tumor的差异分析,当然现在已经是都有了,在GSE62944可以下载所有的表达数据,专门提取假基因的表达数据做差异分析就好了。
但是差异分析的结果是, 没有什么现实意义。所以作者认为normal和tumor这样比较是不科学的,因为tumor本来就不应该按照组织来分类,而是应该按照TCGA的6种数据来分类()
In recent years, various ‘omic’ data, such as mRNA expression, microRNA expression, DNA methylation, somatic copy number alteration and protein expression, have been widely used to classify tumour samples into different molecular subtypes13, 14, 15, 16, 17, 18, 19.
2
第五件事就是把假基因表达数据的分类来跟其它几种分类形式作比较。
那些分类来源于以前的TCGA大文章:
48 in UCEC (endometrioid vs serous)23,
138 in LUSC (basal, classical, primitive and secretory)16,
71 in GBM (classical, mesenchymal, neural and proneural)24 and
547 in BRCA (PAM50 subtypes: luminal A, luminal B, basal-like, Her2-enriched and normal-like)2
文章就是:13141516171819.
3
最后就是做一些生存分析,讲一些好听的故事,比如说这样分类有利于精准医疗。

看起来还不错,值得大家学习一下,数据也都可以下载, 文章中提供了syn编号。

06

所有TCGA的maf格式somatic突变数据均可下载

如果你研究癌症,那么TCGA计划的如此丰富的公共数据你肯定不能错过,一般人只能获取到level3的数据,当然,其实一般人也没办法使用level1和level2的数据,毕竟近万个癌症样本的原始测序数据,还是很恐怖的,而且我们拿到原始数据,再重新跑pipeline,其实并不一定比人家TCGA本身分析的要好,所以我们直接拿到分析结果,就足够啦!

Continue reading

06

突变频谱探究mutation siganures

这也是对TCGA数据的深度挖掘,从而提出的一个统计学概念。文章研究了30种癌症,发现21种不同的mutation signature。如果理解了,就会发现这个其实蛮简单的,他们并不重新测序,只是拿已经有了的TCGA数据进行分析,而且居然是发表在nature上面!

研究了4,938,362 mutations from 7,042 cancers样本,突变频谱的概念只是针对于somatic 的mutation。一般是对癌症病人的肿瘤组织和癌旁组织配对测序,过滤得到的somatic mutation,一般一个样本也就几百个somatic 的mutation。

paper链接是:http://www.nature.com/nature/journal/v500/n7463/full/nature12477.html

从2013年提出到现在,已经有30种mutation siganures,在cosmic数据库有详细记录,更新见:http://cancer.sanger.ac.uk/cosmic/signatures
它的概念就是:根据突变上下文分成96类,然后每类突变的频率不一样画一个条形图,可视化展现。
mutation signature

Each signature is displayed according to the 96 substitution classification defined by the substitution class and sequence context immediately 3′ and 5′ to the mutated base. The probability bars for the six types of substitutions are displayed in different colours.
仔细看paper,还是蛮好理解的,自己写一个脚本就可以做这个分析了,前提是下载各个癌症的somatic mutation文件,一般是maf格式的,很多途径下载。
In principle, all classes of mutation (such as substitutions, indels, rearrangements) and any accessory mutation characteristic, for example, the sequence context of the mutation or the transcriptional strand on which it occurs, can be incorporated into the set of features by which a mutational signature is defined. In the first instance, we extracted mutational signatures using base substitutions and additionally included information on the sequence context of each mutation. Because there are six classes of base substitution—C>A, C>G, C>T, T>A, T>C, T>G (all substitutions are referred to by the pyrimidine of the mutated Watson–Crick base pair)—and as we incorporated information on the bases immediately 5′ and 3′ to each mutated base, there are 96 possible mutations in this classification. This 96 substitution classification is particularly useful for distinguishing mutational signatures that cause the same substitutions but in different sequence contexts.

很多癌症都发现了不止一种mutation signature,甚至高达6种,说明癌症之间差异还是蛮大的!
In most cancer classes at least two mutational signatures were observed, with a maximum of six in cancers of the liver, uterus and stomach. Although these differences may, in part, be attributable to differences in the power to extract signatures, it seems likely that some cancers have a more complex repertoire of mutational processes than others.
Most individual cancer genomes exhibit more than one mutational signature and many different combinations of signatures were observed
但是,我最后也没能绝对的界限是什么,因为总不能用肉眼来看每个突变频谱不一样吧?
The set of signatures will be updated in the future. This will include incorporating additional mutation types (e.g., indels, structural rearrangements, and localized hypermutation such as kataegis) and cancer samples. With more cancer genome sequences and the additional statistical power this will bring, new signatures may be found, the profiles of current signatures may be further refined, signatures may split into component signatures and signatures may be found in cancer types in which they are currently not detected.
分类会持续不断更新,随着更多的cancer type和样本加入,新的signature会被发现,现有的signature也可能会被重新定义,或者被分割成多个小的signature
28

使用可视化工具MutationMapper来看看基因上面突变的分布

how to generate lollipop diagrams ?
这个工具是网页版的,不需要下载,只需上传两列数据即可,就是基因上面的第几个氨基酸是如何突变的这样的信息
网页就可以把这些信息画在基因上面,而且注释到domain和cosmic数据库,挺好用的!
Hugo_Symbol HUGO symbol for the gene TP53
Protein_Change Amino acid change V600E
结果也很容易理解!
我输入了4个基因,网页就会对每个基因画一个图,其实ABCA1这个基因有6个突变,就都画在该基因上面了,基因突变分布在基因的哪个domain也看得一清二楚!
1
我的输入文件是这样的!
2
非常好用:
  • Support mutation data with annotated protein effects
  • Mutation diagram/lollipop view
  • Mutation table view
  • 3D structure view if available
而且pfam数据库本身也提供了这个功能:

Pfam provides an online tool to not only generate the domain information in JSON format, but to draw the lollipop diagram using javascript as well. They have more information here: http://pfam.xfam.org/help#tabview=tab9

IMHO, not as pretty as cBioPortal's but it gets you close to a solution.

EDIT / SHAMELESS PLUG: After seeing the data available and how easy it'd be, I made my own quick tool to fetch the data and draw the diagram for me in a style similar to cBioPortal - feel free to fork it and add features: https://github.com/pbnjay/lollipops

Example output (w/ labels per the comments)

3

甚至还有高手用D3.JS写了一个能实现同样需求的模块

We found ourselves in the same need, we wanted such a plot (JavaScript). Thus, I add our solution, Mutations Needle Plot. The library creates an SVG image (with D3), which then may be downloaded.

You will npm in order to be able to install & run the library.

Examples may be found in the snippets folder or also the index.html - The one displayed here below

类似的高手很多:
My colleague @SolenaLS recently asked me to write something like this: (uniprot+SVG+javascript : )http://lindenb.github.io/pages/uniprot/paintsvg.html
而且ensembl数据库也提供类似的功能
Ensemble also gives similar plot.

 

22

用TCGA数据做cox生存分析的风险因子(比例风险模型)

再次强调一下,R里面实现生存分析非常简单!
用my.surv <- surv(OS_MONTHS,OS_STATUS=='DECEASED')构建生存曲线。
用kmfit2 <- survfit(my.surv~TUMOR_STAGE_2009)来做某一个因子的KM生存曲线。用 survdiff(my.surv~type, data=dat)来看看这个因子的不同水平是否有显著差异,其中默认用是的logrank test 方法。
用coxph(Surv(time, status) ~ ph.ecog + tt(age), data=lung) 来检测自己感兴趣的因子是否受其它因子(age,gender等等)的影响。

Continue reading

14

用RNA-SeQC得到表达矩阵RPKM值

这个软件不仅仅能做QC,而且可以统计各个基因的RPKM值!尤其是TCGA计划里面的都是用它算的
一、程序安装
直接在官网下载java版本软件即可使用:http://www.broadinstitute.org/cancer/cga/tools/rnaseqc/RNA-SeQC_v1.1.8.jar
但是需要下载很多注释数据
clipboard

二、输入数据

clipboard

箭头所指的文件,一个都不少,只有那个rRNA.tar我没有用, 因为这个软件有两种使用方式,我用的是第一种
三、软件使用
软件的官网给力例子,很容易学习:
RNA-SeQC can be run with or without a BWA-based rRNA level estimation mode. To run without (less accurate, but faster) use the command:
java -jar RNASeQC.jar -n 1000 -s "TestId|ThousandReads.bam|TestDesc" -t gencode.v7.annotation_goodContig.gtf -r Homo_sapiens_assembly19.fasta -o ./testReport/ -strat gc -gc gencode.v7.gc.txt 
我用的就是这个例子,这个例子需要的所有文件里面,染色体都是没有chr的,这个非常重要!!!
代码如下:
 java -jar RNA-SeQC_v1.1.8.jar  \
-n 1000 \
-s "TestId|ThousandReads.bam|TestDesc" \
-t gencode.v7.annotation_goodContig.gtf \
-r ~/ref-database/human_g1k_v37/human_g1k_v37.fasta  \
-o ./testReport/ \
-strat gc \
-gc gencode.v7.gc.txt \
To run the more accurate but slower, BWA-based method :
java -jar RNASeQC.jar -n 1000 -s "TestId|ThousandReads.bam|TestDesc" -t gencode.v7.annotation_goodContig.gtf -r Homo_sapiens_assembly19.fasta -o ./testReport/ -strat gc -gc gencode.v7.gc.txt -BWArRNA human_all_rRNA.fasta
Note: this assumes BWA is in your PATH. If this is not the case, use the -bwa flag to specify the path to BWA
四、结果解读
运行要点时间,就那个一千条reads的测试数据都搞了10分钟!
出来一大堆突变,具体解释,官网上面很详细,不过,比较重要的当然是RPKM值咯,还有QC的信息
clipboard
TCGA数据里面都会提供由RNA-SeQC软件处理得到的表达矩阵!
Expression
  • RPKM data are used as produced by RNA-SeQC.
  • Filter on >=10 individuals with >0.1 RPKM and raw read counts greater than 6.
  • Quantile normalization was performed within each tissue to bring the expression profile of each sample onto the same scale.
  • To protect from outliers, inverse quantile normalization was performed for each gene, mapping each set of expression values to a standard normal.
软件的主页是: