于长永, 李淼淼, 赵楚, 马海涛
东北大学秦皇岛分校 计算机与通信工程学院, 河北 秦皇岛 066004
收稿日期：2019-03-13
基金项目：国家自然科学基金资助项目(61772124, 61401080, 61402087);河北省自然科学基金资助项目(F2015501049);河北省教育厅项目(QN2014339);中央高校基本科研业务费专项资金资助项目(N150402002)。
作者简介：于长永(1981-), 男, 辽宁海城人, 东北大学副教授。

摘要：针对相似性确认步骤中编辑距离计算的高复杂性问题, 提出了一种在编辑距离限制下的基于鸽笼原理的字符串相似性确认算法.首先找到满足编辑距离片段映射的片段, 以此片段为基准, 将长度为500 bp的read分段.然后对满足编辑距离片段映射的左右部分递归地进行编辑距离计算, 将各段得到的编辑距离相加即为最后结果.最后根据最长公共子串的下限将需要验证的片段数目降到最低, 得到优化方案.实验结果表明, 基于鸽笼原理的分段递归计算编辑距离的确认算法减少了验证步骤的时间, 并能保证假阳率和假阴率都为零.
关键词：读取映射编辑距离相似性查询鸽笼原理确认算法
A Novel Similarity Verification Algorithm Under Edit Distance Limitation
YU Chang-yong, LI Miao-miao, ZHAO Chu, MA Hai-tao
School of Computer & Communication Engineering, Northeastern University at Qinhuangdao, Qinhuangdao 066004, China
Corresponding author: LI Miao-miao, E-mail: lmm_neu@126.com
Abstract: For the high complexity problem of edit distance calculation in the similarity confirmation step, a string similarity confirmation algorithm based on pigonhole principle under the limitation of edit distance is proposed. Firstly, it found a segment satisfies the edit distance segment mapping, then based on this segment, the read length of 500 bp was segmented. Then the edit distance of right and left parts was calculated and all parts sum to the final result. Finally, based on the longest common substring, the number of segments to be verified is minimized, and an optimization scheme was obtained. The experimental result show that this verification algorithm reduces the time of the verification step, and ensures the zero false positive and false negative rate.
Key words: read mappingedit distancesimilarity searchpigonhole principleverification algorithm
字符串在计算机系统中无处不在, 因此字符串处理得到了广泛的研究.字符串处理中最重要的问题是如何有效地评估两个字符串间的相似性^[1].相似性查询, 即给定一个字符串集合和一个查询字符串, 在这个集合中找到所有与查询串满足相似性度量的字符串^[2].相似性查询的应用场景非常广泛, 在基因序列比对、拼写检查、复制检测等方面都有涉及^[3-6].因此, 字符串的相似性查询成为数据领域的一个关注点.
现有的大多数相似性查询都要求具有一个特定的相似性度量函数, 例如, Jaccard相似性、余弦相似性和编辑距离^[7-11].本文只依据两个字符串间的编辑距离来判断它们是否相似, 但是, 当数据集很大时, 穷举所有字符串对的编辑距离是不现实的.为了解决这一问题, 使用过滤加验证的框架进行相似性查询^[12].在验证阶段, 对候选集合进行编辑距离验证, 得到满足阈值条件的结果集^[13].很显然, 高效的编辑距离计算方法对算法的时间性能有很大的提升.
本文研究基于鸽笼原理的分段递归计算编辑距离的字符串相似性查询方法.为了证明本文所提算法的有效性, 与现有的方法进行了比较.实验结果表明, 该方法能够在时间上达到更优的效果.
1 问题定义及相关工作定义1??(编辑距离限制下的相似性查询问题):设S={s₁, s₂, …, s_n}是一个字符串集合, q表示一个查询字符串, τ代表编辑距离阈值, 式(1)表示在编辑距离限制下的相似性查找结果.

(1)

