使用Pilon对基因组进行polish

简介: 软件安装官方提供了编译好的jar包,方便使用wget https://github.com/broadinstitute/pilon/releases/download/v1.22/pilon-1.22.jarjava -Xmx16G -jar pilon-1.22.jar如果要顺利运行程序,要求JAVA > 1.7, 以及根据基因组大小而定的内存,一般而言是1M大小的基因对应1GB的内存。

软件安装

官方提供了编译好的jar包,方便使用

wget https://github.com/broadinstitute/pilon/releases/download/v1.22/pilon-1.22.jar
java -Xmx16G -jar pilon-1.22.jar

如果要顺利运行程序,要求JAVA > 1.7, 以及根据基因组大小而定的内存,一般而言是1M大小的基因对应1GB的内存。

总览

Pilon有如下作用

  1. 对初步组装进行polish
  2. 寻找同一物种不同株系间的变异,包括结构变异检测

他以FASTA和BAM文件作为输入,根据比对结果对输入的参考基因组进行提高,包括

  • 单碱基差异
  • 小的插入缺失(indels)
  • 较大的插入缺失或者block替换时间
  • 填充参考序列中的N
  • 找到局部的错误组装

最后它输出polish后的FASTA文件, 以及包含变异信息的VCF文件(可选)

分析流程

推荐使用PCR-free建库测序得到的Illumina paired-end数据,这样子避免了PCR-duplication,有效数据更多,也不需要在分析过程中标记重复。

下面步骤,假设你的组装文件为draft.fa, 质量控制后的illumina双端测序数据分别为read_1.fq.gzread_2.fq.gz

第一步:比对

bwa index -p index/draft draft.fa
bwa mem -t 20 index/draft read_1.fq.gz read_2.fq.gz | samtools sort -@ 10 -O bam -o align.bam
samtools index -@ 10 align.bam

第二步:标记重复(非PCR-free建库)

sambamba markdup -t 10 align.bam align_markdup.bam

第三步:过滤高质量比对的read

samtools view -@ 10 -q 30 align_markdup.bam > align_filter.bam
samtools index -@ 10 align_filter.bam

第三步:使用Pilon进行polish

MEMORY= #根据基因组大小而定
java -Xmx${MEMORY}G -jar pilon-1.22.jar --genome draft.fa --frags align_filer.bam \
    --fix snps,indels \
    --output pilon_polished --vcf &> pilon.log

一般Pilon迭代个2到3次就够了,所谓事不过三,过犹不及。

关于Pilon的一些参数说明:

  • --frags表示输入的是1kb以内的paired-end文库,--jumps表示 大于1k以上的mate pair文库, --bam则是让软件自己猜测
  • -vcf: 输出一个vcf文件,包含每个碱基的信息
  • --fix: Pilon将会处理的内容,基本上选snpsindels就够了
  • --variant: 启发式的变异检测,等价于--vcf --fix all,breaks, 如果是polish不要使用该选项
  • minmq: 用于Pilon堆叠的read最低比对质量,默认是0。

阅读日志输出

这个日志文件是标准输出而不是标准错误输出,不过保险起见用&>

最开始,Pilon会输出他的版本信息(如下示例),以及将会对基因组做的调整,

Pilon version 1.14 Sat Oct 31 14:30:00 2015 -0400
Genome: genome.fasta
Fixing snps, indels

其中Fixing后面的含义为:

  • "snps": 单碱基差异
  • "indels":小的indel的差异
  • "amb": 替换原有的N
  • "gaps": 填充基因组的gap
  • "local": 检测和修改错误组装
  • "all": 上述所有
  • "none": 不是上述的任何一种

接着Pilon会分染色体对BAM文件进行处理,根据BAM文件进行堆叠(pileup), 这个时候它会输出有效reads的深度,这里的有效reads包括未成对的read和正确成对的read。

Processing ctg1:1-5414473
frags align_mkdup.bam: coverage 19
Total Reads: 808985, Coverage: 19, minDepth: 5

从Pilon v1.4开始,Pilon还会输出基因组得到确认的碱基比例。

Confirmed 5403864 of 5414473 bases (99.80%)

