纳米孔组装纠错

生信喵

       pilon 组装结果纠错，之前我们有介绍过。顺路就优化下测序数据吧
       pilon 是由 Broad Institute 研究所开发的纠错工具，输入原始拼接结果以及原始测序数据比对到拼接结果的 bam 文件即可。pilon 通过比对后的 bam 文件，可以识别拼接中非一致性的序列，包括单碱基的不同，小的 indel，大的 indel，后者空位 gap，以及错误拼接的区域。
       输入的 bam 可以来自于二代测序数据的比对，也可以来自于三代测序数据比对得到的 bam，bam 文件需要排序并建立索引。
      
       pilon 纠错流程图
       利用二代数据进行纠错

1. #三种拼接的结果
   
2. seqkit stat scaffolds.fasta ../megahit/final.contigs.fa ../../nanopore/flye/assembly.fasta
   
3. file                                format  type  num_seqs      sum_len  min_len    avg_len    max_len
   
4. scaffolds.fasta                     FASTA   DNA    227,948  296,416,270       56    1,300.4  1,188,045
   
5. ../megahit/final.contigs.fa         FASTA   DNA    170,308  275,083,965      200    1,615.2    481,000
   
6. ../../nanopore/flye/assembly.fasta  FASTA   DNA      1,151  203,303,214      510  176,631.8  5,785,557
   
7. 
   
8. #利用二代数据拼接好的16个contig进行纠错
   
9. mv 16.fasta complete.fasta
   
10. ASSEMBLY=complete.fasta
    
11. READS1=/share/home/xiehs/18.mags/2/illumina/ERR4007992_1.fastq.gz
    
12. READS2=/share/home/xiehs/18.mags/2/illumina/ERR4007992_2.fastq.gz
    
13. #对拼接结果建立索引
    
14. bwa-mem2 index $ASSEMBLY
    
15. 
    
16. #illumina比对排序建索引
    
17. echo "time bwa-mem2 mem -t 24 $ASSEMBLY $READS1 $READS2 -o illumina.sam 1>bwa.log 2>bwa.err" >bwa.sh
    
18. bsub -q fat -n 24 -o %J.log -e %J.err sh bwa.sh
    
19. #1.9版本 samtools 1.14-25-ga90ff1f
    
20. #conda install samtools=1.9
    
21. echo "time samtools sort -@ 36 -o illumina.sorted.bam illumina.sam" >samtools.sh
    
22. bsub -q fat -n 36 -o %J.log -e %J.err sh samtools.sh
    
23. samtools index illumina.sorted.bam
    
24. #利用pilon进行纠错
    
25. echo "time java -Xmx32G -jar ~/biosoft/pilon-1.24.jar --genome $ASSEMBLY --fix all --changes --frags illumina.sorted.bam --output pillon --outdir pillon_result --threads 24 --vcf 2> pilon.log" >pilon.sh
    
26. bsub -q fat -n 24 -o %J.log -e %J.err sh pilon.sh

复制代码

      --genome 提供输入参考基因组
       --frags 表示输入 < 1kb 的文库 BAM
       --jumps 输入 > 1kb 的文库 BAM
       -unpaired 输入非配对的 BAM。
       --output 表示输出的前缀
       --outdir 表示输出文件夹
       --changes 列出发生变化的部分,以 FASTA 形式保存
       --vcf 以 VCF 形式保存。
       --fix 声明对参考基因组做哪方面的改进，包括“snps”,”indels”,”gaps”,”local”, 默认是”all”