跳到主要内容

uShuffle:一个有用的工具,可以在保留k-let计数的同时对生物序列进行洗牌

摘要

背景

随机洗牌序列通常用于序列分析,以评估生物序列的统计意义。在许多情况下,生物学家需要复杂的洗牌工具,这些工具不仅可以保存不同字母的计数,还可以保存高阶统计信息,如双重计数、三重计数,通常,k个-让我们数数。

结果

我们提出了一种序列分析工具(称为uShuffle),用于生成生物序列(如DNA、RNA和蛋白质)的均匀随机排列,以保持精确的k个-让我们数数。uShuffle工具实现了Euler算法的最新变体,并在树木生成的关键步骤中使用Wilson算法。它经过精心设计,效率极高。uShuffle工具通过允许任意字母大小和let大小来实现最大的灵活性。它可以用作命令行程序、web应用程序或实用程序库。提供了C、Java和C#的源代码,以及Perl和Python的集成说明。

结论

uShuffle工具在性能和灵活性方面都优于Euler算法的现有实现。它是生物信息学社区的有用工具。

背景

随机洗牌序列通常用于序列分析,以评估生物序列的统计意义。例如,评估RNA序列热力学稳定性的常用方法是将其折叠自由能与随机序列的大样本的折叠自由能进行比较。众所周知,RNA二级结构的稳定性关键取决于相邻碱基对的堆叠;因此,在这种分析中,随机序列中不同偶极子的频率是重要的考虑因素[4,25]. 此外,自然生物序列往往表现出某些近邻模式:真核和原核核酸序列在双频上显示出一致的层次结构;在编码区,密码子的使用也可能明显不均匀。在许多情况下,生物学家需要复杂的洗牌工具,这些工具不仅可以保存不同字母的计数,还可以保存高阶统计信息,如双重计数、三重计数,通常,k个-让我们数数。

随机序列生成方法

生成随机序列通常使用几种方法。这个基本置换法工作原理如下:对于序列S公司[1,n个],选择一个随机数介于1和之间n个,交换这两个元素S公司[]和S公司[n个],然后在子序列上递归S公司[1,n个- 1]. 基本置换方法生成的随机序列保留了字母表中每个不同字母的精确计数,但不保留k个-让我们数数。这个马尔可夫方法[12],它基于马尔可夫链生成随机序列k个-仅平均计数:单个序列的计数可能偏离输入分布。这个交换方法[15]这是一种现在流行的方法,通过反复交换两侧不相交的子序列来生成随机序列(k个-1)-let;它确实保存了k个-让我们精确地计数,但产生的随机序列只是一致的渐近序列,可能需要大量的交换步骤。

欧拉算法保持精确k个-让计数

欧拉算法是一种不太为人所知但非常有效的算法,用于生成真正一致的随机数k个-let-preserving序列[2,12,15]. 我们简要回顾一下它的历史。惠誉(Fitch)[12]首先注意到保偶置换与有向多重图的欧拉游动有关;然而,他提出的算法并不是以等概率生成所有排列。Altschul和Erickson[2]提出了第一种算法(也是基于有向多重图中的欧拉游动),用于生成真正一致的随机序列,该序列既保留二重计数,也保留三重计数,或同时保留二者;然而,他们生成随机树状图的算法的关键步骤依赖于一个试错过程,这是一个潜在的性能瓶颈。坎德尔等人消除了这一瓶颈[15]他用一个基于有向多重图中随机游走的简单而有效的程序取代了试验和错误程序。他们还推广了欧拉算法,以保持k个-让任意计数k个,并建议了一个简单的数据结构来实现。此数据结构基于查找表,需要O(运行)(σ2k个-2个)空间与时间;随着字母表的大小,它很快就会变得效率低下σ和出租大小k个增加。自从坎德尔等人的工作以来[15],Wilson提出了一种更好的算法[19,23]生成随机树状图,这是Kandel等人提出的Euler算法的关键步骤[15]改进了Altschul和Erickson[2]. Wilson的树状图生成算法相对于Altschul和Erickson的前两种算法的优越性[2]Kandel等人[15]两者都是由威尔逊在理论意义上证明的[19,23],并通过将我们的实现与以前的实现进行比较(稍后讨论),从实际意义上进行了演示。

欧拉算法的实现

我们知道欧拉算法早期变体的两个先前实现。这个洗碗Clote等人的项目[6]由Altschul和Erickson实现Euler算法的原始版本[2]. 这个洗牌Coward的节目[11]实现了Kandel等人改进的Euler算法[15]. 在这篇文章中,我们提出了一个序列分析工具(称为uShuffle),用于在保留k个-让我们数数。uShuffle程序基于欧拉算法的最新变体[2,15]并使用Wilson的算法[19,23]在乔木生成的关键步骤。我们的目标是提供尽可能高效和灵活的通用工具:

任意字母大小和let大小

在特定应用中,字母大小σ和出租大小k个通常是固定的:对于生物序列,典型的字母大小是4(对于DNA或RNA)和20(对于蛋白质),典型的let大小是2(对于二核苷酸)和3(对于密码子)。虽然只针对固定的字母表并让大小在手边,实现Euler算法是很诱人的,但我们相信任意字母表和let大小的灵活性是有用的。Clote等人的洗碗节目[6]例如,是硬编码的,用于改变RNA序列,保留二核苷酸计数(字母大小σ=4,让大小k个= 2). 很明显,这种实现不能轻易地用于具有不同字母表和let大小的其他应用程序。

