在生物信息学领域,一个非常有趣而重要的问题就是序列拼接(Sequence Assembly),而从头开始的序列拼接(de novo Assembly)尤为重要,因为很多时候我们并不能使用参考用的序列信息(比如当我们对全新的物种进行测序,或者我们必须不能依赖参考序列以获得新的具有生物学价值的信息的时候)。
为了确保所有的同学(有或没有生物学背景)都可以很好地理解本题,优酱首先先简单介绍一下目前较为通用的二代测序技术(NGS)。
二代测序(高通量测序)
我们都知道DNA分子是一个非常非常长的,由携带四种不同碱基(Adenine、Cytosine、Guanine、Thymine)的脱氧核苷酸组成的高分子长链。理想状态下,如果我们有一种魔法机器,可以令解旋的DNA链从一端(5'端),到另一端(3'端)从这个机器通过,同时读出碱基的类型,我们就可以知道DNA的序列了。然而很不幸的是,即使是当今可能最为先进的纳米孔(Nanopore)测序技术,在非常极端的情况下也只能最多允许1兆碱基的DNA片段通过,而人类基因组的总长高达3234.83兆碱基。而当今较为通用的第二代高通量测序技术,能通过的单个DNA片段仅仅有几百碱基。因此由于技术所限,我们目前做测序的一般流程都是:首先将DNA分子打碎成很短的碎片,然后加入引物、Barcode等人工合成的DNA片段用来辅助扩增和测序,然后进入PCR扩增到充足的DNA量,最后进入测序仪进行高通量测序。
链接:
在这样的工作流程下,注定我们手里的数据只能是若干较短的字符串,而我们知道整个基因组的长度如此巨大。我们要做的就是从较短的字符串中找到规律,将这些短串想办法拼接成可信度较高的长串,这便是序列拼接问题。(当然,如果可以使用参考序列的话,我们可以考虑将短串同参考序列进行比对(Alignment),也就是进行重测序(reSequencing)了,本题不讨论这种情况)
既然我们要对字符串的重叠关系建图,那么一个非常重要的课题就是:我们怎么样才能在海量的较短串中,以非常快的速度查询两个串的最长公共前后缀(因为我们在拼接串的时候一定是前缀和后缀拼接在一起的)?
我们将这个问题形式化如下:
如我们有串:"ACGT"、"AACCGGTT"、“AAAA”,则查询(2,3)的结果为0,因为"AAAA"在右面欲与左面的"AACCGGTT"进行拼接,显然并没有公共的前后缀区域可以使两串拼接。而查询(3,2)的结果为2,因为此时“AAAA”放在了左面欲与右面的"AACCGGTT"拼接,这时有最长的公共区域"AA"。
序列拼接问题的常见解决方案
目前针对二代测序的序列拼接主要有两种策略:一种是通过对字符串之间的重叠关系进行图论建模,得到字符串图,然后在字符串图上尝试找可信度最高的最长路径。另外一种则是通过将序列拆解为彼此重叠的"k-mer",然后对"k-mer"建de Bruijn Graph,在de Bruijn Graph上找欧拉回路或高可信度路径。本题主要关注的是第一种策略(如果对第二种策略感兴趣可以去看看Youhane Doorlock)。既然我们要对字符串的重叠关系建图,那么一个非常重要的课题就是:我们怎么样才能在海量的较短串中,以非常快的速度查询两个串的最长公共前后缀(因为我们在拼接串的时候一定是前缀和后缀拼接在一起的)?
我们将这个问题形式化如下:
最长公共前后缀查询问题
我们给出一个字符串的集合SS,在这个集合中有NN个字符串,编号为ii的字符串SiSi的长度为NiNi。我们给出QQ个查询,每个查询是一个整数对(l,r)(l,r),表示我们假定ll号串在左边,rr号串在右边,我们要查询的是最大的长度lenlen使rr号串的前缀pfx(r)lenpfxlen(r)与ll号串的后缀sfx(l)lensfxlen(l)相等。如我们有串:"ACGT"、"AACCGGTT"、“AAAA”,则查询(2,3)的结果为0,因为"AAAA"在右面欲与左面的"AACCGGTT"进行拼接,显然并没有公共的前后缀区域可以使两串拼接。而查询(3,2)的结果为2,因为此时“AAAA”放在了左面欲与右面的"AACCGGTT"拼接,这时有最长的公共区域"AA"。