其中, ed(s_i, q)表示将q转换为s_i所需的最少操作次数, 包括插入、删除、替换操作^[14].字符串相似性查询问题是字符串研究方向的基础问题.最原始的查询算法是将查询串和数据库中每一个字符串两两计算编辑距离, 很显然, 这是非常耗时的.为了减少相似性查询所用的时间, 采用过滤加验证的框架.大量实验数据表明, 验证占据了大部分的时间.目前可以用来优化验证的方法主要分为两类:①减少候选集合的规模; ②提高编辑距离计算速度.对于方法①, 由于候选集合的规模大小会有一个下限, 当达到一定规模时, 再减小候选集合规模, 会导致过滤成本大大增加.对于方法②, 大部分算法是基于动态规划进行优化.例如, 基于范围的计算^[1]、分层计算、前缀树^{[1, 15-16]}等方法.但是, 对于大规模数据而言, 如何进行验证步骤的优化、实现快速地相似性查询仍然具有挑战性.
1.1 动态规划方法设s和r是两个字符串, 利用|r|+1行|s|+1列的矩阵D计算ed(s, r).其中, |r|和|s|分别表示字符串r和字符串s的长度.s[i, j]表示s中从i到j的子串.ed(s, r)表示将r转换到s所需要的最少操作次数(插入、删除、替换).D[i, j]表示前缀r[1, h]和前缀s[1, j]的编辑距离, 很显然, 对于0≤i≤|r|, D[i，0]=i成立; 对于0≤j≤|s|, D[0，j]=j成立.按照式(2)迭代计算D[i, j]:

(2)

当s[j]=r[i]时, &=0, 否则&=1.该方法时间复杂度为O(|r|×|s|).给定数据集字符串s=series, 查询串r=seraji, 图 1给出了利用原始方法计算编辑距离的距离矩阵.
图 1(Fig. 1)

图 1 动态规划计算编辑距离Fig.1 Computing edit distance by dynamic programming

1.2 范围给定字符串s和r, |s|=m, |r|=n, 编辑距离阈值为τ, 那么只需要计算矩阵中的K条对角线即可, K的值可以由式(3)确定.同时, 当主对角线上的值出现τ-|m-n|+1时, 即可判断字符串s和r不相似.

(3)

例如, s=series, r=seraji, 如果给定阈值τ=2, 那么(τ-(m-n))%2=0, 最多需要计算τ+1条对角线.
1.3 层级令<i, j>表示距离矩阵中第i行第j列的实体, E_x表示矩阵中值为x的实体集合.给定字符串s和r, 编辑距离阈值τ, 为了计算r和s的编辑距离, 需要递归的计算E_x, x=0, 1, …, n.如果<|r|, |s|>∈E_x, 那么ed(s, r)=x.当ed(s, r)≤τ时, 加入结果集中.例如, s=series, r=seraji, 逐步计算编辑距离的过程如下:
首先，E₀={<1, 1>, <2, 2>, <3, 3>}, E₁={<2, 1>, <1, 2>, <3, 2>, <2, 3><4, 3>, <3, 4>, <4, 4>}, <|r|, |s|>不属于E₁, 根据E₁计算E₂, E₂={<3, 1>, <1, 3>, <4, 2>, <5, 3>, <6, 4>, <2, 4>, <3, 5>, <4, 5>, <5, 5>, <6, 4>}, <|r|, |s|>不属于E₂, 根据E₂计算E₃, E₃={<6, 5>, <5, 6>, <6, 6>}, <|r|, |s|>∈E₃, 因此ed(s, r)=3.计算过程见图 2, 时间复杂度为O(τ×min(|r|, |s|)).
图 2(Fig. 2)

图 2 分层计算编辑距离Fig.2 Leveled edit distance calculation

这种方法需要枚举数据集中每个字符串, 计算复杂度也很大.除此之外, 时间复杂度和min(|r|, |s|)以及字符串的长度有关.很显然, 当字符串很长时, 时间性能会有所下降.由于暴力枚举数据集中的每个字符串存在困难, 可以利用不同字符串间的共同前缀共享计算.
1.4 前缀树由于数据库中的字符串类似于基因序列、Web网址等具有高度重复性, 因此, 对集合中的字符串动态构建一个前缀树结构, 有利于解决具有相同前缀的字符串重复性问题.
给定字符串数据集S和查询字符串q, S={s₁, s₂, …, s_n}, 对S中的字符串构建前缀树, 相同字符串具有相同的前缀, 在前缀树上表示为具有相同的节点.例如, s₁=AGGAA, s₂=AGGAC, q=AGGCT, 那么s₁和s₂共享4个前缀树节点.因此, 可以直接利用s₁的子串AGGA与q=AGGCT的编辑距离快速得到s₂和q的编辑距离.
但是, 对候选集合中存在的较长字符串而言, 前缀树结构也不能很好地解决复杂的编辑距离计算问题.
通过对上面三种编辑距离计算方法进行分析, 相对于最初的利用动态规划来计算编辑距离, 虽然在时间上有所提升, 但是, 当ed(s, r)>τ, 或者max (|s|, |r|)很大时, 以上方法并不能达到很好的效果.为了解决该问题, 本文提出一种在片段映射上递归地利用鸽笼原理计算编辑距离的算法, 并对该算法提出了优化方案.
2 Recursive-ed算法2.1 编辑距离片段映射定义2??(编辑距离片段映射):给定两个字符串r和s.P(r)={seg^r(i)}和P(s)={seg^s(i)}, i=1, …, n, 分别表示s和r的分段集合.如果编辑距离ed(r, s)满足式(4), 那么称f(seg^r(i))=seg^s(i)为编辑距离分段映射(ed-mapping).设ed_l，ed_r表示seg^s(i)与seg^r(i)之间左、右两部分的编辑距离, 如果满足式(5), 则称seg^r(i)与seg^s(i)之间是ed-mapping, 记为f^ed(*).