效率

当字母表大小和let大小都是小常数时,Euler算法的运行时间(树状图生成的三种变体中的任何一种[2,15,23])在序列长度上是线性的。因此,洗牌程序的效率似乎不会成为问题,因为对随机数据的任何可想象的下游分析都会比洗牌慢得多。然而,我们注意到,线性运行时间仅在字母表和let大小不变的情况下得到了证明[15]. 对于任意字母表和let大小,欧拉算法的线性性能是否可扩展,目前尚不清楚。如前所述,Kandel等人建议的“标准”数据结构[15]具有时间和空间复杂性O(运行)(σ2k个-2个),当字母表大小为σ和出租大小k个变大,接近序列长度的顺序。事实上,正如我们稍后将要讨论的那样,我们有理由相信Coward在shufflet程序中使用了这种数据结构[11].

此外,欧拉算法的实现(特别是树图生成的关键步骤)是不平凡的,因为它大量使用了图论概念,如有向多重图和欧拉行走。尽管威尔逊的著名算法[19,23]可追溯到1996年,并在理论计算机科学界广为人知,科沃德于1999年实现了shufflet[11]仍然使用Kandel等人的旧树形算法[15]. 我们不知道Wilson算法在生物信息学应用中有任何实现。通过仔细选择算法和数据结构,并通过严谨的算法工程,我们努力实现最高效的实现。

多种形式和编程语言

Clote等人的洗碗节目[6]用Python编写;科沃德的洗牌节目[11]是一个C语言的web应用程序。为了达到最广泛的受众,我们以多种形式提供uShuffle程序。它可以用作命令行程序、web应用程序或实用程序库。提供了C、Java和C#的源代码,以及Perl和Python的集成说明。

实施

本节由四个小节组成。在前两小节中,我们在概念层面上讨论了Euler算法及其树木生成的关键步骤,为实现细节的讨论做准备。在第三小节中,我们介绍了实现的算法工程细节。在第四小节中,我们描述了uShuffle工具的软件组织和用户界面。为了证明我们的算法选择并解释我们的优化技术,前三小节中的讨论必然是技术性的。对图形算法的理论讨论或算法工程的技术细节不感兴趣的读者可以安全地跳到软件组织和用户界面的第四小节。

欧拉算法

在本小节中,我们将回顾欧拉算法的一些基本概念。

有向多重图

A类k个-下面是k个序列中的连续元素。S公司是要排列的序列。T型 k个 是一个保持k个-让计数S公司(例如,T型1是的简单排列S公司,以及T型2是的排列S公司具有相同的二核苷酸计数。)要生成T型 k个 对于k个≥2,欧拉算法[2,15]首先构造一个有向多重图G公司我们参考图1例如。对于每个不同的(k个-1)-引入S公司,G公司有一个顶点。对于每个k个-出租L(左)在里面S公司,其中包含两个(k个-1)-出租L(左)1L(左)2这样的话L(左)1先于L(左)2,G公司具有来自顶点的定向边L(左)1到的顶点L(左)2.的副本k个-lets可能存在于S公司,因此顶点之间可能有多条边。

图1
图1

序列AATAT的定向多重图。

排列与欧拉游动的对应

当我们扫描k个-让我们进去S公司一个接一个,我们也走在有向多重图中G公司从顶点到顶点。当所有k个-让我们扫描,每个边G公司只参观一次:步行欧拉学派另一方面,让一个欧拉式的人进来G公司,我们可以通过拼写出(k个-1)沿着行走的顶点的let(并放弃重叠)。自每个k个-让进来S公司对应于中的边G公司,每一次欧拉式的走进G公司对应于具有相同k个-让我们算作S公司坎德尔等人[15]表明,只要欧拉行走在相同的两个顶点开始和结束t吨与开头和结尾相对应(k个-1)-出租S公司,的-让所有1的计数≤k个被保存了下来。因此,生成统一的随机序列T型 k个 简化为生成均匀随机欧拉步入G公司t吨.

欧拉步行与树木群落的对应关系

对于一个欧洲人来说G公司,每个顶点v(v)属于G公司除了结束顶点t吨有一个最后边缘e v(v) 从中退出的v(v)最后一次。除以下顶点之外的所有顶点的最后一条边集t吨形成一个乔木植根于t吨:所有顶点都可以到达的有向生成树t吨.树木繁茂A类植根于t吨,从t吨最后边缘符合A类可以很容易地生成[2,15]:

  1. 1

    对于每个顶点v(v),收集边缘列表E类 v(v) 退出v(v).排列每个边缘列表E类 v(v) 分开保存e(电子) v(v) 列表上的最后一条边。

  2. 2

    浏览图表G公司根据边缘列表{E类 v(v) }:开始于(套u个),取第一个未标记的边(u、 v(v))从列表中E类 u个 ,标记边,然后移动到下一个顶点v(v)(套u个v(v)); 继续进行,直到所有边缘都标记好,步行结束于t吨.

