搜索
查看: 2235|回复: 1

【直播】我的基因组 45:SNV突变(6种)频谱的制作

[复制链接]

103

主题

133

帖子

860

积分

版主

Rank: 7Rank: 7Rank: 7

积分
860
发表于 2017-2-7 19:55:52 | 显示全部楼层 |阅读模式
【直播】我的基因组 45:SNV突变(6种)频谱的制作

受热心读者的委托,特意为他讲解一下SNV突变(6种)频谱的制作,同时欢迎大家留言其它需求!

突变频谱呢,就是对含有SNV的VCF格式的文件进行一个统计。



全基因组SNP突变可以分成6类(C>A, C>G, C>T, A>C, A>G, A>T)。肯定会有人问为什么是六类?

以A:T>C:G为例,此种类型SNP突变包括A>C和T>G。由于测序数据既可比对到参考基因组的正链,也可比对到参考基因组的负链,当T>C类型突变出现在参考基因组正链上,A>G类型突变即在参考基因组负链的相同位置,所以将T>C和A>G划分成一类,换句话说我们只考虑正链的突变形式,参考碱基只允许有C或者T,因为它们等价于G或者A。所以全基因组SNP突变可以分成这6类。

很明显,我们只需要考虑VCF文件的第4,5行即可!

cat realign.vcf |grep -v INDEL |grep -v "^#" |cut -f 1-5 |head


f45a5023432d4d6bc3d3f6b6b71e3a54.png


cat realign.vcf |grep -v INDEL |grep -v "^#" |cut -f 4,5|sort |uniq -c |grep -v ","

我们过滤掉了多种变异形式的SNV,比如T,突变成G或者C!最后的结果如下:


d6cc515288b0c91877e74aa8a160475e.png



一般来说,是要可视化一下的,我用R语言的ggplot来画一个呗~

  • dat <- data.frame(type=c('C>A(G>T)','C>T(G>A)','C>G(G>C)','T>A(A>T)','T>G(A>C)','T>C(A>G)'),
  • counts=c(180515+181567,698322+697568,184176+185144,148387+148580,177215+177415,676816+675821)
  • )
  • library(ggplot2)
  • p=ggplot(dat,aes( x=type,y=counts))+geom_bar(stat="identity")
  • print(p)



当然,mutation spectrum这个画图代码只能出一个最简单的条形图,如果你想达到下面的效果,需要学习ggplot啦!


画条形图请参考:[url=]http://docs.ggplot2.org/0.9.3.1/geom_bar.html[/url]


c96234fbabe4a2aa52587c15d10df3a2.png



如果要区分染色体,可以重新考虑第1行,拿去可视化!
[url=]http://www.bio-info-trainee.com/1619.html[/url]

文:Jimmy

图文编辑:吃瓜群众






上一篇:【直播】我的基因组 44:比对文件画profile和heatmap图
下一篇:【直播】我的基因组46:SNV突变(96种)频谱的制作
基因组,转绿组,肿瘤信息,生物统计,Python, Linux.
回复

使用道具 举报

0

主题

3

帖子

40

积分

新手上路

Rank: 1

积分
40
发表于 2017-2-12 16:26:01 | 显示全部楼层
顺带学习了一下linux命令!
回复 支持 反对

使用道具 举报

您需要登录后才可以回帖 登录 | 立即注册

本版积分规则

QQ|手机版|小黑屋|生信技能树 ( 粤ICP备15016384号  

GMT+8, 2019-11-20 13:28 , Processed in 0.032001 second(s), 27 queries .

Powered by Discuz! X3.2

© 2001-2013 Comsenz Inc.