本文首發(fā)于“生信補(bǔ)給站” https://mp.weixin.qq.com/s/k7yzk9hPX3Bi-ohAo83ZYw
還有其他 R統(tǒng)計(jì) 繪圖 生信的干貨，也許有需要的呢？
Mutational Signatures 首次出現(xiàn)在2013年的nature文章Signatures of mutational processes in human cancer中（https://www.nature.com/articles/nature12477）。**將mutation位置加上前后一個(gè)堿基，構(gòu)成三堿基模式，然后統(tǒng)計(jì)96（6 * 4 * 4）種突變組合的情況。
好奇為什么是96種的，可以查一下文獻(xiàn)。
本文介紹如何利用deconstructSigs-R包進(jìn)行mutation signature分析。
一 準(zhǔn)備R包，數(shù)據(jù)
#install.packages("deconstructSigs")
library(deconstructSigs)
#讀入數(shù)據(jù)
head(sample.mut.ref)
Sample chr pos ref alt
1 1 chr1 905907 A T
2 1 chr1 1192480 C A
3 1 chr1 1854885 G C
4 1 chr1 9713992 G A
5 1 chr1 12908093 C A
6 1 chr1 17257855 C T
class(sample.mut.ref)
## [1] "data.frame"
只需要將自己的數(shù)據(jù)整理成以上五列（ID，chr，pos，ref，alt ）信息即可，如果是TCGA中的MAF文件也是很好提取的。
二 mut.to.sigs.input構(gòu)建輸入文件
使用 mut.to.sigs.input 函數(shù)，構(gòu)建計(jì)算signature的輸入文件，得到每個(gè)樣本的96種三堿基類型。
# Convert to deconstructSigs input
sigs.input <- mut.to.sigs.input(mut.ref = sample.mut.ref,
sample.id = "Sample",
chr = "chr",
pos = "pos",
ref = "ref",
alt = "alt")
注：這一步也許會(huì)提示沒(méi)有XX包，按照要求下載指定R包即可（也許是數(shù)據(jù)庫(kù)，耐心安裝）。
#查看結(jié)果信息
dim(sigs.input)
#[1] 2 96
head(t(sigs.input)) #只有兩個(gè)sample：“1”和“2”
1 2
A[C>A]A 9 1
A[C>A]C 7 1
A[C>A]G 5 0
A[C>A]T 7 0
C[C>A]A 10 3
C[C>A]C 18 2
以上就得到了sample.mut.ref文件中的每一個(gè)sample的96種三堿基類型的結(jié)果了。
三 推斷signature的組成
# Determine the signatures contributing to the two example samples
sample_1 = whichSignatures(tumor.ref = sigs.input,
signatures.ref = signatures.cosmic,
sample.id = 1,
contexts.needed = TRUE,
tri.counts.method = 'default')
其中：
tumor.ref：每個(gè)sample的96種三堿基突變序列
signatures.ref：已知的signatures參考文件，可選signatures.nature2013和signatures.cosmic
sample.id：對(duì)應(yīng)tumor.ref文件中的樣本名
contexts.needed ：是否需要突變上下文
tri.counts.method：三核酸序列標(biāo)準(zhǔn)化方式，默認(rèn)“default” 不進(jìn)行標(biāo)準(zhǔn)化 ；或者選擇exome,genome,exome2genome,genome2exome 來(lái)限定區(qū)域。
3.2 查看返回結(jié)果
#查看結(jié)果
class(sample_1)
#查看權(quán)重結(jié)果
sample_1$weights
#輸出tumor的三堿基序列百分比
sample_1$tumor
#三堿基序列百分比 * 權(quán)重
sample_1$product
whichSignatures會(huì)輸出5個(gè)元素的list文件：
weights -- data frame containing the weights assigned to each of the k signatures of the input signatures matrix
tumor -- matrix of the trinucleotide contexts for the tumor sample used as input
product -- matrix obtained when the tumor matrix is **multiplied by the assigned weights **
diff -- matrix representing the difference between the tumor matrix and product matrix
unknown -- numeric weight not assigned to any of the input signatures
3.3 指定signature權(quán)重
通過(guò)associated參數(shù)指定參與計(jì)算的signature
sample_1.associate = whichSignatures(tumor.ref = sigs.input,
signatures.ref = signatures.cosmic,
sample.id = 1,
associated = c("Signature.1","Signature.22"),
contexts.needed = TRUE,
tri.counts.method = 'default')
sample_1.associate$weights
3.4 設(shè)定signature的閾值
通過(guò)signature.cutoff設(shè)定閾值，小于此值的為0
sample_1.cutoff = whichSignatures(tumor.ref = sigs.input,
signatures.ref = signatures.cosmic,
sample.id = 1,
contexts.needed = TRUE,
signature.cutoff = 0.08 ,
tri.counts.method = 'default')
sample_1.cutoff$weights
四 plotSignatures 可視化
# Plot example
plot_example <- whichSignatures(tumor.ref = sigs.input,
signatures.ref = signatures.cosmic,
sample.id = 1)
# Plot output
plotSignatures(plot_example, sub = 'example')
查看sample1的signature的組成情況，就是上面plot_exampletumor , plot_example$product 的結(jié)果可視化。
參考資料：
https://github.com/raerose01/deconstructSigs
◆ ◆ ◆ ◆ ◆
精心整理（含圖版）|你要的全拿走！有備無(wú)患 （R統(tǒng)計(jì)，ggplot2繪圖，生信圖形可視化匯總）
【覺(jué)得不錯(cuò)，右下角點(diǎn)個(gè)“在看”，期待您的轉(zhuǎn)發(fā)，謝謝！】
