三代纳米孔宏基因组拼接

生信喵

背景
（1）纳米孔宏基因组测序的优点：纳米孔测序技术有着独一无二的优势，具体有8点，
       1.DNA/RNA直接测序
       2.可扩展性:便携式或台式测序仪
       3.真正的实时性
       4.10分钟文库制备
       5.无需对测序设备的资金投入
       6.高保真度
       7.超长读长
       8.对大基因组的高数据量测序
（2）利用Nanopore测序拼接的优势
       1、超长读长，解决重复序列问题；
       2、超长读长解决相位问题；
       3、超高通量测序，单张芯片最大可产出300G以上数据；

一、纳米孔测序宏基因组拼接下载数据
       从宏基因组测序中拼接细菌完成图，《Complete and validated genomes from a metagenome》2012 年加拿大 Northern Alberta 油砂尾矿池附近海藻细菌培养分离样品，2019 年重新培养提取。
      
       纳米孔测序宏基因组拼接案例
二、利用 metaFlye 拼接基因组
       flye 软件原来叫做 abruijn，是一款适合单分子测序的拼接软件，可以用于 Pacbio 与 nanopore测序数据。flye 适合拼接多种数据集，从小的细菌基因组到很大的哺乳动物基因组。flye 其实也是一个完整的拼接流程，首先原始测序数据作为输入，最后就可以直接输出经过纠错校正的基因组序列，同时给出拼接结果统计。flye 也有一个独立的模块 metaFlye 可以用来拼接宏基因组。
       软件官网：https://github.com/fenderglass/Flye
       查看数据：

1. seqkit stat ERR3994080_1.fastq.gz
   
2. file                   format  type  num_seqs        sum_len  min_len  avg_len  max_len
   
3. ERR3994080_1.fastq.gz  FASTQ   DNA    724,734  6,557,524,714        2  9,048.2  195,618

复制代码

      之前我们用二代数据megahit拼接完，最长才48k多，这个纳米孔原始数据最长就已经19k多了。
       使用案例
       就像它简介的名字一样，flye 使用起来也非常简单。直接输入 pacbio 或者 ont 的测序数据即可，可以是 fasta 格式，也可以是 fastq 格式文件，支持压缩格式。可以是原始数据，也可以是纠错后的数据，但需要通过选项参数分别设定。

1. conda install -y flye
   
2. echo "time flye --nano-raw ERR3994080_1.fastq.gz --genome-size 200M --out-dir flye --meta --threads 24 1>flye.log 2>flye.err" >flye.sh
   
3. bsub -q fat -n 24 -o %J.log -e %J.err sh flye.sh

复制代码

常用选项参数：
       --pacbio-raw ：输入原始 pacbio 数据
       --pacbio-corr ：输入纠错后的 pacbio 数据
       --nano-raw：输入原始 nanopore 数据
       --nano-corr ：输入原始 nanopore 数据
       --genome-size：预估基因组大小，用于评估覆盖深度
       --out-dir：输出结果文件路径
       --threads：cpu 线程数据
       --iterations：纠错迭代次数
       --min-overlap：最小 overlap 连接大小
       --meta: 拼接宏基因组数据
       --plasmids: 拼接质粒数据
输出结果

1. total 404M
   
2. drwxrwxr-x 2 xiehs xiehs 4.0K Dec 28 17:07 00-assembly/
   
3. drwxrwxr-x 2 xiehs xiehs 4.0K Dec 28 17:23 10-consensus/
   
4. drwxrwxr-x 2 xiehs xiehs 4.0K Dec 28 17:34 20-repeat/
   
5. drwxrwxr-x 2 xiehs xiehs 4.0K Dec 28 17:36 30-contigger/
   
6. drwxrwxr-x 2 xiehs xiehs 4.0K Dec 28 18:13 40-polishing/
   
7. -rw-rw-r-- 1 xiehs xiehs   92 Dec 28 18:13 params.json
   
8. -rw-rw-r-- 1 xiehs xiehs 328K Dec 28 18:13 assembly_graph.gv
   
9. -rw-rw-r-- 1 xiehs xiehs 202M Dec 28 18:13 assembly_graph.gfa
   
10. -rw-rw-r-- 1 xiehs xiehs 198M Dec 28 18:13 assembly.fasta
    
11. -rw-rw-r-- 1 xiehs xiehs  52K Dec 28 18:13 assembly_info.txt
    
12. -rw-rw-r-- 1 xiehs xiehs 2.8M Dec 28 18:13 flye.log
    

复制代码

最后结果目录中有三个文件比较重要。
       1、assembly.fasta ：最终拼接得到的基因组序列，fasta 格式。
       2、assembly_graph.{gfa|gv} ：拼接过程中用到的 repeat graph。
       3、assembly_info.txt：拼接结果统计信息，也可以自己单独使用 seqkit 工具统计。
总结
       1、软件使用简单；
       2、自带纠错功能；
       3、可以拼接宏基因组和质粒；

1. grep '\bY\b' assembly_info.txt
   
2. contig_1389     5785557 17      Y       N       1       *       1389
   
3. contig_765      5143764 118     Y       N       4       *       765
   
4. contig_308      4678198 21      Y       N       1       *       308
   
5. contig_227      4471145 116     Y       N       4       *       227
   
6. contig_777      4139087 13      Y       N       1       *       777
   
7. contig_37       3976090 20      Y       N       1       *       37
   
8. contig_1311     3954350 69      Y       N       2       *       1311
   
9. contig_307      3834242 83      Y       N       3       *       307
   
10. contig_1258     3764263 15      Y       N       1       *       1258
    
11. contig_1287     3740429 47      Y       N       2       *       1287
    
12. contig_1074     3729233 38      Y       N       1       *       1074
    
13. contig_282      3220921 18      Y       N       1       *       282
    
14. contig_268      3174875 237     Y       N       9       *       268
    
15. contig_517      3147655 429     Y       N       16      *       517
    
16. contig_1403     2881730 43      Y       N       1       *       1403
    
17. 以上15个可能是拼好的环状基因组
    
18. samtools faidx assembly.fasta contig_1389 >contig_1389.fasta
    
19. samtools faidx assembly.fasta contig_765 >contig_765.fasta
    
20. samtools faidx assembly.fasta contig_308 >contig_308.fasta
    
21. 可以与nt库比较出是哪个物种，网页ncbi比较也可以。
    
22. head contig_1389.fasta
    
23. 序列复制去 https://blast.ncbi.nlm.nih.gov/Blast.cgi 鉴定。
    
24. 当然也可以用centrifuge软件物种鉴定。

复制代码

生信喵

1. awk -F "\t" '{if ($4 == "Y" && $2 >= 300000) print $0}' assembly_info.txt | wc
   
2. 拼接结果筛选，原文大于300k，排除线粒体。
   
3. 目前16个，原文13个。
   
4. awk -F "\t" '{if ($4 == "Y" && $2 >= 300000) print $1}' assembly_info.txt | xargs
   
5.         contig_1389 contig_765 contig_308 contig_227 contig_777 contig_37 contig_1311 contig_307 contig_1258 contig_1287 contig_1074 contig_282 contig_268 contig_517 contig_1403 contig_795
   
6. samtools faidx assembly.fasta contig_1389 contig_765 contig_308 contig_227 contig_777 contig_37 contig_1311 contig_307 contig_1258 contig_1287 contig_1074 contig_282 contig_268 contig_517 contig_1403 contig_795 >16.fasta
   
7. 再用centrifuge软件物种鉴定。也可以复制去ncbi比对。

复制代码