算法设计Project:宏基因组组装

背景

DNA序列

DNA序列为字符集为Σ={A, C, G, T, N} 的字符串。

其中N表示未知,可能为A, C, G, T中任一个。(也可以理解为通配符)

测序技术

常规测序

对于普通的测序,我们先通过PCR技术将DNA复制若干份。再将它们随机切割为长度接近的小段,分别进行测序。
这样我们将得到大量互相交叠的小段,可通过重叠部分将其拼接起来得到完整的DNA序列。

简而言之,测序后得到的数据将是大量长度接近,均匀分布的原字符串的子串。

序列的方向

方向 说明 示例
原始方向 一般任意选择一条链为序列的原始方向 ACAGT
相反方向 与原始方向方向相反 TGACA
互补方向 DNA中另一条链上的序列方向,即原串按照A-T,C-G转换后得到的序列 TGTCA
互补相反方向 既互补又相反 ACTGT

Illumina Paired end短序列

在测序时,若复制次数过少即子串数量过少,则交叠部分不足,可能无法拼接出完整的序列,而增加复制次数则会成倍增加成本。

因此提出了paired end技术。即先将原串切割为长度接近的较长片段,再从长段两侧分别测序之前小段的长度。
这样我们将会得知某一对字串之间的大致距离,方便拼接出完整的序列。
由于paired-end测序时从两条链的两端向中间进行测序,因此会得到两条有间隔的互为互补反向的短序列。

Pacbio Single end长序列

Pacbio使用了实时单分子DNA测序。
将脱氧核苷酸用荧光标记,使用显微镜实时记录荧光的强度变化来对DNA进行测序。

相比于Illumina Paired end技术得到的短序列,该方法得到的序列长度有显著增加,使得能从中获得更多关于长度关系的序列信息。
但是受限于测序原理,该方法得到的序列测序错误率远高于Illumina方法得到的短序列。

测序错误

由于技术原因,DNA测序得到的子串会有一定几率错误(对于Illumina测序方法几率小于每个字符1%,对于Pacbio测序方法,几率约为15%~20%)
测序错误体现为

  1. 替换:某一位置字符替换为其他字符。

  2. 插入:两字符间插入一字符。

  3. 删除:某字符未被测出而缺失。

宏基因组测序

宏基因组学通过通过直接从环境样品中提取全部微生物的DNA,构建宏基因组文库,利用基因组学的研究策略研究环境样品所包含的全部微生物的遗传组成及其群落功能。
宏基因组测序即对环境样品的所有微生物进行测序,测序得到的结果为所有微生物测序序列的混合。因此要分析这类数据较为困难。

基因组组装

基因组组装的任务是将测序所得的序列片段按照某种算法进行拼接,从而得到更长的片段(contig)用于后续分析。
最理想的拼接结果为还原完整的目的参考基因组序列。

拼接的主要依据是参考不同序列中连续出现的部分(overlap),下面展示了将两条序列拼接成一条长序列的操作(实际设计的算法不一定要采用这样的方法)。

问题描述

在本次课程Project中,你将设计一个宏基因组拼接的算法,输入数据为包含多个物种的测序数据(包含pacbio长序列和illumina短序列),你可以使用所有数据进行组装,也可以仅使用其中的部分。输出结果为拼接后得到的长片段。
衡量拼接结果好坏的依据主要为长度和准确率。

输入与输出

输入

1. 短序列:"short_1.fasta" "short_2.fasta"

文件中每两行代表一个序列。

相同行数的序列成对,互为反向互补方向(paired-end),但short_1.fasta中不一定全是原始方向。

fasta格式每两行代表一个序列,其中第一行为">"加上序列名称,第二行为对应序列。
例:example.fasta:

>Read_name
AAATTCCGGG

2. 长序列:"long.fasta"

文件中每两行代表一个序列,序列方向可能为原始方向或者反向互补方向。

3. "param.json"

包含测序参数的json文件

其中包含:
  • name:测试数据名称

  • short_read_length:短序列

  • pair_distance:paired-end短序列间平均距离,这里的距离包含两条read的长度以及它们之间的间隔

  • standard_deviation:paired-end短序列间距离的标准差

  • short_read_error_rate:短序列错误率

  • long_read_length:长序列长度

  • long_read_error_rate:长序列错误率

其中非整数的参数均用小数表示。

输出

1. 拼接结果:"contig.fasta"

每两行输出一个拼接好的片段

评测方式

评测采用提交答案的形式,对提交的组装结果使用软件metaquast(Mikheenko A, Saveliev V, Gurevich A. MetaQUAST: evaluation of metagenome assemblies[J]. Bioinformatics, 2015, 32(7): 1088-1090.) 其中评测的主要参考指标为:

  • Genome fraction(%):组装好的contig比对到参考参考基因组后占参考基因组总长度的比例。

  • NGA50: 将组装好的contig比对按长度从大到小匹配到到参考基因组,当匹配总占比超过原基因组50%时该contig的长度
    对不同物种分别计算后取平均值。

  • Misassemblies:组装好的contig中连续的片段映射到参考基因组中非连续的片段上,则判定该contig为misassemlby(组装错误),最后统计所有misassembly的数目。

  • Mismataches per 100kbp: 匹配回参考基因组后每100000个字符中与参考基因组不同的字符数(可以理解为拼接结果的错误率,参考测序的错误率描述)

测试数据下载

数据名称 数据描述 输入文件
Toydata1 参考基因组为五个物种,每个物种genome长度为10000。短序列平均复制次数34x,长序列平均复制次数5x。 下载

联系我们

如果你对题目有疑问,欢迎联系:

黄剑铮 14307130264@fudan.edu.cn