(4)

(5)

例如, 给定数据集中的一个字符串s=abcd, 查询串r=acxd, 根据式(4), 可以列出从字符串s到字符串r的编辑距离片段映射，ε表示空字符.
1) f(a)=a, f(b)=ε, f(c)=c, f(ε)=x, f(d)=d,
2)f(ab)=a, f(c)=c, f(d)=xd,
3)f(ab)=a, f(c)=cx, f(d)=d,
4) f(abc)=ac, f(d)=xd.
根据定义2, 如果知道字符串r和s之间的1个编辑距离片段映射, 那么2个字符串间复杂的编辑距离的计算可以转换为多个子字符串编辑距离的求和计算, 减少长字符串编辑距离计算的复杂度.但是使用定义2去确定两个字符串中的分段是否满足ed-mapping是一个复杂的过程, 下面通过定理1和定理3来介绍如何巧妙地运用ed-mapping解决加速计算编辑距离的问题.
定理1??给定两个字符串s和r, 将r分成n段, 如果s与r在给定的编辑阈值τ下相似, 那么至少存在一个seg^r(i)(1≤i≤n)和一个s中的子串s_m, 满足:
1) ；
2) f^ed(seg^r(i))=s_m.
证明??由于ed(s, r)≤τ, 即对r进行不多于τ次的编辑操作就可以将r变为s, 因此, 将r分为n段后, 至少有一段在变换的过程中发生的编辑操作次数小于, 设该段为seg^r(i), 同时设该段在变换中变为s的子串s_m, 则有seg^r(i)和s_m满足条件1和条件2, 证毕.
在定理1中, 给出了快速判断两个字符串分段中是否存在ed-mapping, 同时以下推论也成立:
1) 如果τ<n, 那么seg^r(i)=s_m.
2) seg^r(i)在字符串r中的位置与s_m在s中的位置差小于等于τ.
2.2 Recursive-ed算法定理2??给定两个字符串s和r, 将r分成n段, 如果s与r在给定的编辑阈值τ下相似, 那么可以推断出定理1中的结论.并且以下结论也成立:
3) 如果找到编辑距离映射ed(seg_i^r, s_m), 那么左半部分的r_l和s_l必须满足定理1, 其中左半部分的编辑距离阈值为τ_l=τ-(|r_r|-|s_r|+ed(seg_i^r, s_m)).
4) 右半部分同理, r_r和s_r必须满足定理1, 此时的编辑距离阈值为τ_r=τ-(|r_l|-|s_l|+ ed(seg_i^r, s_m)).
其中，seg^r(i)和s_m满足f^ed(seg^r(i))=s_m, r_l, r_r和s_l, s_r分别表示seg^r(i)和s_m的左、右两部分, τ_l, τ_r分别表示r_l和s_l, r_r和s_r之间的编辑距离阈值.通过上式可递归的计算编辑距离.
定理3??如果两个字符串r和s满足ed(s, r)≤τ, r的分段集合为P(r)={seg^r(i)}, 那么有.其中，f(seg^r(i))=s_mⁱ表示seg^r(i)和s_mⁱ在编辑距离限制下是分段映射.值得注意的是, 这样的映射是否为编辑距离映射是不必要的.只要满足片段seg^r(i)能够映射到s_mⁱ, 那么就可以用ed_i(s, r)=ed_l+ed_r+ed(seg_i^r, s_m)来表示总的编辑距离.根据定理1, 定理2和定理3显然成立.
2.3 Recursive-ed算法的优化定义3??(最小覆盖字符串):给定字符串s和r, 如果c(s, r)满足:
1) ed(s, r)=ed(s, c)+ed(c, r);
2) s→c过程中只有插入和替换操作, c→r过程中只有删除和替换操作.
那么称c(s, r)为s和r的最小覆盖字符串.很显然, c(s, r)一定存在并且不唯一, 并且|c(s, r)|>|s|, |c(s, r)|>|r|.|c(s, r)|可以通过以下定理得到:
定理4??给定字符串s和r, 令F是从s到r编辑距离映射的集合.那么,

