请选择 进入手机版 | 继续访问电脑版

生信人

找回密码
立即注册
搜索
热搜: 活动 交友 discuz
发新帖

0

收听

12

听众

246

主题
发表于 2022-12-29 22:33:33 | 查看: 8574| 回复: 2
本帖最后由 生信喵 于 2022-12-29 23:19 编辑

       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. #利用二代数据拼接好的16个contig进行纠错
  8. mv 16.fasta complete.fasta
  9. ASSEMBLY=complete.fasta
  10. READS1=/share/home/xiehs/18.mags/2/illumina/ERR4007992_1.fastq.gz
  11. READS2=/share/home/xiehs/18.mags/2/illumina/ERR4007992_2.fastq.gz
  12. #对拼接结果建立索引
  13. bwa-mem2 index $ASSEMBLY

  14. #illumina比对排序建索引
  15. echo "time bwa-mem2 mem -t 24 $ASSEMBLY $READS1 $READS2 -o illumina.sam 1>bwa.log 2>bwa.err" >bwa.sh
  16. bsub -q fat -n 24 -o %J.log -e %J.err sh bwa.sh
  17. #1.9版本 samtools 1.14-25-ga90ff1f
  18. #conda install samtools=1.9
  19. echo "time samtools sort -@ 36 -o illumina.sorted.bam illumina.sam" >samtools.sh
  20. bsub -q fat -n 36 -o %J.log -e %J.err sh samtools.sh
  21. samtools index illumina.sorted.bam
  22. #利用pilon进行纠错
  23. 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
  24. 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”

本帖子中包含更多资源

您需要 登录 才可以下载或查看,没有账号?立即注册

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

QQ|Archiver|手机版|小黑屋|生信人

GMT+8, 2024-3-28 22:25 , Processed in 0.138556 second(s), 21 queries .

Powered by Discuz! X3.5

© 2001-2024 Discuz! Team.

快速回复 返回顶部 返回列表