（概率论与数理统计专业论文）多重序列比对和基因芯片数据分析.pdf

摘要 随着人类基因组( 测序) 计划( h u m a ng e n o m ep r o j e c t ) 的实现以及分子生物 学相关学科的迅猛发展，越来越多的动植物、微生物基因组序列得以测定，基因序列 数据正在以前所未有的速度迅速增长与此同时综合了分子生物学，半导体微电子， 激光，化学染料等领域的最新科学技术的基因芯片，也在以其无可比拟的信息量、高 通量、快速、准确地分析基因的能力，也在基因功能研究、临床诊断及新药开发等方 面显示了巨大的潜能而生物信息学分析是大规模基因信息( 序列信息，表达信息) 分析的主要方法本文主要结果由两部分组成； 第一部分是针对基因序列的快速多重比对算法及应用 论文的第一和第二章针对基因序列的多重序列的比对问题，根据序列突变与比 对的”模代数”结构理论 1 ，并在序列两两比对的基础上，应用系统聚类等方法， 给出了同源多重序列的超级快速比对算法( 简称s m a ) 该算法可适用于大规模( 如m 5 0 0 ，n 1 0kb p ) 的同源多重序列的比对计算，s m a 的程序和测试数据 我们已在网上公布，并提供比对计算服务我们分别对8 3 3 0kb p 的s a r s 序列 与7 0 6 l okb p 的v - l 序列进行比对，主频3 ogp c 机上完成多重比对的计 算时间分别为2 1 分钟与3 4 小时，s m a 算法在速度上对h m m e r 有明显优势，而 且根据相似比等各项优化指标测试，结果不差于h m m e r 多重序列在完成比对后，对序列的结构分析是多重序列比对后的主要与关键问 题论文的第三章给出了利用突变网络理论来研究多重比对后基因的突变分析利用 该理论和正交化方法来研究基因突变，我们可以得到基因组多重序列比对的突变过 程更清楚的描述我们以s a r s 病毒基因组为例，说明它们的突变网络系统模型与 正交化运算，绘制了s a r s 基因突变图，并由此得到s a r s 病毒从早期传播到爆发 的基因突变过程 第二部分是针对基因芯片表达数据的聚类算法和网络调控分析 第四章针对乳腺癌在转移过程中的基因表达特征的聚类分析法分析，我们改进 了k - m e a n s 聚类算法，给出了k r - m e a l l s 算法，使之具有自动搜索聚类数的功能，并 且有助于改善k - m e a i l s 算法的聚类结果陷入最小值的状况通过对平均聚类误差指 标的比较，k r m e a 丑s 要明显优于k - m e a n s 算法本文所得到的结果可供乳腺癌诊断 与病变分析参考，同时可以应用于小型基因检测芯片的制备，也可以用于构建基因网 摘要 络调控图 第五章针对表达数据提取基因调控矩阵从而构建基因网络的问题，我们通过线 性微分方程模型可以初步构建基因网络，了解网络结构，提取最显著的信息然而由 于分子生物学的条件限制或者数据来源的限制，导致实验数据不充分，使方程组无 解我们使用三次样条方法，对2 6 例临床、病理资料完备的具有淋巴结转移的乳腺 癌基因表达数据进行插值处理，使表达数据满秩，从而使用最小二乘法解出加权矩 阵，构建初步的表达基因调控网络通过对构建的基因网络的初步生物学和医学分析 表明；乳腺癌转移的形成是由多基因异常引起多条传导通路异常，致使细胞恶性转化 的结果，这与生物学上公认的看法是相一致的利用此线性模型方法对基因表达谱进 行分析具有一定可行陛，在认识乳腺癌转移机制，乳腺癌诊断和治疗方面具有一定的 理论和应用价值 关键词；多重比对算法，大规模基因组比对，突变网络分析，基因芯片表达数据 分析，聚类分析，调控网络分析 a b s t r a c i t bc o m p 甜el a r g en u m b e r 8o fg e n o m i cs e q u e n c e so fr e l a t e d 、，i r u s ，s u c ha sh i v ， b j o l o g i s t sh a v ea ni n c r e 崩n gn e e d 妇am e t h o dt h a tc a ne 伍c j e n t l yh a n d l eh l l n d r e d s ，e v e nt h o u 8 a n d s ，o fg e n o l i l i cs e q u e n c e sa c c u r a t e l ye n o u g ht oc o r r e c t l ya l i g i l t h e s ec o n s e r v e df e a t u r e si nt h el s t2 n dc h a p t e r ，w ei r l t r o d u c ean e wa n de m c i e n t t o o in 虹1 e ds m at h a tc a ne 8 s i i ya c c o m m o d a t ei 盯g e _ s c a i ev i r u 8g e n o m i cs e q u e n c e s a h i g h - t k o u g h p u tt 豁t0 n7 0 6h 一1g e n o i n i cs e q u e n e e s8 岍唱t h a ts m a i 8 m c h f a s t e rt h a nt h ea v a i l a b l ep r o 罂a m s 谢t ha tl e a s tt h es 啪ep e r f o r m a n c e s m ai sa g o o di m p r o v e m e n to fe ) ( i s t i n ga l g o r i t h m sf o rh i g h - v o l u m em u l t i p l e8 e q u e n c ea l i g n - m e n t i to f f e r sa n0 p t i o nt h a tp r o v i d e si m p r o v e ds p e e da n da c c u r a c yc o m p a r e dw i t h c u r r e n t l ya v a i l a b l ep r o g r a j l l s b a s e do ni n u l t i p l es e q u e n c e 啦i g n m e n t ，t h ea n a l y s i 8o ft h eb i o l o 盯8 e q u e n c e m u t a t i o ns t r u c t u r en e e d sas e r i e so fc o m p u t a t i o 璐t h es y 8 t e m a t i ct r e ea n dk a s t d i s t a n c et r e ec a nb ec o n s t r u c t e dw h i c hc o u l db er e p r e s e n t e do n l yb yat r e et o p o l o 够 i nt h e3 r dc h a p t e rw ep r e s e n tan e wt o p o l o g ys t r u c t u r et oa n a l y z et h em u t a t i o n 8 t r u c t u r e t h en e wm o d e ll sn a m e da sm l l t a t i o nn e 七w o r k 科8 t e m t h em a j o rp r 如 1 e mo fm u t a t i o nn e 七w d r ks v s t e mi 8t od e t e r 血n et h er e l a t i o 璐o ft h ed i 船r e n tm u t a - t i o n 帅e s t h es i i n p l e s tr e l a t i o ni sd e 矗n e da so r t h o g o n a 】i z a t i o n ，w h i c hi n e a n st h e r e i s oo v e r l 印b e t w e e n 栅om u t a t i o nr e g i o n s t h es i m p l m c a t i o no fo r t h o g o n a l i z a t i o n i ni n u t a t i o nn e t w o r ks y s t e mi 8t h ef o c u si nt h i sc h a p t e r w bp r e 8 e n tt h eb a s i ct h 盼 r e m0 fo r t h o g o n a l i z a t i o nt h e o r yi n m t a t i o nn e t w o r ka n d 也ea p p h c a t i o ni na n a l y s i s o fs a r sg e n o n l 呈cs e q u e n c e s t h ec o m p u t a t i o no fo r t h o g o n a l i z a t j o nd i s c o v e r st h e i r m t a t i o np r o c e s so fs a r sf r o mt h ee a r l yp h a s et om i d d l ep h a s e i nt h e4 t hc h a p t e r ，a1 i n e a r 出任e r e n t i 出m o d e li se m p l o ”dt oc o n s t r u c tr o u g h g e n e t l cn e t w o r ko fm a i n m a r yc a r c i n o m ag e n er e l a t e dt om e t a s t 出s i s c u b i cs p l i n e i n t e r p o l a t i o ni su s e dt oi n t e r p 0 1 a t et h er e 口e s 8d a t ai n t oe x p e r i m e n tg e n ee x p r e s s i o n d a t a h e n c et h el i n e a rm o d e lc a nb ee 8 t i m a t e db yl e s ts q u a r e sm e t h o d w i t ht h e a n a l y s i so ft h eg e n e t i cn e 七w 0 r kc o n s t r u c t e db yt h el i n e a rm o d e l ，w ef o u n dp o t e n t i a l t r a 璐血t i o np a t l l w a y sc a l l s e db ym u l t i p l eg e n ev a r i a n c e a b s t r a c t ，ep r o p o s ea na m e l i o r a t e dk - m e a n sa l g o r i t h mi nt h e5 t hc h a p t e r ，n a m e dk r m e a n 8 ar e c u r s i v ei t e r a t i o ns t e pi 8i n t r o d u c e dt ot h ek r m e a i 玛a l g o r i t h m t h i 8 s t 印a _ v o i d sg u e s s t i m a t i n gt h e 后i nk - m e a n sa l g o r i t h ma n d p u t e 七a u t o m a t i c a l l y a c c o r d i n gt h ea v e r a g ee r r o ro ft h ed u s t e r i n g ，1 【r - m e a n ss h o w sb e t t e rp e r f o r m a n c e c o m p a r i n g 耐t hk _ m e a n s k e yw o r d s ：m u l t i p l ea l g m r n e n t ，m u t a t i o nn e t w o r ks y s t e m ，g e n ee x p r e 8 s i o n d a t ac l u s t e r i n ga n a l y s i s ，g e n ee x p r e s s i o nr e g u l a t o r yn e 怕帕r ka n a l y s i s y9 6 8 3 1 2 南开大学学位论文版权使用授权书 本人完全了解南开大学关于收集、保存、使用学位论文的规定， 同意如下各项内容：按照学校要求提交学位论文的印刷本和电子版 本；学校有权保存学位论文的印刷本和电子版，并采用影印、缩印、 扫描、数字化或其它手段保存论文；学校有权提供目录检索以及提供 本学位论文全文或者部分的阅览服务：学校有权按有关规定向国家有 关部门或者机构送交论文的复印件和电子版；在不以赢利为目的的前 提下，学校可以适当复制论文的部分或全部内容用于学术活动。 学位论文作者签名：鍪写参吹 6 一卵丰l 月f 日 经指导教师同意，本学位论文属于保密，在年解密后适用 本授权书。 指导教师签名：学位论文作者签名： 解密时间：年月日 各密级的最长保密年限及书写格式规定如下： 南开大学学位论文原创性声明 本人郑重声明：所呈交的学位论文，是本人在导师指导下，进行 研究工作所取得的成果。除文中已经注明引用的内容外，本学位论文 的研究成果不包含任何他人创作的、已公开发表或者没有公开发表的 作品的内容。对本论文所涉及的研究工作做出贡献的其他个人和集 体，均已在文中以明确方式标明。本学位论文原创性声明的法律责任 由本人承担。 学位论文作者签名：喜为锔像 ) 口。厂年7 乙月j 日 第一章引言 1 1 多重比对算法目前的进展和存在的问题 随着高通量的d n a 序列数据的从测序仪器中源源不断的产生，早期通过人眼的 目测比对法得到的最优，甚至局部最优的多重比对对于现代的海量数据是不可能的 只有i 酿生大规模基于计算机运行的多重序列序列比对算法，例如本文提出的s m a 算 法，才能解决随着数据量和复杂性增加带来的问题 多重序列比对作为一个重要的生物信息学工具，有着广泛的应用从分子生物学 到遗传学，多重序列比对在分析同源性，预测新基因，发现突变位点，预测蛋白质空 间结构，乃至推测进化树中，都有它的用武之地 分子生物学揭示了在不同的组织间，甚至不同的物种间，其d n a 序列仍然具有 很大的相似性，而相同的序列基本都有相似的功能通过各物种同源物序列的多重序 列比对，可以发现进化中的保守区域，保守位点及建立种系进化关系，并将有助于对 其未知功能做进一步推测 研究相似的序列之间的突变和重组，也是相当重要的课题如抑癌基因通过点突 变，易位，插入，染色体部分甚至整个臂的缺损，而失活，使细胞增殖失控，导致细 胞的恶性转化【2 】利用d n a 分析技术检测肿瘤的染色体定位及丢失情况，可以得 到较为肯定的抑癌基因比如：p 5 3 基因位于人第1 7 号染色体短臂上( 1 7 p ) ，c d n a 全长1 6 k b 2 0 k b ，有u 个外显子，突变集中在第5 8 号外显子，p 5 3 基因编 码分子量为5 3 k d ，p 5 3 的突变，使d n a 损伤的细胞可能逃避正常凋亡过程 3 】通 过对大量正常和突变序列的多重序列比对分析，可以发现基因的突变位点，推测染色 体重组情况 与两两比对比较，多重序列比对还可以解决一些在两两比对中解决不了的问题， 比如：寻找多序列的共同保守区域，以及共同保守区域中的局部的子保守区域在多 序列中的共同保守区域中，可能会有局部上的差异，从而形成若干个分类，这种差异 只有在多重序列比对下才能得到全面的观察和分析 由此可见，多重序列比对在生物信息学的各个领域中都有应用，并且随着生物序 列数据量的爆炸性增长，作为其主要分析工具的多重序列比对的重要性，也越来越收 】 第一章引言2 到大家的重视 多重序列的比对的最优解问题在一些文献中被列为是个非易计算问题，如【4 】 文称之为生物计算中的第一未解决问题，文献 5 ，6 ，7 】等文把它归结为n p - h a r d 与 m a x s n ph a r d 问题，把它的计算复杂确定在d ( 2 ”m “) 范围( 其中m ，n 分别为多 重序列的重数与长度) ，因此在理论上实现多重序列的比对最优解的易计算是一个非 常困难的问题 另一方面，由于多重序列的比对问题的重要性，近几年内有许多多重序列的比对 问题十分活跃，许多算法，软件包与比对结果大量出现，比如m u l a u n 【8 ，c l u s t a l w 【9 ，h m m e r1 10 ，t c o f f e e 【1 1 】，p o a 【1 2 ，m u l t i l a g a n 【13 ，c h a o s d i a u g n 2 【1 4 ，1 5 ，m a d 【1 6 ，m u s c l e 【l7 和m u m m e r 【l8 】 但是，大多数算法和工具是针对少量的长序列，或者大量的短序列比如m u l t i _ l a g a n 和c h a o s d i a l i g n 2 就是针对全基因组多重序列比对，但是只能容纳3 0 条左右的序列而象经典的c 1 喊出w ，h m m e r ，m u l a l i n 和最近出现的p o a ，t c o 丘e e 则是针对大量的在l k b 左右的短基因序列如果序列长度达到l o k b ，基本上计算时 间长的不能忍受( 大于一周) 针对大量的较长的多重序列比对算法目前在国际上相 当少而针对这个问题的多重比对算法也有相当的需求例如：据u n a i d s 组织报 道，目前有三千九百万的艾滋病患者因此有大量病例的爱滋病毒基因组序列已经被 或者将被测序，比如现在至少有7 0 6 条h i v - l 的基因组序列已经被测序，而且被测 序的艾滋病毒序列增长迅速目前几乎没有算法可以容纳7 0 0 多条的长度1 0 k b 左 右的艾滋病病毒序列组，而本文第一部分的目的就是要提出个新的算法，可以较优 的对这样大规划的序列组进行多重比对 1 2 核酸( n u c k i ca c i d ) 和脱氧核醣核酸( d n a ) 序列 生物序列包括核酸序列( d n a 与r n a ) 与蛋白质序列核酸的多重序列比对和 蛋白质的序列的多重比对有相似的地方，也有不同的地方，但是，本文主要讨论核酸 序列的多重比对以下文中所用的”序列”，如无特别指出，均指核酸d n a 序列 核酸是以核苷酸为一系列小分子单元所按照前后顺序排列而成的巨分子，是圾 胞分子量最大的功能性分子，包括d n a 及r n a ；其主要功能为遗传信息的存 第一章；l 窘 3 贮，传递与波达，是现代分子生物学的研究主要对象 核苷酸( n u c l e o t i d e ) 是由磷酸、五碳糖及碱基三部份组成，每一部份又可有几 种可能：例受玎赢碳糖部分，可能怒拔醣( 磁b o s e ) 或是聪戴棱醣( d e o x y r i b o s e ) ，矗f ：种 差造或了n a 与r n a 在毪蒺与褥遗主懿不建 由四种藏本化学物质的d n a 分子形成双螺旋络梅缀令箭成了包含大量遗僚溶慧 的染色体，逡暇种基本化学物质分照腺嘌呤a ( a d e n i n e ) ，胸腺嘧啶t ( t h y 戚n e ) ， 鸟嘌呤g ( g u a n i n e ) 和胞嘧啶c ( e y t o s i n e ) 而按照一定顺序排列的a ，c ，g ，_ r 的 d n a 以痔捌的形势表达出来，也艇通常所指的d n a 序列的一级结构，我们镶称为 d n a ，掌嬲。 图1 1 ：1 9 5 3 年w 醯s o n 和e r 触发现由a ，c ，g ，t 形成的d n a 双螺旋结构 1 3 多黧序列比对的基本记号 为讨论d n a 序列的结构问题，我们记 a = 江1 ，如，t ，) ，g = e l ，现，- ，岛。)( 1 1 ) 分另是在z 4 与z 5 集合中取德的序列，其中，n 。分舅楚序列a ，g 的长度，黼 z 4 = o ，c ，p ，t ) = o ，l ，2 ，3 ) ，z 5 一 n ，c ，9 ，t ， ) = o ，l ，2 ，3 ，4 ) 分别表示碱撼与碱 基+ 虚拟擒入的符号的集合，这嬲称”o 或4 为虚拟符姆以下记 甄= l ，2 ， ，怒= l ，2 ，瓤 ，( 1 2 ) 第一章引言4 分别是序列a ，g 的下标或位点集合，那么 0 ，c 中的元分别是序列a ，g 的下标或 位点例如如下序列a ，g ： 序列a g 酷c t c t c c g 武t a g a c c a g a t t t g a g c c t g g c a g c 序列g ： g 或c t c t c c c c t t t a g a c c a g a t t t g 地g c c t g a t g c c 我们记为： a = ( 0 1 ，0 2 ，n 3 ，o n 。) ，o i 历， g = ( c l ，c 2 ，c 3 ，。c n 。) ，c t 磊 n 。和n 。是a 和g 的序列长度这里我们称a 是比对前的原始序列，e 是比 对后的序列( a l i g n e ds e q u e n c e ) 我们记 = a 。，s = l ，2 ，m ，( 1 3 ) 为一多重序列组，m 是序列的个数( 重数) ，其中a 。，s = 1 ，2 ，m 为单个序列 a s = ( 山函，o s ，n 。) ，o ”毛，( 1 4 ) 其中是序列的长度一4 是直接观测的原始序列的数据，也称观测数据 我们记 c = g ，= 1 ，2 ，m ) ，( 1 5 ) 为比对后的多重序列组( m l l l t i p l ea l i g i l m e n ts e q u e n c e ) ，也称为m s a ，其中m 是序 列的个数( 重数) ，g t ，= l ，2 ，m 为单个比对后的序列 g = ( o “，o 蚰，n s ，n 。) ，n 。，j 磊， ( 1 6 ) 其中n 。是序列的长度，并且经过多重比对后各个序列的长度相同c 是经过多重 比对后的数据 例如，多重序列组a 如下： a 1 = 豇t a g a a g a a g c c a a c a a a g g a g a g a a c ， 第一章引言 a 2 = a c a g g a a g a a g c c a a c a g g a g a g a a c ， a 3 = c a g t t a g a g c c a a c a a a g g 艇赡a 比对后的序列组c 如下： g = 料戥t a g a a g a a g c c a a c a a a g g a g a g a a c ， 岛= a c a g 料+ g a a g a a g c c a a 料c a g g a g a g a a c ， g = c a g t t a g $ a g c c a a c a a a g g a g a g a a c 则上面c 称为a 的一个多重比对( m s a ) 1 4 序列的突变 5 一般来说，序列突变可以分为四种； 1 i 型突变：又称为变异突变，比如由c 突变为t ，a 突变为g 例如：a a a c 酵一 8 a a g 昏这里第四位的c 突变为g 2 i i 型突变：又称为换位突变，它是由核苷酸的片段交换位置而产生的例如： a a a c 醇一c 酗a a a 这里前面三个核苷酸和后面三个交换了位置 3 i i i 型突变：又称为插入突变，它是由核苷酸被插入其他的片段而产生的例如； a a a c 酤一a a a t t c 酵这里在第四位插入了7 r t 这个核苷酸片段 4 i v 型突变：又称为删除突变，它是由核苷酸序列的某个片段丢失而产生的例 如：a a a c 酵一a a a 蓟这里原来序列里面第四位的e 被删除 i i i 型突变和i v 突变也被通称为插入删除型突变般来说，多重比对算法对 于i ，i i i ，型突变比较容易识别，也就是比较容易得到比较好的比对结果而对于换 位突变，又称为重组( r e a r r a n g e m e n t ) 突变，尤其是在染色体上的基因组范围内的大 规模重组( m a i v eg e n o h l er e a r r a n g e m e n t ) 的比对，是多重比对算法的个难点 1 5 多重比对的s p 打分矩阵 在多重比对的结果评估中，各种数学上的指标被用来对多重比对的结果做效果 的估计目前最常用的就是s p 打分矩阵其定义如下： 第一章引言 ( c ) 6 ( 1 7 ) 这里c 是一个多重比对，其中g ，g c 并且他们在磊中取值，n 是g ，a 的长度，。”g ，岛，g m 是c 包含序列的条数d ( 龟，j ，c 。，j ) 是序列间的打分 函数： ac气，岛。=主：2岛：f妻曼i岛j c - s ， 1 6 多重比对的优化准则 序列的多重比对是对( 1 3 ) 中的序列插入一些虚拟符号”+ ”使得得到的c 的得 分最大，或者罚分最小重序列比对有多种优化准贝i j ，如s p 一准则，s h a n n o n 熵准 则等，这些准则在不同意义下多重序列比对所达到的”对齐”程度，但它们存在的 个共同缺点是无法判定多重序列比对最优解的位置依据常见的多种多重序列比 对算法，我们提出多重序列比对的相似比准则这个相似比准则有三个优化指标， 即r = ( r 1 ，见，凰) ，它们的定义如下： ( 1 ) r 1 是实际比对的数据率如果记n 是多重序列4 的原始长度，n 是实际 采用的多重序列爿的长度，那么蜀= 等 在许多文献中，在作的多重序列比对时，实际比对的序列长度往往要比原始序列 长度要短，这就回大大降低多重序列比对的难度，因此r l 是多重序列比对的个重 要指标 ( 2 ) 岛是多重序列比对插入符号的比例记为多重序列4 的碱基总数， c = q ，岛，c k ) 是a 的多重比对序列，记c 中插入”一”的总数为7 ，那么等 是多重序列比对插入符号的比例，那么定义疡= 1 一等在多重序列比对中，要求 插入符号的比例等应尽可能少，因此要求恐尽可能大 ( 3 ) 凡是多重序列比对的相似比如w 是( 3 ) 式中所定义的两两比对的得分 啪 。州一 。闰 第一章引言 7 矩阵，w = ( ：t ) “：而是多重比对序列c 中各序列的相似率，也就是”幺是g 与g 序列的相似率，那么定义 耻高。聂。甓 ( 1 9 ) 对相似比指标矗显然有r 1 ，r 2 ，忍sl 成立，在多重序列比对是则要求r 1 ，忌，兄3 这三个指标尽可能大如果恐十分接近l ，或大于1 ，那么我们可以认为这个多重 序列比对实现或接近最优比对 1 7 主要结果 我们针对基因序列的多重序列的比对问题，根据序列突变与比对的”模代数” 结构理论，并在序列两两比对的基础上，应用系统聚类等方法，给出了同源多重序列 的超级快速比对算法( 简称s m a ) ，同时我们提出了相似比指标，来衡量多重比对 的性能该算法可适用于大规模( 如m 5 0 0 ，n l okb p ) 的同源多重序列的比 对计算s m a 的程序和测试数据我们已在网上公布，并提供比对计算服务我们分 别对8 3 3 0 kb p 的s a r s 序列与7 0 6 1 0 kb p 的h r u l 序列进行比对，主频 3 ogp c 机上完成多重比对的计算时间分别为2 1 分钟与3 4 小时s m a 算法在 速度上对h m m e r 有明显优势，而且根据相似比等各项优化指标测试，结果不差于 h m m e r 第二章基于模结构的多重比对快速算法s m a 沈世镒教授创造性的引进模代数理论来研究序列突变和比对的问题，并提出的 序列的模结构理论f l 】 为多重比对问题开创了一个新的研究方法在沈教授的指导 下，我们运用这个理论，并在序列两两比对的基础上，应用系统聚类等方法，给出了 同源多重序列的多重比对快速算法s m a ( s u p e rm 、l l t i p l ea l i g n m e n t ) 在动态规划和 逐步迭代法的基础上，该算法的计算复杂度为d ( m n 2 ) ，在s p a 1 8 1 算法基础上，该 算法的计算复杂度为0 ( m 2 佗) ，可适用于大规模( 如m 5 0 0 ，n 1 0k b ) 的同源多 重序列的比对计算，在比对规模、计算速度与优化指标等性质明显高与现有的多重比 对算法与软件包同时我们已在网上公布，并提供比对计算服务我们分别对8 3 3 0 k b 的s 址玛序列与7 0 6 1 0k b 的h i v - l 序列进行比对，在主频1 2gp c 机上 完成多重比对的计算时间分别为2 0 分钟与4 小时，其中对s a r s 序列实现了最优 比对，而对h i v 1 的比对结果也接近最优解，由s m a 计算的h i v _ l 结果的各项优 化指标明显高于h 数据库所公布的多重比对结果 2 1 序列的模结构基本记号和定义 我们记下面的比对后的序列： g = ( 岛，岛。，岛。) ， c ；= ( 龟，岛。，岛。) g 。和q 是第s 条和第条多重比对的序列下面我们引进模结构来描述序列 的比对： f 乙= ( t 。，l 。，) ，k = 1 ，2 ，k 。) 这里( t 。，k ，z 。，k ) 和k 表示；在比对后的序列z 的第t 。，k 位置有f 础个插入，并 且序列z 共有k 插入区域 比如对于如下原始序列： a 1 = g t t r 培a a g a a g c c a a c a 嘲g a g a g a 8 c ， a 2 = a c a g g a 唔a a g c c a a c a g g a g a g a 舭， a 3 = c a 卧t a g a g c c a a c 雒嶝a g a g a a c 8 第二章基于模结构的多重比对快速算法s m a 9 经过多重比对后的序列如下； a = 木+ + 酗t a g a a g a a g c c a a c a a a g g a g a g a a c ， q = a c a g ”+ g a a g 蛤g c c a a “c a g g a g a g a a c ， 岛= 4 c a 酵t a g 料料a g c c a a c a a a g g a g a g a a c 我们可以用如下数据结构来表示： 日c 。= ( o ，3 ) ) ， 日凸= ( 4 ，3 ) ，( 1 5 ，2 ) ) ， 凰= ( o ，1 ) ，( 7 ，4 ) ) 我们称这如。，比t ，日凸序列比对模结构，s a m s ( s e q u e n c ea l i g n m e n tm o d _ u l u 8s t r u c t u r e ) ，简称模结构 2 2序列扩张的模结构与模代数理论 2 2 1序列扩张的模结构 为讨论d n a 或r n a 序列的结构问题，我们记 4 = ( 0 1 ，0 2 ，- 一，n 。) ，g = ( c 1 ，c 2 ，一，。)( 2 1 ) 分别是在z 4 与z 5 集合中取值的序列，其中n 。，n 。分别是序列a ，c 的长度，而 z 4 = 口，c ，g ，t ) = o ，1 ，2 ，3 ，z 5 = o ，c ，g ，t ， ) = o ，1 ，2 ，3 ，4 ) 分别表示碱基与碱 基+ 虚拟插入的符号集合，这里称”+ ”或4 为虚拟符号以下记 虬= l ，2 ，) ，肌= 1 ，2 ，咐，( 2 2 ) 分别是序列a ，g 的下标或位点集合，那么蜕，c 中的元分别是序列a ，g 的下标或 位点 定义3 1称序列g 是序列a 的扩张，如果在序列g 中删除取值为4 的分 量，那么e 就变成a 这时称序列a 是序列c 的压缩 如果序列c 是序列a 的扩张，那么我们记 日= ( t k ，k ) ，自= 1 ，2 ，。) ，g = = ( ，y l ，7 2 ，h 。) ，( 2 3 ) 第二章基于模结构的多重比对快速算法s m a l o 是序列g 关于序列a 扩张的二种”模结构”表示，其中： ( 1 ) ( 砝，靠) 表示序列a 在位点缸后插入个长度为如虚拟插入的符号向量 而k 是从a 到g 的插入的符号向量的总数，这时有o e 1 2 o ) 为序列a 的插入位点的集合，那么i k 就是l 中的点，且靠= m 。因 此由g 唯一确定日反之，如果日给定，那么厶= i 1 ，t 2 ，h 确定，而 m ：靠，如果= 缸厶，这时由日唯一确定g lo ， 否则 定理3 1 ( 1 ) 如果序列a 与它的扩张模结构日或g 给定，那么a 的扩张序列c 唯 一确定，这时n 。= n 。+ 皂。靠，且吼在序列g 中的位置是； 五= i + 厶，其中 厶= “矧如 ( 2 ) 定理3 1 一( 1 ) 的逆命题成立，这就是：如果g 是序列a 的扩张，那么模结 构口或g 唯一确定 对此定理的证明我们作以下概要说明如果序列g 是a 的扩展，那么它必有分 解式： = d 1 ，如，一，靠。) ，( 2 4 ) 其中= 2 自。+ 1 ，且0sj 1 j 2 j 。一1 1 4 这里如，1 ，如，2 ，如，l ，如。2 分别是区间以，2 ，以，4 ，矗2 ，如，t 的长度，它们分别是序列e ，d 关 于a ，b 的插入片段 2 2 5序列突变的模结构 序列突变有四种类型，这就是i - 型突变：口一9 ； i l 型突变： ( o c c ) ( 9 讹) 一 ( g “) ( o c c ) ；i i l 型突变：。o 一8 ( u ) n ；与i v - 型突变：n ( “) 口一n 我们称 i ，i i - 型突变为非位移突变，而称i i i ，一型突变为位移型突变对位移突变我们可采 用突变的模结构表示如果 a = ( 。1 ，0 2 ，一，o 。) ，b = ( 6 l ，6 2 ，一，b ) ，m ，z j ， ( 2 1 5 ) 第二章基于模结构的多重比对快速算法s 似 1 5 是二个核苷酸序列，如果b 是a 的位移突变序列，那么它们的位移突变的模结构可 表示为： t = ( 蟊，玖) ，= l ，2 ，一，女。，( 2 1 6 ) 其中如是第一次位移突变的起始位点，0sf 1 t 2 2 ，3 ，那么系统聚类树为 1 7 r ： ( 1 ) ，( 2 ) ：号( 4 ) ，( 3 ) ，( 4 ) 些字( 5 ) ， 其中( 1 ) ，( 2 ) ，( 3 ) 为由a 1 ，a 2 ，a 3 所确定的原始聚类点，( 4 ) ，( 5 ) 为由系统聚类法所产 生的新聚类点，( 3 ) ，( 4 ) ! 弩( 5 ) 未! 示在( 3 ) ，( 4 ) 点的聚类中，( 1 ) ，( 3 ) 点的距离撮：近 ( 或得分最多) 由系统聚类树t 就可得到最大得分树死为：( 2 ) 一( 1 ) 一( 3 ) ，其中( 2 ) 一 ( 1 ) 表示由( 1 ) ，( 2 ) 二点所组成的无向弧 3 由最大得分树蜀与模结构矩阵“构造三重比对的模结构为： 日1 = 研，2 v 凰，3 醌= 凰，1 v ( 日l 一风，2 ) ，日3 = 风，l v ( 日1 风，3 ) 4 由三重比对的模结构宙= ( 皿，也，风) 构造三重序列4 = 1 ，a 2 ，a 3 ) 的比 对序列c = ( 白，q ，g ) ，其中g 为凡的扩张序列，具有模结构为风，s = 1 ，2 ，3 由此实现三重序列4 的多重比对对于m 3 的情形可作类似推广，但计算程 序要复杂得多 2 4基于逐步迭代法的多重比对算法s m a 基于逐步迭代法的多重比对s m a 比基于系统聚类树的多重比对算法s m a 的计 算步骤要简单，在保证了比较好的精确度的同时，也增加了计算速度 s t e p 1 ：与上文相似，对给定的有m 条序列的序列组a ，设j 是平均长度先 计算a 的两两比对相似率矩阵s 这里我们选择s p a 算法 第二章基于模结构的多重比对快速算法s m a1 8 s t e p _ 2 ：搜索相似率矩阵s ，首先选择最相似的两条序列，设为： a ，和a 。然 后，对a a n da 2 进行两两比对，这里我们选择动态规划作为比对算法我们记比 对后的序列a l 和a 2 为a 和q ，并且由q 和q 组成个新的比对的序列组g 然后将g 和q 转换到s a m s ，并记为三b 和b s t e p 3 ：继续搜索相似率矩阵s ，在a 选择条序列，设为ra ，与中任意 一条序列，设为a ，最相似然后，对a 和a 使用动态规划进行两两比对，结果 记为g 和g ，然后将g 和0 转换到s a m s ，并记为日q 和，使用模结构运 算将三b 加入三b ： e 0 r 洲y 日仉，伉c ，七i ， d o ：矗诜= 日吼u 点k 。 最后将g 抽入c s t e p - 4 ：重复s t 印3 直到a 中所以的序列都被比对后加入到e 中 s t e p - 5 ：根据模结构( s a m s ) 日c ，将a 中的序列转换成多重比对( m s a ) 2 5使用s m a 进行大规模病毒基因组比对 2 5 1 测试数据集 为了分析s m a 的计算能力和运算效果，我们使用了2 个病毒基因组数据集来 对其进行测试： 1 我们从g e n b a l l k 【2 l j 选取的8 3 条高质量的s a r s 病毒基因组序列平均长度 为3 0 b 左右由于这8 3 条基因组相似率比较高，因此比较适合对其多重比对 做进一步的分析，比如下文将要讨论的基因组突变的网络结构分析 2 我们从h i v d a t a b a 8 ef 2 2 】中选取了7 0 6 条h i v l 的病毒基因组序列其中分歧 度最大的两条序列的相似率为5 7 4 序列的平均长度是l o kb 我们把测试的 重点放在h i v - l 的数据库上，原因是，目前基本上很少有多重比度算法或者工 具能够处理7 0 6 l o b 的比较高的分歧度的数据 2 5 2用于比较的对比工具 我们希望尽可能多的选择用于比较参照的其他算法和工具，我们选择了当前流 第二章基于模结构的多重比对快速算法s m a 1 9 行的和最新的算法和工具，但是经过运行发现，只有h m m 腿能够参与h i v _ l 数据 库的测试由于下面的各种原因其他算法和工具都没有用于测试； 1 c l u s t a l 、运行速度太慢，在h i v 1 数据集上运行了一周还没有结果 2 t c o 压如和p o a 处理不了超过1 0 髟b 的序列 3 m u l a l i n ，m 1 l l t i - l a g a n ，和m a d 不能处理超过5 0 0 条的序列 4 c h a o s d i a l i g n 和m u 8 c l e 在处理序列数据时出现内存异常的错误 显然这些问题的产生主要是因为这些软件不能处理象h 一l 成百上千条1 0 b 这样规模的病毒基因组数据所有测试我们使用了同一台p 43 o gp c 所有的运算 时也基于算法的默认设置 2 5 3 测试结果 根据上文提到的多重比对最优化标准，s m a 和h m m e r 的运算结果和性能在 表3 1 中表示 表3 1 多重s a r s 与m v 1 序列比对的指标分析表 序列算法比对规模 s p 计算时间主频相似比插入率 名称 名称 m nb p 打分3 0 gp c 机 h m m e r 8 3 3 0 k 1 0 0 1 0 8 3 8 1 分钟 9 9 9 3 04 7 s a r s系统s d a8 3 3 0 k1 0 0 1 0 8 3 0 分钟 9 9 ，9 9 0 ，5 l 渐进s m a 8 3 3 0 k 1 0 0 1