《生物序列数据K-mer频次统计问题的算法①》由会员分享，可在线阅读，更多相关《生物序列数据K-mer频次统计问题的算法①（5页珍藏版）》请在金锄头文库上搜索。
1、2014 年 第 23 卷 第 4 期 http:/www.c-s- 计 算 机 系 统 应 用 Software TechniqueAlgorithm 软件技术算法 121生物序列数据K-mer 频次统计问题的算法 张鑫鑫1,2, 陈 波1,2, 何继凌1, 徐 云1,2 1(中国科学技术大学计算机科学与技术学院, 合肥 230027) 2(安徽省高性能计算重点实验室, 合肥 230027) 摘 要: 生物序列的 k-mer 频次统计是生物信息处理中一个非常基础且重要的问题. 本文针对多序列在对齐模式下, 不同偏移处一段长度范围内的 k-mer 频次统计问题进行了研究. 提出了一种逆向遍历
2、k-mer 计数算法 BTKC. 该算法能够充分利用长度的k-mer统计信息, 快速得到长度的k-mer统计信息, 从而避免了统计任意长度的k-mer频次信息时都需要对所有序列进行遍历. 算法的时间复杂度分析及实验结果表明, 相比于传统的前向遍历 FTKC算法, BTKC 算法性能提升非常明显, 且其时间复杂度与 k-mer 长度的变化范围无关, 非常适合于在 k-mer 长度变化范围较大的情况下使用. 关键词: k-mer; k-mer 计数; 频次统计; 逆向遍历; 生物信息处理 Algorithms for Biological Sequence K-mer Frequency Coun
3、ting Problem ZHANG Xin-Xin1,2, CHEN Bo1,2, HE Ji-Ling1, XU Yun1,2 1(School of Computer Science and Technology, University of Science and Technology of China, Hefei 230027, China) 2(Key Laboratory of High Performance Computing of Anhui Province, Hefei 230027, China) Abstract: K-mer counting of biolog
4、ical sequence is a fundamental and very important problem in biological information processing. This paper focuses on counting k-mers at each position of multiple sequences within aligned mode. We present a new backward traverse k-mer counting algorithm called BTKC. BTKC algorithm takes full advanta
5、ge of the k+1-mers statistic information to obtain k-mers statistic information quickly. Thus, its no need to traverse the whole sequences when counting each single k-mer. Both the algorithms time complexity and experiment results show that BTKC gets an obvious improvement compared with forward trav
6、erse k-mer counting algorithm FTKC, and its time complexity was found not to be realted with the range of k-mer length. Key words: k-mer; k-mer counting; frequency statistic; backward traverse; biological sequence processing 生物序列主要指RNA,DNA及蛋白质序列等. 自人类基因组计划实施以来, 特别是随着第二代测序技术大规模的使用与推广, 现今一次深度测序技术可以产生
7、TB 级的序列数据. 而在对这些海量数据的分析和应用中, k-mer 频次统计是一个非常基础且重要的问题. k-mer 频次统计信息可以用来揭示生物序列中各种子序列的分布规律, 它是一种衡量序列相似性的重要工具. 因而其在单体分型, 模体发现, 物种识别, 宏基因组分类, 序列拼接, 多序列比对, 变异探测及序 列纠错等众多的生物学问题上都有着重要且广泛的应用. 特别是在单体分型及模体发现等一些需要探究序列中块属性的问题上, 常常还需要对多序列中不同偏移处一段长度范围内的 k-mer 进行频次统计. 鉴于此, 本文针对多序列在对齐模式下, 不同偏移处一段长度范围内的 k-mer 频次统计问题进
8、行了研究, 并提出一种逆向遍历 k-mer 计数算法 BTKC, 相比于传统的前向遍历算法算法FTKC, 该算法显著降低了k-mer频次统计算法的时间复杂度. 基金项目:国家自然科学基金(60970085) 收稿时间:2013-08-29;收到修改稿时间:2013-09-26 计 算 机 系 统 应 用 http:/www.c-s- 2014 年 第 23 卷 第 4 期 122软件技术算法 Software TechniqueAlgorithm 1 相关工作 郝柏林院士等在2000年提出了一种使用2D分形图像的可视化方法来反映序列的 k-mer 的频次分布, 并在此基础上设计了基于B树的快速
9、k-mer频次统计算法1,2. 王树林等在2007年设计并实现了k-长DNA子序列内部计数算法和外部计数算法, 在该算法中, 利用一个巧妙构造的hash函数把k-长DNA子序列映射为整数, 从而把 k-长 DNA 子序列的频次统计问题转化为整数关键字的重复计数问题, 使得能够利用经典 B 树算法来解决 k-长 DNA 子序列的出现频次计数问题3. Marcais 等人在 2011 年提出了并实现了并行的 Jellyfish 算法对 k-mer 的频次进行统计, 该算法在多线程环境下利用了前缀数组和优化的无锁哈希表, 该算法不仅效率高并且是内存有效的4. Melsted等人则在 2011 年使用
10、 Bloom filter(布隆过滤器, 一种概率数据结构)来存储 k-mer 频次统计过程中出现的各种 k-mer 子串, 然后再进行二次扫描获取具体的k-mer 频次, 试验数据表明, 该算法相对于其他基于hash 的算法, 能够显著减小内存的使用5. 而在 2013年, Rizk 等人又提出了一种分块的 k-mer 频次统计算法 DSK. 该算法综合利用了 hash 映射和分块缓存, 确保在任意规模的数据集上, 该算法均可以在内存和磁盘空间受限情况下正常运行6. 但对于单体分型及模体发现等一系列需要对序列中不同偏移位置处一段长度范围内的k-mer进行频次统计从而探究序列的块属性的问题,
11、目前鲜有研究. 2 问题定义及描述 我们给出 k-mer 的定义如下: 设 s 为长为 m 的生物序列, 则12msq qq=?, 其中, ( 为生物序 列 的 符 号 空 间 , 比 如 对 于 基 因 序 列 , 则有 ). 一个长为 k 的子串是指序列 s 中从任意位置处开始的个连续符号, 称之为 k-mer. 对于多序列在对齐模式下, 不同偏移处一段长度范围内的k-mer频次统计问题, 我们给出如下定义: 记 为n条长为m的基因序列在对齐模式下构成的序列集合, 则 , 对于任意的1, i jn, 要求is和 的起始位置对齐. 对于给定的 k-mer 长度变化范围1k和212()kkk,
12、 分别对于各个 k-mer 长度12()k kkk在各个偏移位置(0)llmk 上, 计算由 n 条序列中每一条序列从偏移l处取长为 k 的k-mer 所构成的序列块中不同 k-mer 子串的出现频次. 图 1 给出了一个简单的示例, 所有序列的起始位置对齐, k-mer 长度的变化范围为 35, 首先计算偏移为0时, 长度为3的红色序列块中各个k-mer的出现次数, 然后计算长度为4的黑色序列块中各个k-mer的出现次数, 直到最大值5. 然后偏移位置向右滑动一位变为 2, 再依次统计各 k-mer 长度下各个子串的出现频次. 不断的重复上述过程, 直至到达序列的尾部. 图 1 k-mer
13、示例图 3 算法描述及复杂度分析 对于 2 中描述的问题, 目前的常规算法是前向遍历计数算法, 我们简称其为 FTKC. 在 FTKC 算法中, 对于每一个偏移位置l, 依次从k-mer长度变化范围的最小值1k遍历到最大值2k. 在每一次遍历时, 均会扫描所有 n 条序列, 取每条序列当前偏移处, 长度为当前 k-mer 长度的子串, 使用 hash 表来来统计和存储各子串的出现频次. 对于 FTKC 算法, 由算法的运行过程可知, 在由n条长度为 m 的序列构成的序列空间中, 对于任意k-mer长度12()k kkk共构成了1mk+个从左向右滑动的块. 对每一个块进行k-mer频次统计时,
14、均需要依次访问 n 条序列, 进行 n 次 hash 操作, 因而该算法的时间复杂度为21()O kk mn. 注意到上述 FTKC 算法每次在块中进行 k-mer 统计时, 均需要遍历所有 n 条序列, 而没有考虑充分利用已有的统计结果. 但我们仔细观察图 1 可以看出, 对于同一偏移处, 长度为 k 的块中每一条序列均是长度为1k +的块中每一条序列的长为 k 的前缀字符串. 比如, 对于偏移 0 处的第一条序列, ATT 为 ATTG 的前缀字符串. 更进一步考虑, 比如对于 80 条序列构成的序列空间, 如果我们已经得到了偏移为 0, k-mer 长度为 4 时的统计结果为ATTA:1
15、2, ATTG:8, ATTC:9, 1q , ,A T G C =12,ns ss =?js2014 年 第 23 卷 第 4 期 http:/www.c-s- 计 算 机 系 统 应 用 Software TechniqueAlgorithm 软件技术算法 123ATTT:10, ATCA:8, ATCG:7, ATCC:15, ATCT:11, 那么下一步我们仅仅只需要对上述得到的 k-mer 频次结果集合进行遍历, 对提取的每一条记录取其长为 3 的子串并使用 hash 表进行计数, 则只需要 8 次遍历操作, 就可以得到偏移为 0, k-mer 长度为 3 时的结果为 ATT:39, ATC:41, 然后再遍历 k-mer 长度为 3 时的结果集合, 得到 k-mer 长度为 2 时的结果为AT:40. 因此在我们提出的 BTKC 算法中, 对于每一偏移处, 首先遍历所有 n 条序列, 取每条序列中当前偏移处长度为2k(k-mer 长度最大值)的子串, 将该子串做为 key, 该子串的出现频次做为 value, 使用hash 表进行统计和存储, 得到结果为 . 下一步我们首先构建一个空的 hash 表 , 然后对上一步得到的