(6)

(7)

由于|c|不能直接计算出来, 因此用|c|*=max {|s|, |r|}代替|c|, 并随着假设的映射片段的改变, |c|*逐步改变.令作为|c|的近似值.
证明??c(s, r)表示r变为s的过程中, 先对s进行所有的插入操作, 得到字符串c, c的长度为|c(s, r)|.因此, 将s分段后, 每段变为r的一段, 这变换中每个对应段的表示每段都只进行插入, 得到字符串c_i, c_i的长度和就是.因此整体变换和分段变换的结果是相同的, 所以式(6)成立.
取s和r的分段映射, 每段(s_i, r_i)至多包含一个插入或一个删除, 那么有c(s_i, r_i)=max{|s_i|, |r_i|}, 所以式(7)成立.
引理1??(共享q-gram的下限):给定字符串s和r, 满足ed(s, r)≤τ.那么s和r至少共享N_q-gram=|c(s, r)|-q+1-τ×q个公共q-grams.当s→r只有删除或插入操作时, LB=max{|s|, |r|}-q+1-q×τ是共享q-grams的最少数量.
引理2??(最长公共子串的下限):给定字符串s和r, 满足ed(s, r)≤τ.那么必然存在长度至少为.并且, 当共享q-gram的数量大于N_q-gram时, L^*=.
在优化算法中, 通过减少被验证片段的数目, 达到优化目的.由引理1，引理2可知, 共享q-gram的下限和最长公共子串的下限都与|c(s, r)|有关, 因此, 当一个新的片段映射被假设之后, |c|*, N_q-gram, L^*的值都被更新.
由鸽笼原理可知, 字符串r被分为τ+1个片段, 每一个片段都有可能ed-mapping到s的子串, 因此, 为了验证ed(s, r)≤τ, 最坏的情况下需要对τ+1种情况进行验证.实际上, 根据引理2, 只有较长片段才会被验证.
定理5??给定字符串s和r, 如果ed(s, r)≤τ, 那么必定存在s的子串sub^s和r的子串sub^r满足:
1) sub^s=sub^r;
2) F_ed(sub^r)=sub^s;
3) |sub^r|≥L^*(L₂).
根据定理5, 只需要考虑映射片段长度大于L₂的情况.为了找到长度大于L₂的所有映射片段, 需要将字符串r分成个相同长度的片段.同时, 利用q-gram的计数过滤, q=L₂.
3 实验结果与讨论3.1 数据集和实验环境本文实验在Linux操作系统下, g++编程环境上利用C++语言实现, 本文实验使用真实的测序read数据集, 具体见下表 1, 该数据集来自于文献[14], 从该文献的实验数据集中选择两个数据作为本文的实验数据, read长度为100 bp, 在实验过程中, 将长度为100 bp的read扩展到500 bp.
表 1(Table 1)

表 1 数据集Table 1 Database 数据集 大小/Mb 长度/bp ERR240726_1 984 100 ERR240726_2 880 100

表 1 数据集 Table 1 Database

3.2 实验结果根据定理1~4, 本文提出的算法的假阳性和假阴率均为0.当两个字符串s和r的编辑距离小于τ时, Recursive-ed返回真实的编辑距离, 否则Recursive-ed返回τ+1.因此, 本文提出的算法的准确性是100%.
对于上述两个数据集, 分别随机产生编辑距离τ≤x(x=3, 5, 7)的10 000条和30 000条基因序列, 序列数目可以根据实验要求进行更改.
表 2(Table 2)

表 2 验证候选集Table 2 Candidate database 候选集 大小/Mb 基因序列/条 Pairdata_1_3 784.9 10 000 Pairdata_1_5 783.5 10 000 Pairdata_1_7 782.3 10 000 Pairdata_2_3 784.9 30 000 Pairdata_2_5 783.5 30 000 Pairdata_2_7 782.3 30 000

表 2 验证候选集 Table 2 Candidate database

在实验中, 为了验证本文所提算法的有效性, 将序列的长度定义为500 bp.对每个候选集利用本文列举的几种不同编辑距离计算方法, 实验结果见表 3.
表 3(Table 3)

表 3 Pairdata_1_x各函数所用时间Table 3 Time spent on Pairdata_1_x functions s τ 原始 层级 范围 Recursive-ed 3 31.37 2.37 2.59 0.14 5 31.84 2.43 2.68 0.22 7 32.11 2.51 2.74 0.32

表 3 Pairdata_1_x各函数所用时间 Table 3 Time spent on Pairdata_1_x functions