在有向多重图中,欧拉游动和树状图之间有很好的对应关系:每个树状图都以t吨对应于完全相同数量的欧拉行走[,15]. 因此,生成均匀随机的欧拉走进去G公司t吨简化为在G公司植根于t吨在下一小节中,我们将讨论生成随机树状图的算法,其中一些算法又基于随机漫步,这很有趣。

生成随机树丛

在本小节中,我们回顾了现有的树木生成算法,并解释了我们对Wilson算法的选择[19,23]. 生成随机树状图和生成树的主要方法有两种:行列式算法和随机遍历算法。

行列式算法

行列式算法基于矩阵树定理[3,第二章,定理14]。对于图形G公司,特定边缘的概率e(电子)出现在均匀随机生成树中的是两个数字的比率:包含边的生成树的数目e(电子),以及生成树的总数。矩阵树定理允许人们通过计算图的组合拉普拉斯矩阵(或基尔霍夫矩阵)的行列式来计算图的生成树的确切数目。根据边的概率,通过反复收缩或删除边,可以生成随机生成树。

第一个行列式算法是由Guénoche给出的[14]和库尔卡尼[16]:对于的图形n个顶点和边,可以在中生成随机生成树O(运行)(n个)时间。这个运行时间后来被改进为O(运行)(n个) [7]. 科尔伯恩、梅尔沃尔德和纽菲尔德[8]简化了O(运行)(n个)时间算法,并表明运行时间可以进一步减少到O(运行)(n个2.376),二乘的最佳上界n个×n个矩阵[9].

随机遍历算法

随机遍历算法使用一种完全不同的方法来生成随机生成树。奥尔德斯[1]和Broder[5](在与Diaconis讨论矩阵树定理之后)独立地发现了随机生成树和随机游动之间的有趣联系:

模拟图中的均匀随机游动G公司从任意顶点开始直到访问了所有顶点。对于每个顶点v(v),收集边缘{u、 v(v)}对应于第一个入口v(v).集合T型边的一致随机生成树G公司.

对于图形G公司和一个顶点v(v)在其中,定义覆盖时间C v(v) (G公司)作为从开始的随机行走的预期步数v(v)需要访问的所有顶点G公司.Aldous-Broder算法的运行时间[1,5]在覆盖时间内明显呈线性。在改变生物序列的背景下,Kandel等人[15]扩展Aldous-Broder算法[1,5]生成覆盖时间内欧拉有向图的均匀随机树。威尔逊和普罗普[24]然后提出了一种在18个覆盖时间内生成一般有向图的均匀随机树状图的算法。

威尔逊算法

威尔逊[19,23]结果表明,使用循环弹出算法(cycle-popping algorithm)可以比覆盖时间更快地生成随机树形图和生成生成树,该算法模拟了循环增强的随机行走。对于图形G公司和两个顶点u个v(v)在其中,定义击球时间hu个,v(v)(G公司)作为随机行走的预期步数u个v(v)威尔逊算法的运行时间[19,23]在相应随机图的最大或平均命中时间内是线性的。作为威尔逊[19,23]值得注意的是,平均和最大命中时间总是小于覆盖时间,并且在某些图形中差异可能非常显著。因此,为了生成均匀随机树状图,Wilson的算法[19,23]优于Kandel等人的算法[15].

为了表示的完整性,我们在下面包含了Wilson算法的伪代码[19,23]:

具有根的RandomTree(第页)

1个用于1至n个

2      InTree中[]

3个下一步[第页]

4   InTree中[第页] 真的

5用于1至n个

6      u个

7而不是InTree中[u个]

8个下一步[u个] 随机后续任务(u个)

9         u个下一步[u个]

10      u个

11而不是InTree中[u个]

12         InTree中[u个] 真的

2013年u个下一步[u个]

14返回下一步

E类 u个 是从顶点退出的有向边集u个.函数RandomSuccessor(u个)选择均匀随机边(u、 v(v))来自E类 u个 ,然后返回顶点v(v).

与Aldous-Broder算法不同[1,5],它模拟从根到访问所有顶点的单个随机行走,Wilson的算法[19,23]模拟多个随机行走:从每个未访问的顶点开始,随机行走一直持续到它加入一个最初只包含根的生长乔木。随机行走遵循下一步[·]指针;每当再次遇到以前访问过的顶点时,就会形成一个循环并立即删除,因为下一步[·]指针被覆盖(在第一个while循环中)。当行走到达正在生长的树丛时,行走中的所有顶点都会作为一个分支加入树丛。

两种方法的比较

我们现在比较两种生成随机树状图的方法。Kandel等人[15]证明了欧拉有向多重图的覆盖时间n个顶点和边缘是O(运行)(n个2). 根据我们前面对覆盖时间和命中时间的讨论,可以得出Wilson算法的预期运行时间[19,23]最多在同一个多重图上O(运行)(n个2)同样,忽略了对数n因素。

对于多重图,数字边的数量可以任意大于n个个顶点。因此,看起来Colbourn等人的行列式算法[8],以确定性运行O(运行)(n个)时间或甚至O(运行)(n个2.376)时间,将是比随机行走算法更好的选择[15,19,23]. 然而,我们注意到,当行列式计算的中间值也可以很大。在当今典型的计算机系统中,浮点数的算术运算没有足够的精度来保证行列式算法数值计算的准确性和稳定性。随机行走算法[15,19,23]另一方面,只需要对小整数进行基本运算,而不需要这些数值问题。因此,我们决定实现Wilson的随机遍历算法[19,23]用于乔木生成。

