跟着PNAS学数据分析:MUM&Co软件基于基因组做结构变异检测

2023-10-20 14:09:39 浏览数 (1)

MUM&Co软件对应的论文

MUM&Co: accurate detection of all SV types through whole-genome alignment

https://academic.oup.com/bioinformatics/article/36/10/3242/5756209?login=false

工具的github主页 https://github.com/SAMtoBAM/MUMandCo

工具是一个shell脚本,需要安装MUMmer4

能够检测的变异类型有 Deletions, insertions, tandem duplications and tandem contractions (>=50bp & <=150kb) Inversions (>=1kb) , translocations (>=10kb)

可可树基因组结构变异对应的论文

Genomic structural variants constrain and facilitate adaptation in natural populations of Theobroma cacao, the chocolate tree

https://www.pnas.org/doi/10.1073/pnas.2102914118

巧克力树结构变异PNAS.pdf

这篇论文里有对应的分析代码,github主页是

https://github.com/thamala/cacaoSV/tree/main

31个基因组,二倍体,做的是单倍型的组装,每个单倍型分别与参考基因组比对做结构变异检测,论文里提供了工具可以报两个单倍型分别检测sv的结果合并成一个二倍体的结果。之前看人类的结构变异分析的论文里也提供了类似的工具

用拟南芥的数据尝试一下这个工具

首先是比对做结构变异检测

代码语言:javascript复制
time bash ../MUMandCo-master/mumandco_v3.8.sh -r An1.fa -q C24.fa -g 120000000 -o C24 -t 8
time bash ../MUMandCo-master/mumandco_v3.8.sh -r An1.fa -q Eri.fa -g 120000000 -o Eri -t 8
time bash ../MUMandCo-master/mumandco_v3.8.sh -r An1.fa -q Kyo.fa -g 120000000 -o Kyo -t 8
time bash ../MUMandCo-master/mumandco_v3.8.sh -r An1.fa -q Sha.fa -g 120000000 -o Sha -t 8

用拟南芥的数据,基因组大小在120M左右,5分钟一个样本

输出的文件内容

image.png

这里有一个vcf文件,有一个tsv文件,tsv文件是下面合并多个样本需要用到的数据

合并多个样本的工具是mumco2vcf.c

首先是编译

代码语言:javascript复制
gcc mumco2vcf.c -o mumco2vcf -lm

mumco2vcf 这个是可执行文件

将需要合并的样本名整理到一个文本文件里,一行一个样本,如果是一个样本的两个单倍型,放在同一行,用制表符分隔

比如我这里是四个样本,准备的示例文件

代码语言:javascript复制
ls */*all.tsv > vcf.list


C24_output/C24.SVs_all.tsv
Eri_output/Eri.SVs_all.tsv
Kyo_output/Kyo.SVs_all.tsv
Sha_output/Sha.SVs_all.tsv

合并运行的命令

代码语言:javascript复制
../cacaoSV-main/mumco2vcf -i vcf.list > out.vcf

还有好几个参数

Usage: -i [file] File listing samples to combine. Haplotypes are in one row, separated by a tab. -t [str] either INS, DEL, DUP, TRA, or INV. (optional) -l [int]-[int] minimum and maximum length of the SVs. (optional) -o [int] overlap in base pairs to combine variants. (optional)

如果是单倍型

vcf.list的文件内容

代码语言:javascript复制
C24_output/C24.SVs_all.tsv      Eri_output/Eri.SVs_all.tsv
Kyo_output/Kyo.SVs_all.tsv      Sha_output/Sha.SVs_all.tsv

合并的命令是一样的

欢迎大家关注我的公众号

小明的数据分析笔记本

小明的数据分析笔记本 公众号 主要分享:1、R语言和python做数据分析和数据可视化的简单小例子;2、园艺植物相关转录组学、基因组学、群体遗传学文献阅读笔记;3、生物信息学入门学习资料及自己的学习笔记!

image.png

0 人点赞