为了测试候选序列数目对算法的影响, 对验证集Pairdata_2_x选取Number=30 000的情况, 实验结果见表 4.
表 4(Table 4)

表 4 Pairdata_2_x各函数所用时间Table 4 Time spent on Pairdata_2_x functions s τ 原始 层级 范围 Recursive-ed 3 95.96 7.28 7.88 0.40 5 96.73 7.59 8.14 0.69 7 95.79 7.66 7.94 0.43

表 4 Pairdata_2_x各函数所用时间 Table 4 Time spent on Pairdata_2_x functions

4 结语本文提出的Recursive-ed算法用来快速计算编辑距离.在Recursive-ed算法中, 利用鸽笼原理, 通过将序列分段递归计算编辑距离, 有效地完成相似性验证步骤.和已有的其他4种算法相比较, 本文提出的方法更加适用于较长字符串的相似性比较, 能够克服因max(|s|, |r|)很大而造成的计算复杂的问题.现如今序列数据呈指数模式增长, 使得序列间的相似性查询越来越困难.本文提出的算法可以作为该问题的潜在解决方案.
参考文献

[1]	Deng D, Li G L, Feng J H, et al.Top-k string similarity search with edit-distance constraints[C]// IEEE 29th International Conference on Data Engineering(ICDE).Washington D C, 2013: 925-936.
[2]	Zhang Z J, Hadjieleftheriou M, Ooi B C, et al.Bed-tree; an all-purpose index structure for string similarity search based on edit distance[C]// Proceedings of the ACM SIGMOD International Conference on Management of Data.Indianapolis, 2010: 915-926.
[3]	Deng D, Li G L, Feng J H.A pivotal prefix based filtering algorithm for string similarity search[C]// SIGMOD Conference.Snowbird, 2014: 673-684.
[4]	Hadjieleftheriou M, Koudas N, Srivastava D.Incremental maintenance of length normalized indexes for approximate string matching[C]// SIGMOD Conference.Rhode Island, 2009: 429-440.
[5]	Li C, Lu J H, Lu Y M, et al. Efficient merging and filtering algorithms for approximate string searches[J].International Conference on Data Engineering, 2008, 7(12): 257–266.
[6]	Wandelt S, Wang J, Leser U, et al. State-of-the-art in string similarity search and join[J].ACM Sigmod Record, 2014, 43(1): 64–76.DOI:10.1145/2627692.2627706
[7]	Jiang Y, Li G, Feng J, et al. String similarity joins:an experimental evaluation[J].International Conference on Very Large Data Bases, 2014, 7(8): 625–636.
[8]	Xiao C, Wang W, Lin X M, et al. Efficient similarity joins for near-duplicate detection[J].ACM Transaction, 2011, 36(3): 131–140.
[9]	Jiang Y, Deng D, Wang J N, et al.Efficient parallel partition-based algorithms for similarity search and join with edit distance constraints[C]//EDBT/ICDT.New York, 2013: 341-348.
[10]	Bi F, Chang L, Zhang W, et al.Efficient string similarity search: a cross pivotal based approach[C]// International Conference on Database Systems for Advanced Applications.Hanoi, 2015: 545-564.
[11]	Levenstein V. Binary codes capable of correcting spurious insertions and deletions of ones[J].Problems of Information Transmission, 1965, 1(1): 8–17.
[12]	Yu M, Wang J, Li G L, et al. A unified framework for string similarity search with edit-distance constraint[J].International Conference on Very Large Data Bases Journal, 2017, 26(2): 249–274.DOI:10.1007/s00778-016-0449-y
[13]	米琳. 基于q-gram的字符串相似性查询研究[J].现代计算机(专业版), 2014(6): 12–16. ( Mi Lin. Research on string similarity query based on q-gram[J].Modern Computer(Professional Edition), 2014(6): 12–16.)
[14]	Mishra S, Gandhi T, Arora A, et al.Efficient edit distance based string similarity search using deletion neighborhoods[C]//EDBT/ICDT.New York, 2013: 375-383.
[15]	Rheinl?nder A, Knobloch M, HochmuthN, et al.Prefix tree indexing for similarity search and similarity joins on genomic data[C]//International Conference on Scientific and Statistical Database Management.Berlin, 2010: 519-536.
[16]	Xin H, Greth J, Emmons J, et al. Shifted hamming distance:a fast and accurate SIMD-friendly filter to accelerate alignment verification in read mapping[J].Bioinformatics, 2015, 31(10): 1553–1560.DOI:10.1093/bioinformatics/btu856