实施详细信息

在本小节中,我们描述了欧拉算法实现的细节[2,15,19,23]用于生成k个-let表示随机序列。

Kandel等人的数据结构

正如Kandel等人建议的那样[15],Euler算法的简单实现[2,15]可以使用大小为的查找表σk个-1尽一切可能(k个-1)-lets作为有向多重图中的顶点G公司,然后构建大小为的邻接矩阵σk个-1×σk个-1对于中的边G公司。当两者同时存在时σk个是小常数,这是这种简单方法的空间要求,σ2k个-2个,可能看起来并不严重。然而,计算表明,即使对于σ=20(蛋白质的字母大小)和k个=3(典型的选择k个),所需空间总计

σ2k个-2个=204= 160, 000.

另一方面,蛋白质序列的典型长度低于1000。尽管序列本身可能仅存储在1KB中,但置换算法仍然需要数百倍的空间。k个进一步增加了:即使是看起来很天真的参数σ=20和k个=5,空间要求

σ2k个-2个=208> 168= 232

超过了32位计算机所能容纳的全部4GB内存!我们注意到Coward的两组参数[11]用于实验他的洗牌程序的只有

σ= 4,k个= 6,σ2k个-2个= 1, 048, 576

σ= 20,k个= 3,σ2k个-2个= 160, 000.

我们将在结果和讨论部分对uShuffle和shufflet进行比较,进一步讨论这一点。

线性空间中有向多重图的表示

为了使uShuffle程序具有可伸缩性,很明显,在实现过程中需要仔细的算法工程。正如我们在前面关于欧拉算法的小节中所讨论的那样,有向多重图G公司包含每个不同的顶点(k个-1)-引入S公司.由于(k个-1)-进入S公司确实是-k个+ 2,G公司最多有个-k个+2个顶点,因此-k个+连续之间有1条定向边(k个-1)-出租。这意味着G公司实际上长度是线性的序列的S公司待置换。有了合适的数据结构,uShuffle只需要线性空间。

在下面,我们首先解释有向多重图的构造和表示G公司然后解释了图构造后随机序列的生成。图的构造包括两个步骤:确定顶点集,然后添加有向边。

确定顶点

我们使用哈希表来确定顶点集。哈希表由一个大小为b条=-k个+2,数量(k个-1)进入S公司,以及每个bucket处的链接列表,以避免通过链接发生冲突[10]. 每个(k个-1)-出租x个=x个1x个2x个k个-1具有多项式哈希代码

小时 ( x个 ) = x个 1 k个 1 + x个 2 k个 2 + + x个 k个 2 2 + x个 k个 1 = ( ( ( x个 1 + x个 2 ) + + x个 k个 2 ) + x个 k个 1 ) , 数学类型@MTEF@5@5@@=feaafart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xI8qiVKYPFjYdHaVhbbf9v8qaqFr0xc9vqFj0dXdba91qpepeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaaeGaciGaaiaabeqaqaqabibiWaaaGcb aqbaeaabiWaaaqaaaiabdIgaOjabcIcaOiabdIha4jabcMcaPaqaaiabg2da9aqaaiabbdIha4 naaBaaalaaacqaXaqmaeqaaOGaemyyae2aaWbaaSqabeacqWGRbWAcqGHsislcqaIXaqmaaGccqGHRaWkcqWG4baEdaWgaaWcbaGaeGOmaidabeaakiabdggaHnaaCaaaleqabaGaem4AaSMaeyOeI0IaeGOmaidaaOGaey4kaSIaeS47IWKaey4kaSiaemiEaG3aaSbaaSbaaQaaabdUgaRjabgkHiTiabikdaYaqabaGccq WGHbqydaahaaWcbeqaaiabikdaYaaabgUaabGqaUcabUcabRiabdIha4naaBaaaaaacqaRbWAaqaGHsislc qaIX aqmaeqaaOGaemyyaegabaabaaGaeyypa0dabaGaeiikaGIaeiikaGIaS47IWKaeiakaGIaemiEaG3aaSbaaSqaaiabigdaXaqabaGccqWGHbqycqGHRaWkcqWG4baEdaWgaaWcbaGaeGOmaidabakiabcMcaPiabdggaHjabgUcaRiabl+UimjabgUcaRiabdIha4naaBaaaleaacqWGRbWAcqGHsislcqaIYaGmaeqaaOGaeiykaKIaemyyaeMaey4kaSIaemiEaG3aaSbaaSqaaabdUgaRjabgHiTiabigdaXaqabaGqaGGPaqkcqWGHbqycqGGSaaaaaaaaa@762E@

哪里= ( 5 1 ) 数学类型@MTEF@5@5@+=feaafart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacPC6xNi=xH8viVGI8Gi=hEeeu0xXdba9frFj0xb9qqpG0dXdb9aspeI8k8fiI+fsY=rqGqVepae9pg0db9vqaiVgFr0xfr=xfr=xc9adbaqaaeGaciGaaiaabeqaaeqabiWaaaGcbaGaeiikaGYAaOaaaaacqaI1aqnaSqabaGccqGHsislcqaIXaqmcqGGPaqkaaa@3081@ /2是黄金比率的倒数;的索引x个到桶数组是