后续是Pilon将会对原参考基因组做的一些调整的总体情况,如下表示纠正2个snp, 2个小的插入,4个缺失。

Corrected 2 snps; 0 ambiguous bases; corrected 2 small insertions totaling 12 bases, 4 small deletions totaling 6 bases

最后声明当前部分处理结束

Finished processing ctg1:1-5414473

如果,在--fix中选了gaps, 那么输出的内容还有如下内容。其中82048 -0 +276解释为在坐标82428处移除0个碱基,插入276个碱基。

# Attempting to fill gaps
fix gap: scaffold00001:82428-93547 82428 -0 +276 ClosedGap

参考资料

目录
相关文章
|
SQL 监控 IDE
实时开发IDE!数据开发效率开挂
杭州奥零数据科技有限公司成立于2023年,专注于数据中台业务,维护开源项目AllData并提供商业版解决方案。AllData提供数据集成、存储、开发、治理及BI展示等一站式服务,支持AI大模型应用,助力企业高效利用数据价值。
实时开发IDE!数据开发效率开挂
|
人工智能 算法 新能源
【2023高教社杯】A题 定日镜场的优化设计 问题分析及数学模型
本文介绍了2023年高教社杯数学建模竞赛A题的定日镜场优化设计问题,涉及问题分析和数学模型构建,旨在提高太阳能光热发电效率并实现电力系统的新能源转型。
806 1
【2023高教社杯】A题 定日镜场的优化设计 问题分析及数学模型
|
人工智能 自然语言处理 关系型数据库
不写一行代码,用MCP+魔搭API-Inference 搭建一个本地数据助手! 附所有工具和清单
还在为大模型开发的复杂技术栈、框架不兼容和工具调用问题头疼吗?MCP(Model Context Protocol servers)来拯救你了!它用统一的技术栈、兼容主流框架和简化工具调用的方式,让大模型开发变得简单高效。
2838 1
|
数据安全/隐私保护 Python
Python3,10行代码,给pdf文件去水印,再也不用花费冤枉钱了。
Python3,10行代码,给pdf文件去水印,再也不用花费冤枉钱了。
34794 1
Python3,10行代码,给pdf文件去水印,再也不用花费冤枉钱了。
|
搜索推荐 算法 C语言
【排序算法】八大排序(上)(c语言实现)(附源码)
本文介绍了四种常见的排序算法:冒泡排序、选择排序、插入排序和希尔排序。通过具体的代码实现和测试数据,详细解释了每种算法的工作原理和性能特点。冒泡排序通过不断交换相邻元素来排序,选择排序通过选择最小元素进行交换,插入排序通过逐步插入元素到已排序部分,而希尔排序则是插入排序的改进版,通过预排序使数据更接近有序,从而提高效率。文章最后总结了这四种算法的空间和时间复杂度,以及它们的稳定性。
1056 8
|
机器学习/深度学习 弹性计算 人工智能
滴滴派单算法_从算法模型思路到评估方案 - 详解
导读:说到滴滴的派单算法,大家可能感觉到既神秘又好奇,从出租车扬召到司机在滴滴平台抢单最后到平台派单,大家今天的出行体验已经发生了翻天覆地的变化,面对着每天数千万的呼叫,滴滴的派单算法一直在持续努力让更多人打到车,本篇文章会着重介绍我们是如何分析和建模这个问题,并且这其中面临了怎样的算法挑战,以及介绍一些我们常用的派单算法,这些算法能够让我们不断的提升用户的打车确定性。
11158 123
滴滴派单算法_从算法模型思路到评估方案 - 详解
|
C语言 C++ Linux
查看gcc/g++默认include路径
转自:http://gcc.gnu.org/ml/gcc-help/2007-09/msg00205.html   `gcc -print-prog-name=cc1plus` -v `g++ -print-prog-name=cc1plus` -v   例如,CentOS 6.
2366 0
|
数据安全/隐私保护
|
异构计算 Windows
使用MATLAB连接USRP
使用MATLAB连接USRP
1019 2
技术笔记:MOS管常用效应
技术笔记:MOS管常用效应
1802 0