06

用SomaticSignatures包来解析maf突变数据获得mutation signature

mutation signature这个概念提出来还不久,我看了看文献,最早见于2013年的一篇nature文章,主要是用来描述癌症患者的somatic mutation情况的。

首先要自己分析癌症样本数据,拿到somatic mutation,TCGA计划发展到现在已经有非常多的somatic mutation结果啦,大家可以自行选择感兴趣的癌症数据拿来研究,解析一下mutation signature 。

我这里给大家推荐一个工具,是R语言的Bioconductor系列包中的一个,SomaticSignatures

其实它的说明书写的非常详细了已经,如果你理解了mutation signature的概念,很容易用那个包,其实你自己写一个脚本也是非常任意的,就是根据mutation的位置在基因组中找到它的前后一个碱基,然后组成三碱基突变模式,最后统计一下那96种突变模式的分布状况!

我这里简单讲一讲这个包如何用吧!

首先下载并加载几个必须的包:

library(SomaticSignatures)  ## 程序
library(SomaticCancerAlterations) ## 自带测试数据
library(BSgenome.Hsapiens.1000genomes.hs37d5)  ## 我们的参考基因组
library(VariantAnnotation)
## 这个对象很重要: GRanges class of the GenomicRanges package
##其中SomaticCancerAlterations这个包提供了测试数据,来自于8个不同癌症的外显子测序的项目。
sca_metadata = scaMetadata()
###可以查看关于这8个项目的介绍,每个项目都测了好几百个样本。但是我们只关心突变数据,而且只关心somatic的突变数据。
sca_data = unlist(scaLoadDatasets())

然后根据突变数据做好一个GRanges对象,这个可以看我以前的博客

sca_data$study = factor(gsub("(.*)_(.*)", "\\1", toupper(names(sca_data))))
sca_data = unname(subset(sca_data, Variant_Type %in% "SNP"))
sca_data = keepSeqlevels(sca_data, hsAutosomes())
## 这个对象就是我们软件的输入数据
sca_vr = VRanges(
    seqnames = seqnames(sca_data),
    ranges = ranges(sca_data),
    ref = sca_data$Reference_Allele,
    alt = sca_data$Tumor_Seq_Allele2,
    sampleNames = sca_data$Patient_ID,
    seqinfo = seqinfo(sca_data),
    study = sca_data$study
)
## 这里还可以直接用readVcf或者readMutect 来读取本地somatic mutation文件
## 提取突变数据,并且构造成一个Range对象。
sca_vr
###可以简单看看每个study都有多少somatic mutation
sort(table(sca_vr$study), decreasing = TRUE)
    LUAD   SKCM   HNSC   LUSC   KIRC    GBM   THCA     OV
   208724 200589  67125  61485  24158  19938   6716   5872
##用mutationContext函数来根据Range对象和下载好的参考基因组文件来获取突变的上下文信息。
sca_motifs = mutationContext(sca_vr, BSgenome.Hsapiens.1000genomes.hs37d5)
head(sca_motifs)
##可以看到Range对象,增加了两列:alteration        context
## 接下来根据做好的上下文突变数据矩阵来构建 the matrix MM of the form {motifs × studies}
sca_mm = motifMatrix(sca_motifs, group = "study", normalize = TRUE)
## 根据96种突变的频率,而不是次数来构造矩阵
head(round(sca_mm, 4))
## 然后直接画出每个study的Mutation spectrum 图
plotMutationSpectrum(sca_motifs, "study")
 mutation spectrum
## 还要把spectrum分解成signature!!
## 这个包提供了两种方法,分别是NMF和PCA
n_sigs = 5
sigs_nmf = identifySignatures(sca_mm, n_sigs, nmfDecomposition)
sigs_pca = identifySignatures(sca_mm, n_sigs, pcaDecomposition)
##还提供了很多函数来探索:signatures, samples, observed and fitted.
需要我们掌握的是assessNumberSignatures,用来探索我们到底应该把spectrum分成多少个signature
n_sigs = 2:8
gof_nmf = assessNumberSignatures(sca_mm, n_sigs, nReplicates = 5)
gof_pca = assessNumberSignatures(sca_mm, n_sigs, pcaDecomposition)
plotNumberSignatures(gof_nmf) ## 可视化展现
## 接下来可视化展现具体每个cancer type里面的各个个体在各个signature的占比
library(ggplot2)
plotSignatureMap(sigs_nmf) + ggtitle("Somatic Signatures: NMF - Heatmap")
plotSignatures(sigs_nmf) + ggtitle("Somatic Signatures: NMF - Barchart")
plotObservedSpectrum(sigs_nmf)
plotFittedSpectrum(sigs_nmf)
plotSampleMap(sigs_nmf)
plotSamples(sigs_nmf)
同理,PCA的结果也可以同样的可视化展现:
plotSignatureMap(sigs_pca) + ggtitle("Somatic Signatures: PCA - Heatmap")
plotSignatures(sigs_pca) + ggtitle("Somatic Signatures: PCA - Barchart")
plotFittedSpectrum(sigs_pca)
plotObservedSpectrum(sigs_pca)
mutation signature NMF
值得一提的是,所有的plot系列函数,都是基于ggplot的,所以可以继续深度定制化绘图细节。
p = plotSamples(sigs_nmf)
## (re)move the legend
p = p + theme(legend.position = "none")
## (re)label the axis
p = p + xlab("Studies")
## add a title
p = p + ggtitle("Somatic Signatures in TGCA WES Data")
## change the color scale
p = p + scale_fill_brewer(palette = "Blues")
## decrease the size of x-axis labels
p = p + theme(axis.text.x = element_text(size = 9))
###当然,对上下文突变数据矩阵也可以进行聚类分析
clu_motif = clusterSpectrum(sca_mm, "motif")
library(ggdendro)
p = ggdendrogram(clu_motif, rotate = TRUE)
p
## 最后,由于我们综合了8个不同的study,所以必然会有批次影响,如果可以,也需要去除。
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