(x个) =小时(x个b条国防部b条.

将哈希表初始化为空,然后尝试插入(k个-1)-逐个进入哈希表。如果(k个-1)-让它是同类中的第一个,它被分配一个新的顶点编号,然后插入到哈希表中;它的序列起始索引S公司也会被记录。如果(k个-1)-let之前已插入,但未插入到哈希表:它的顶点数和序列的索引S公司是从第一个(k个-1)-让它的种类。插入后,我们可以从指定的最大顶点数推断出有向多重图中的顶点总数。然后分配顶点的内存。

添加定向边

为了添加有向边,我们使用邻接列表表示法来避免邻接矩阵的过度内存需求。在邻接列表表达法中,每个顶点需要维护两个边列表:传入边列表和传出边列表。传出边列表是生成欧拉行走所必需的[2]. 当Kandel等人的算法时,传入边列表是生成树状图所必需的[15]使用(如Coward的实现[11])。我们使用威尔逊算法[19,23]用于生成树状结构。正如我们在上一节中所讨论的,Wilson的算法[19,23]比Kandel等人的算法更快[15]. 此外,我们在这里注意到Wilson的算法[19,23]与Kandel等人的算法相比还有一个优点[15]在易于实施方面。而不是一个向后的从终点随机行走t吨达到Kandel等人算法中的所有其他顶点[15],威尔逊算法[19,23]使用多个向前地从每个未访问的顶点随机行走以加入植根于t吨:仅输出边列表就足以生成欧拉行走和树丛。

表示边缘列表和管理内存

为了获得最大效率,我们将每个边列表实现为一个顶点数组。传出边的数量因顶点而异;如果我们为每个顶点分配一个固定大小的数组,那么在最坏的情况下,我们必须使每个数组足够大,以容纳所有边,由此产生的空间需求将变为长度的二次方序列的S公司。我们当然可以首先计算每个顶点的传出边的数量,然后为每个顶点分配一个足够大的单独数组。然而,这需要我们为每个顶点调用一次相对昂贵的内存分配函数。

在我们的实现中,我们为所有边分配一个大数组(边的总数为-k个+1),然后将块打包到各个顶点。为了实现这一点,我们首先扫描序列S公司要计算每个顶点的传出边数,请将每个顶点的数组(传出边列表)指向大数组的连续偏移。通过这种优化,内存分配的数量减少到只有4个:一个用于hashtable bucket数组,一个用于(k个-1)-将一个顶点数组和一个边数组作为散列表条目。一旦构造了定向多重图,就可以释放桶数组和哈希表项的内存。

图构造后的序列生成

在构造了有向多重图之后,我们可以分三步生成随机序列。如前一节所述,我们需要首先模拟循环增强的随机游动[19,23]为了生成树状图,接下来排列各个边列表,同时保持最后一条边,然后模拟由边列表引导的欧拉行走,并沿行走输出序列。由于每个边列表都是作为数组实现的,因此可以非常有效地执行置换。沿着行走输出随机序列也很容易,因为每个顶点都保持其在输入序列中第一次出现的起始索引。

uShuffle工具的软件组织和用户界面

在本小节中,我们将描述uShuffle工具的软件组织和用户界面。

C库和命令行工具

uShuffle的初始实现是用C编程语言实现的。uShuffle的C版本由两个组件组成:一个uShuffle库(uShuffle.C和uShuffle.h)和一个命令行工具(main.C)。

在典型场景中k个-为每个输入序列生成let-preserving随机序列。uShuffle程序的图形构建阶段只需对多个输出序列执行一次。为了给用户提供优化选项,我们在uShuffle库中导出了三个接口函数:

void shuffle(const char*s,char*t,int l,int k);

void shuffle1(const char*s,int l,int k);

无效洗牌2(字符*t);

函数shuffle接受四个参数:s是要置换的序列,t是输出随机序列,l是s的长度,k是let大小k个函数shuffle只需首先调用shuffe1,然后调用shuwe2:shuffle1实现有向多重图的构造;shuffle2实现了有向多重图中的循环随机游动和随机序列的生成。随机排列的统计行为在很大程度上取决于随机数生成器。

胆小鬼[11]注意到随机数生成器在各种平台上的默认实现通常不令人满意,所以他使用一种可以说更好的算法实现了自己的生成器。我们注意到,有许多随机数生成算法,并且不断有新算法被提出:一种算法是否优于另一种算法可能非常主观。我们没有限制用户使用特定的实现,而是将默认生成器设置为标准C库中的随机函数,然后导出一个接口函数,以便高级用户自定义生成器:

typedef长(*randfunc_t)();

无效集randfunc(randfunc_t randfunch);

命令行uShuffle工具是uShuff库的最小前端,它演示了库的典型用法。它有以下四个选项:

  • s<string>指定输入序列,

  • n<number>指定要生成的随机序列数,

  • k<number>指定let大小,

  • seed<number>指定随机数生成器的种子。

Java小程序

uShuffle程序被移植到Java编程语言。除了有一个库和命令行工具外,uShuffle程序的Java版本还可以在web浏览器中作为小程序运行。我们参考图2uShuffle Java小程序的屏幕截图1:小程序的界面很小,由三部分组成:顶部的输入文本区域、底部的输出文本区域和中间的控制面板。控制面板包含两个文本字段和一个按钮。最大出租规模k个和数字n个的输出序列可以在两个文本字段中设置。单击“Shuffle”按钮时,小程序从输入文本字段获取输入序列,去掉空白,生成n个保留k个-让我们计数,然后在输出文本区域中输出序列。n个>1:每个输出序列前面都有一个注释行,其中包含从1到n个.

图2
图2

uShuffle Java Applet的屏幕截图。

uShuffle Java小程序将所有输出序列保存在内存中,以便在输出文本区域中显示。当数字n个输出序列和输入序列长度例如,过大,n个=100000000和=100,则保存输出序列所需的总内存可能超过Java虚拟机(JVM)的最大堆大小,并且小程序可能会挂起。这不是我们程序中的错误,而是由于JVM的限制;尽管如此,我们还是准备了一个网页来指导用户如何增加JVM的最大堆大小。

C#/Perl/Python版本

uShuffle程序也被移植到C#编程语言。Perl和Python是用于生物信息学的流行编程语言;它们允许与用C编写的程序轻松集成。我们没有在源代码级别将uShuffle程序移植到Perl和Python,而是准备了两个网页来指导用户如何使用uShuff库扩展Perl和Cython环境。

结果和讨论

我们进行了两组实验来测试uShuffle工具的两种主要形式的性能:首先对uShuff C库的性能进行基准测试,然后将uShufle Java applet的性能与Coward的shuffle程序进行比较[11].

uShuffle C库的性能

我们在台式PC上测试了uShuffle C库2测试数据由真实的生物序列和人工生成的随机序列组成。

真实生物序列实验

真正的生物序列是从两个来源获得的:首先,从人类蛋白质参考数据库中采集了152个蛋白质序列(共91262个氨基酸),152个分子类别中的每一个都有一个序列;其次,从补充数据中提取了家鼠的69个microRNA前体序列(共4773个核苷酸)4Bonnet等人[4].

我们对这些真实生物序列的实验表明,uShuffle库非常有效:在一秒钟内,它可以为152个蛋白质序列中的每个序列生成700个双重保存的随机序列,或者(ii)为69个RNA序列中的每一个序列生成12000个双重保存的随机序列。

人工生成随机序列的实验

为了分析uShuffle在各种参数集下的性能,我们还对uShuff在人工生成的随机序列上进行了系统测试。为了简单起见,序列长度是2的精确幂12至224也就是说,从4000到16000000。这些数字有些随意;只要计算机有足够的内存来存储输入序列,并且有我们实现所需的一些额外(如实现部分中所讨论的,非常小的)内存,就没有什么可以阻止用户在很长的序列上运行uShuffle,即使是在基因组范围内。

对于每个序列长度,在英文字母[a-z]上生成64个均匀随机序列作为测试序列;对于每个测试序列,64k个-然后通过uShuffle生成let-preserving随机序列。uShuffle生成64×64=4096的总运行时间k个-记录每个序列长度的let-preserving随机序列。在测试程序中放置了两个getrusage系统调用,以夹住正在进行基准测试的代码区域;两个时间戳的差异用于计算运行时间。

我们参考图uShuffle程序总运行时间的对数图k个=2和k个在不同序列长度下=3。该图显示,uShuffle程序的运行时间在待洗牌序列的长度上基本上是线性的。

图3
图3

uShuffle的运行时间 k个 =2和 k个 = 3uShuffle在不同序列长度下的运行时间(毫秒)k个=2和k个= 3.

绝对运行时间并不能很有效地证明uShuffle程序的极端效率。我们参考图4对于更具说明性的比率图。对于k个=2和k个=3,对于每个序列长度,图中显示的不是uShuffle程序的绝对运行时间,而是两个运行时间的比率:

图4
图4

运行时间比率 k个 =2和 k个 = 3uShuffle和简单置换法在不同序列长度下的运行时间比率k个=2和k个= 3.

  1. 1

    uShuffle程序生成k个-let-preserving序列,以及

  2. 2

    简单置换法的运行时间[10](在“背景”部分中查看)来洗牌相同数量的随机序列没有保留k个-让我们数数。

比率图显示uShuffle程序的平均运行时间仅为简单置换方法的1.5倍k个=2,仅2次k个= 3. 简单的置换方法是最小的:对于输入序列的每个位置,它只执行一个随机函数调用和一个交换。另一方面,uShuffle程序执行了更多的工作;尽管64k个-每个测试序列的let表示随机序列是由一个shuffle1和64个shuff函数调用生成的,以避免冗余的多图构造,每个shuffle函数调用仍然包括通过循环加速随机游动生成树状图[19,23]以及由简单排列的单个边缘列表引导的欧拉行走的生成。鉴于uShuffle程序和简单置换方法的复杂性的对比,它们的运行时间的小比率是显著的。细心的读者会从图中注意到一个有趣的事实4,当序列长度增加到224(约1600万),uShuffle程序的运行时间k个=2为偶数小于简单的置换方法!“奇怪”的现象让我们困惑了很长一段时间,直到我们最终确信这不是一个bug,而是一个特性。我们注意到,在每个步骤中,简单置换方法随机交换分散在2个大数组中的两个元素24元素。另一方面,uShuffle程序在小型多重图中执行随机漫步(最多26个顶点k个=2并且在[a-z]字母表上)并且排列各个边缘列表(每个边缘列表具有大约224/26个元素)。uShuffle程序的内存引用更多地方的而不是简单排列法。具有现代内存体系结构的计算机通过复杂的缓存方案利用本地内存引用积极优化代码,从而提高uShuffle程序的性能。

我们参考图5uShuffle程序在不同参数值下的运行时间k个,其中测试序列长度固定为1024。uShuffle程序的运行时间峰值为k个=4,约为其运行时间的三倍k个=2,然后逐渐减小k个增加,最后在k个=1024,因为在序列长度为1024的情况下,唯一保留1024个let的随机序列是输入序列本身。此图显示uShuffle程序对以下所有可能的值都是有效的k个.

图5
图5

不同Let大小uShuffle的运行时间 k个.uShuffle在不同let大小的固定序列长度下的运行时间(毫秒)k个.

uShuffle Java Applet与shuffle的比较

Euler算法还有另外两种实现。Clote等人的洗碗节目[6]由Altschul和Erickson实现Euler算法的原始版本[2]. 硬编码的洗牌RNA序列保留二核苷酸计数,洗牌不是任意字母表和大小的通用工具。另一个节目,Coward的shufflet[11]实现了Kandel等人对Euler算法的改进版本[15]对于任意let大小k个正如我们在“实现”一节中所解释的,Kandel等人[15]虽然优于Altschul和Erickson的算法[2],仍然不如Wilson的算法[19,23];此外,它的查找表数据结构对于大字母表和let大小来说效率很低。

就功能而言,shufflet实现[11]更接近我们的uShuffle实现。Shufflet是用C编程语言编写的,并作为web应用程序托管(但已离线)。我们无法对uShuffle和shufflet进行全面比较。然而,Coward[11]提到了在Digital DEC/Alpha 2100 web服务器上进行的两个实验:

  1. 1

    10000个核苷酸的DNA序列的100次洗牌k个=6大约需要2.5秒;

  2. 2

    1000个氨基酸的蛋白质序列的100次洗牌k个=3不到1秒。

我们在Apple iMac计算机(运行MacOS 10.4.9、Firefox 2.0.0.3和Java 1.5.0的2 GHz PowerPC G5)上使用uShuffle Java applet进行了类似的实验:

  1. 1

    10000个核苷酸的DNA序列的1000次洗牌k个=6大约需要1.5秒;

  2. 2

    1000个氨基酸的蛋白质序列的4000次洗牌k个=3不到1秒。

假设这两台计算机的性能相当,我们估计在核苷酸实验中,我们的uShuffle Java小程序比shufflet快大约15-20倍(k个=4),在氨基酸实验中,其速度约为shufflet的40倍(k个= 20).

我们当然理解这样一个比较的困难:1999年的网络服务器与2005年的台式计算机;本机中的C程序与虚拟机中的Java小程序。然而,这一比较说明了我们的uShuffle Java applet对于大let大小具有更好的可伸缩性。这两个性能比率之间的差异(15–20与40)表明,uShuffle即使对于较大的let大小也仍然有效,而随着let大小的增加,shuffle变得更加低效,这很可能是因为Kandel等人使用了低效的查找表数据结构[15].

结论

uShuffle工具基于高级图形算法,经过精心设计,效率极高。它通过允许任意字母大小和let大小来实现最大的灵活性,并以多种形式提供给不同类型的用户。我们相信uShuffle对生物信息学界来说是一个有用的工具。

可用性和要求

项目名称:uShuffle。

项目主页:网址:http://www.cs.usu.edu/~mjiang/u缓冲/

操作系统:独立于平台。

编程语言:C、Java、C#、Perl、Python。

其他要求:无。

许可证:FreeBSD。

非学者使用的任何限制:无。

注释

1MacOS 10.4.9、Firefox 2.0.0.3、Java 1.5.0。

2Dell XPS M1710:2.1 GHz Intel双核2处理器,2 GB双列直插式RAM;Microsoft Windows Vista Business Edition,Cygwin 1.5.24,gcc 3.4.4,带-O3选项。

http://www.hprd.org/moleculeClass网站.

4http://bioinformatics.psb.ugent.be/supplementary_data/erbon/nov2003/.

工具书类

  1. Aldous DJ:均匀生成树和均匀标记树的随机行走构造。SIAM离散数学杂志1990, 3(4):450–465.

    第条 谷歌学者 

  2. Altschul SF,Erickson BW:核苷酸序列比对的重要性:一种保留二核苷酸和密码子用法的随机序列置换方法。分子生物学与进化1985, 2(6):526–538.

    中国科学院 公共医学 谷歌学者 

  3. Bollobás B:现代图论。施普林格;2002

    谷歌学者 

  4. Bonnet E,Wuts J,RouzéP,van de Peer Y:微RNA前体与其他非编码RNA不同,其折叠自由能低于随机序列的证据。生物信息学2004, 20(17):2911–2917.

    第条 中国科学院 公共医学 谷歌学者 

  5. Broder A:生成随机生成树。第30届计算机科学基础年会论文集(FOCS’89)1989, 442–447.

    第章 谷歌学者 

  6. Clote P,Ferr’e F,Kranakis e,Krizanc D:结构RNA的折叠能量低于相同二核苷酸频率的随机RNA。核糖核酸2005, 11: 578–591.

    第条 公共医学中心 中国科学院 公共医学 谷歌学者 

  7. Colbourn CJ,Day RPJ,Nel LD:对图的生成树进行排名和排序。算法杂志1989, 10: 271–286.

    第条 谷歌学者 

  8. Colbourn CJ,Myrvold WJ,Neufeld E:两种用于解排树丛的算法。算法杂志1996, 20(2):268–281.

    第条 谷歌学者 

  9. Coppersmith D,Winograd S:通过算术级数进行矩阵乘法。符号计算杂志1990, 9: 251–280.

    第条 谷歌学者 

  10. Cormen TH、Leiserson CE、Rivest RL、Stein C:算法简介第二版。麻省理工学院出版社;2001

    谷歌学者 

  11. Coward E:洗牌:洗牌序列,同时保留k-let计数。生物信息学1999, 15(2):1058–1059.

    第条 中国科学院 公共医学 谷歌学者 

  12. 惠誉WM:随机序列。分子生物学杂志1983, 163: 171–176.

    第条 中国科学院 公共医学 谷歌学者 

  13. Freyhult E,Gardner PP,Moulton V:RNA折叠措施的比较。BMC生物信息学2005, 6: 241.

    第条 公共医学中心 公共医学 谷歌学者 

  14. 盖诺切A:随机生成树。算法杂志1983, 4: 214–220.

    第条 谷歌学者 

  15. 坎德尔D、马蒂亚斯Y、昂格R、温克P:混乱的生物序列。离散应用数学1996, 71(1–3):171–185.

    第条 谷歌学者 

  16. Kulkarni VG:生成随机组合对象。算法杂志1990, 11(2):185–207.

    第条 谷歌学者 

  17. Le S-Y,Chen J-H,Currey KM,Maizel JV Jr:预测重要RNA二级结构的程序。计算机在生物科学中的应用1988, 4(1):153–159.

    中国科学院 公共医学 谷歌学者 

  18. Nussinov R:Sone规则DNA中核苷酸的排序。核酸研究1980, 8(19):4545–4562.

    第条 公共医学中心 中国科学院 公共医学 谷歌学者 

  19. Propp JG,Wilson DB:如何从一般马尔可夫链中获得完全随机样本,并生成有向图的随机生成树。算法杂志1998, 27(2):170–217.

    第条 谷歌学者 

  20. Rivas E,Eddy SR:在检测非编码RNA时,仅二级结构通常没有统计学意义。生物信息学2000, 16(7):583–605.

    第条 中国科学院 公共医学 谷歌学者 

  21. Seffens W,Digby D:mRNAs的负折叠自由能大于洗牌或密码子选择随机序列。核酸研究1999, 27(7):1578–1584.

    第条 公共医学中心 中国科学院 公共医学 谷歌学者 

  22. Smith TF,Waterman MS,Sadler JR:核酸序列功能域的统计表征。核酸研究1983, 11(7):2205–2220.

    第条 公共医学中心 中国科学院 公共医学 谷歌学者 

  23. Wilson DB:生成随机生成树的速度比覆盖时间更快。第28届ACM计算理论年会(STOC'96)会议记录1996, 296–303.

    谷歌学者 

  24. Wilson DB,Propp JG:如何在覆盖时间内从通用马尔可夫链中获取精确样本,并从有向图中获取随机生成树样本。第七届ACM-SIAM离散算法年会论文集(SODA'96)1996, 448–457.

    谷歌学者 

  25. Workman C,Krogh A:没有证据表明mRNA的折叠自由能低于具有相同二核苷酸分布的随机序列。核酸研究1999, 27(24):4816–4822.

    第条 公共医学中心 中国科学院 公共医学 谷歌学者 

下载参考资料

鸣谢

我们感谢匿名评论员的宝贵意见。这项研究得到了国家科学基金会拨款DBI-0743670和犹他州立大学拨款A13501的部分支持。

作者信息

作者和附属机构

作者

通讯作者

与的通信江明辉.

其他信息

作者的贡献

MJ设计了软件和实验,实现了uShuffle程序的C和Java版本,并编写了技术报告。JA将Java程序移植到C#,研究软件许可证,并设计uShuffle徽标。JG编写了测试程序,执行了实验,并编写了Perl集成的指令。MM编写了Java小程序接口和Python集成说明。所有作者都审查了源代码,并为主页的构建做出了贡献。

作者提交的原始图像文件

权利和权限

本文由BioMed Central Ltd.授权发布。这是一篇根据知识共享署名许可条款发布的开放存取文章(http://creativecommons.org/licenses/by/2.0)它允许在任何介质中不受限制地使用、分发和复制原始作品,前提是正确引用了原始作品。

转载和许可

关于本文

引用这篇文章

Jiang,M.、Anderson,J.、Gillespie,J。等。uShuffle:在保留k-let计数的同时,对生物序列进行洗牌的有用工具。BMC生物信息学 9, 192 (2008). https://doi.org/10.1186/1471-2105-9-192

下载引文

  • 收到以下为:

  • 认可的以下为:

  • 出版以下为:

  • 内政部以下为:https://doi.org/10.1186/1471-2105-9-12

关键词