sequence assembly

1. DNA测序技术

第一代：sanger

第二代：高通量 illumina， 454测序，等

2.5代：lon torrent/proton 测序

第三代：pacbio

2.mate-pair（大长度） 文库的构建

≥2k 片段→加接头→环化→打断→短片段→回收带标记的片段对较短的片段进行测序

3. 各测序平台特性

Illumina

reads：49-300bp
错误率 1%左右，3‘端会高一点
pair-end文库（170-40kb）

数据格式 4行：（fastq）

1. reads id，@开头
2. 序列信息，N代表不确定
3. 质量值的id，一般和 reads 的 id 相同
4. 质量值的信息（ASCII 码-64 = 质量值）

成对的 reads 对应的行号是一样的。

454

reads：700bp（peak），平均 500
error and indel 错误
pair-end 1kb-20kb

数据格式：

1. 序列信息：fasta格式，> 开头 id
2. 质量信息：存放于另一个文件中，> 号开头，直接用数值表示

4. 组装的基本思路

策略1：

分级组装：人类基因组
缺点：实验环节较慢，需要图谱，消耗高
优点：算法简单，质量高

策略2：

鸟枪法（whole-genome shotgun）：
优点：实验简单，花费低
缺点：计算能力要求，质量低

5. 基本概念

overlap：重合部分
contig：组装的初级序列
pair-end：中间有未知的区域
scaffold：

6. 基因组大小的估计

流式细胞仪

Kmer 分析

根据近源物种估计，得到基础值，得到一些测序数据，利用Kmer估计

1. 泊松分布 poisson distribution

测序深度：基因组中每个碱基被测序的频率。一般情况下，个体测序：30 层

假设基因组大小为G，随机测序，read 的长度为L，当生成一个read 的时候，一个固定碱基被测到的概率是 L/G， reads 的总量为n_r.

L/G 非常小，而 n_r 非常大，每个碱基被测到的频率的期望：
d_n = L/G * n_r
它服从泊松分布（n 很大而 p 很小的二项分布）。

Kmer：固定长度的核苷酸序列，奇数（避免正反链一致，中间的肯定不一样）

假设基因组中 K 是独一无二的，则 G 大小的基因组中的 Kmer 的数量也为G，reads 的长度为 L，则一个固定的 Kmer 被测到的概率为（L- K + 1）/G（因为L上可以有L-K+1个Kmer）。

L-K+1/G 很小，而 n_r 非常大，每个Kmer被测到的频率的期望和总个数：
d_k = (L-K+1)/G * n_r
n_k = (L-K+1) * n_r
可推导出：
G = n_k/d_k
d_n/d_k = L/(L-K+1)

image.png

图片为，reads为100时，Kmer长度为25，在d_n 为 40时， d_k 为30.4。
100/(100-25+1) * 40 = 30.4
基因组大小 G = n_k/d_k

7. Kmer 的其他应用

下图为一个比较简单的基因组的Kmer分布

image.png

下面较复杂：

最前面：测序错误造成。
两个峰值：一个峰高，一个峰低（高杂合基因组造成，深度正好差一倍）。
尾部：由于重复片段造成的。

image.png

下面：

同样两个峰，80处有个小峰，基因组带有高重复率，但是杂合度不高。

image.png

Kmer在组装上的应用

组装基于 Kmer 而非 reads

长度为17的Kmer，相邻有16bp是相同的。

read 上的Kmer数： L-K+1

8. 基于二代测序结果的组装

一些名词：

image.png

Illumina

原理：根据Kmer（Kmer 图）

软件：SOAPdenovo（华大）、velvet ABySS、ALLPATH-LG

以SOAPdenovo为例：

1. 构建图（De bruijn graph）

   
   image.png
   

   中间的数字显示出现的次数

   
   image.png
2. 简化图，并获得contigs
   去掉低频kmer链接（测序错误造成的）
   利用各种方法去掉其他的错误
   最终使 reads 连接成为 contigs
3. 将 pair-end mapping 到 contigs 上
   a. 将pair-end 定位到 contigs（完整的序列（长于Kmer），存在于reads和contig上）
   b. 测序方向的问题。
   c. 估计gap的大小，用 N 填充。
4. 构建 scaffolds
5. pair-end 中间 gap 的填充（延伸contig）
   a. 利用pair-end 的数据
   b. 长reads可以跨过洞

sanger/454测序

方法：overlap-layout-consensus

代表软件：Newber、Celera、Phrap

高杂合基因组

fosmid 文库 + 鸟枪法（牡蛎为例），花费高，周期长
BAC克隆（华大） + 鸟枪法，花费高，周期长
Ultra-Deep de novo，高深度测序 （下一节），花费少，效果差，周期短

9. 高深度测序原理及方法

利用Kmer分布图（峰），识别出 unique Kmer 和 repeat Kmer

image.png

10. 影响组装的因素

1. 测序错误

N太多（超过 10%，5%，3%，视数据多少而定），去掉

测序质量

接头序列，reads和接头比对

小的 insert size

去除 PCR duplication

测序错误引起的问题

Kmer 峰的前移

纠正
包含错误位点的Kmer出现的频率低，识别出，一次替换为其他三个碱基，直至Kmer的出现频率变高。

为了不人为改变造成错误，一般一个reads只改一个

先纠错 再进行 contigs 的组装

2. 重复序列

N50和N90 推测影响

3. 杂合

N50和N90 推测影响

如果杂合度超过千分之5，就要考虑是否考虑 高深度测序等。

11. 辅助组装的技术

1. Fosmid-end/ BAC-end
2. 遗传图谱和物理图谱
3. optical mapping，酶切，获得片段长度，辅助组装
4. The Irys System (BioNano)，密度高
   可以侦测 gap 的大小及错误的插入

尽量将组装结果，连接到假染色体上。

12. 评价基因组

基于长度

N50：contigs从长到短排序，累加到总长50及以上是，那个conig的长度。
N90：同理

genome 的覆盖度（80%以上）

基于准确度

与EST 或 转录组数据 比对
golden standard （比如 fosmid 等）
同源方法，与同源物种或模式生物比较基因模式

GC 深度和 测序深度的分析

下图，红圈部分，可能是y染色体

image.png

下图，红圈部分，左右两端的小圈区域，GC高的部分可能是细菌污染（可以过滤掉相关reads）

image.png image.png

gap 可能的产生原因：repeat 造成的深度高

基因颠倒，证明是否正确

基因结构预测，与近源物种比较

13 SOAPdenovo 组装流程

1. 下机数据准备

后有index，用于区分样品

质量值：Q = -10 log₁₀P，ASCII码，减64
PE（pair-end） reads：成对的reads（paired-end reads 和 mate-pair reads）
SE（single-end） reads

2. 过滤，数据质控

Q20，占全部碱基的80以上。

过滤 adaptor：fastq转化为fasta，利用nucmer将接头序列与reads相比较，过滤掉adaptor（SOAPfilter）

过滤低质量的reads：

N 占太多的 过滤掉（10%）

过滤掉 small insert 的 reads

两端的碱基准确率可能较低（两个reads，A,T不相等），可以剪掉

3. 下机数据纠错

KmerFraq_AR_ 软件
-K 设置成17，可以估计出小于16G的基因组
关注最后一个文件，标记各种kmer选项。

4. 基因组大小和杂合度的估计

需要软件
命令：./KmerFreq_AR_v2.0 -k 17 -t 4 -p species fq.lst

5. SOAPdenovo

配制文件，格式在官网上有介绍

avg_ins： 插入片段的长度

reverse_seq：是否环化，0 或 1

... ...

6. SOAPaligner