混合拼接

生信喵

背景
      在当前测序市场条件下，一代二代三代测序平台共存，这是因为每一个测序平台都不是完美的，都有优缺点，例如一代测序准确性高，读长长，但是相应的通量低，价格贵。二代测序显著特点是通量大，准确性高，价格便宜，但是读长短，而三代测序读长长，通量高，但是错误率高，这就需要我们能够根据每种数据的特点，充分利用每一种测序平台数据，得到最优的方案。例如当前基因组拼接中，尤其对于一些大型的基因组，往往混合多种测序平台数据进行拼接，达到最佳的拼接效果。
一、混合拼接方案
      我们对同一样品，分别进行了二代 illumina 双末端测序，三代单分子 pacbio 测序以及三代纳米孔 nanopore 测序数据。这些数据有多种分析方案，例如以下方案：
      方案一：二代 illumina 数据单独拼接
      方案二：三代数据单独拼接
      方案三：二代 illumina 数据为主+三代测序连接+补洞；
      方案四：三代测序单独拼接，利用二代数据连接+纠错；
      由于另种数据拼接结果肯定要优于只使用一种数据，因此，这里只需要比较方案三与方案四。也就是以哪种测序数据为主，哪种方案为辅。
当前大型基因组最佳方案：
      nanopore搭建框架
      illumina补洞
      Hi-C挂载基因组

二、Hi-C
      高通量染色体构象捕获（High-throughput chromosome conformation capture, Hi-C）则是一种可用于研究基因组三维结构的创新技术。自2009年美国Job Dekker 研究组开发出了Hi-C 技术，现已广泛应用于细菌、酵母、拟南芥、棉花、水稻、玉米、小鼠及人类基因组研究中。
      参考文献：PMID：24185095
      
Hi-C实验步骤
      

实验步骤示意图（Lieberman-Aiden E. et al. Science. 2009）
      1、使用多聚甲醛处理细胞，固定DNA的构象；
      2、裂解细胞后，使用限制性内切酶处理交联的 DNA，产生粘性末端；
      3、DNA末端补平修复，并同时引入生物素，标记寡核苷酸末端；
      4、使用DNA连接酶连接DNA 片段；
      5、蛋白酶消化解除与 DNA 的交联状态，纯化 DNA并随机打断至长度为 300~500bp 的片段；
      6、使用亲和素磁珠捕获标记的DNA ，进行二代建库测序。
      采用高通量测序获得的互作数据具有两个规律，一是染色体内的互作强度大于染色体间的互作；二是同一染色体内互作强度随线性距离的增加而减弱。Hi-C技术可研究基因组空间调控机制，还可将Hi-C数据用于基因组组装中，对三代测序组装获得的contig序列进行scaffolding，达到染色体水平。

Hi-C辅助染色体定位
      
      基于Hi-C数据上述规律，首先，基于Hi-C互作强度对Hi-C read聚类，判断scaffolds/Contigs分组；针对group内部互作强度，确定scaffolds/Contigs前后相邻关系；最后进一步确定相邻scaffolds/Contigs方向，最后获得染色体水平参考基因组序列。

三、不同基因组拼接方案比较
1 单独使用 illumina 拼接

1. spades.py --isolate -o ill -t 24 -1 /share/home/xiehs/05.assembly/data/clean.1.fq.gz -2 /share/home/xiehs/05.assembly/data/clean.2.fq.gz 1>spades.log 2>spades.err

复制代码

2 利用 illumina 数据+pacbio 数据拼接

1. spades.py --isolate -o ill_pac -t 24 -1 /share/home/xiehs/05.assembly/data/clean.1.fq.gz -2 /share/home/xiehs/05.assembly/data/clean.2.fq.gz --pacbio /share/home/xiehs/05.assembly/data/pacbio.fastq.gz 1>ill_pac.log 2>ill_pac.err

复制代码

3 利用 illumina 数据+nanopore 数据拼接

1. spades.py --isolate -o ill_ont -t 24 -1 /share/home/xiehs/05.assembly/data/clean.1.fq.gz -2 /share/home/xiehs/05.assembly/data/clean.2.fq.gz --nanopore /share/home/xiehs/05.assembly/data/nanopore.fastq.gz 1>ill_ont.log 2>ill_ont.err

复制代码

     比较前三种

1. seqkit seq -m 500 ill/scaffolds.fasta | seqkit stat
   
2. seqkit seq -m 500 ill_ont/scaffolds.fasta | seqkit stat
   
3. seqkit seq -m 500 ill_pac/scaffolds.fasta | seqkit stat

复制代码

4 unicycler 混合拼接

1. conda activate unicycler
   
2. echo "unicycler -1 /share/home/xiehs/05.assembly/data/clean.1.fq.gz -2 /share/home/xiehs/05.assembly/data/clean.2.fq.gz --long /share/home/xiehs/05.assembly/data/nanopore.fastq.gz -o unicycler -t 12 >unicycler.log 2>unicycler.err" > unicycler.sh

复制代码

     结果统计

1. seqkit stat assembly.fasta
   
2. file            format  type  num_seqs    sum_len  min_len    avg_len    max_len
   
3. assembly.fasta  FASTA   DNA          6  5,691,015    3,478  948,502.5  5,311,635

复制代码

四、quast 比较不同拼接方案
      利用 quast 比较不同拼接方案的拼接结果，选择最优的拼接结果。

1. quast.py -r MGH78578.fasta ill.fa ill_pac.fa ill_ont.fa uni.fa -o quast
   

复制代码

     结果report.html
      
      报告可以看出unicycler的结果最接近参考序列。