进化树的研究
1 引言 生物信息学是生物技术的核心,是在分子生物学和信息科学共同发展的基础上产生的一门交叉学科,包含对生物数据的获取、处理、存储、分发、分析、挖掘等方面的研究内容。生物信息学的研究对于最终改善人类自身生活质量,解决健康问题等也有重大的作用。随着分子生物学的不断发展,人们惊奇地发现DNA 的双螺旋结构中蕴涵着生命的密码,四种核苷酸的排列、变化包含着许多遗传、进化信息。人类基因组计划以来,有关核酸(或蛋白质)序列和结构的数据成指数增长,而面对如此复杂的数据,计算机在此方面的应用必不可少。因此,生物信息学研究的目的就在于,人们通过数学、计算机科学等各种工具,可以阐明和理解大量数据包含的生物学意义。
由于深度测序和基因芯片技术的不断完善和发展,表达谱、转录组、基因组等数据不断增长。到目前为止,已被测序的昆虫基因至少有10个,被报道的转录组数据也有30多个。生物信息学在昆虫学研究中的应用价值随着昆虫学研究的不断深入和昆虫生物数据的大量积累越来越明显。大量医学昆虫、经济昆虫和农业昆虫的基因组在模式昆虫果蝇的基因组测序成功之后也相继被测序。昆虫种类繁多、进化关系复杂、个体发育系统多样对于生物的多样性组成也占有举足轻重的地位。此外,昆虫与人类的日常生活和生产亦有密切的关系。例如,家蚕、蜜蜂等经济类益虫能够为人类提供日常生产资料和生活资源,害虫能给人类带来巨大的损失。对昆虫基因组进行深入研究不仅能为传统昆虫学科的发展提供崭新的机遇,而且对深入了解昆虫的多样性及其生物学特征与本质具有重大意义。 所有生物都可以追溯到共同的祖先,生物的产生和分化就像树一样的生长,分叉,因此以树的形式来表示生物间的进化关系是非常合理的。根据各类生物间的亲缘关系的远近,把生物安置在树状图表上,简明地表示生物的进化历程和亲缘关系的树状结构就是进化树。在进化树上每个叶子结点代表一个物种,每一条边都被赋予一个适当的权值的话,两个物种之间的差异程度就可以用两个叶子结点间的最短距离来表示。
2 生物信息学
2.1 生物信息学的诞生 1953年Watson 和Crick 发现了DNA 的双螺旋结构,开辟了现代分子生物学的新纪元。遗传中心法则的提出大大推动了分子生物学的发展。随后,限制性内切酶的发现、重组DNA 克隆技术的实现是新一代深度测序的基础,同时也是海量生物数据产生的重要推动力[1]。1961年,计算科学首次被应用于基因和蛋白质的进化分析中。1990年人类基因组计划正式启动。此后更多的生物基因组测序计划开始进行,生物数据不断积累。然而缓慢的信息挖掘速度与成指数级增长的生物数据形成了巨大差距,这种差距进一步推动了生物信息学的发展。
2.2 生物信息学的概念
生物信息学以数学、信息学和计算机科学为主要手段,运用计算机软、硬件和计算机网络,对核酸、蛋白质等生物大分子数据库进行研究;对数据库中的原始数据进行存储、管理、注释、加工等, 形成具有明确生物意义的生物信息;通过信息查询、搜索、分析、比较, 获取基因编码、核酸和蛋白质结构功能及其相互关系等知识;在大量信息和知识基础上, 探索生命起源、进化以及细胞、器官和个体的发生、发育、病变、衰亡等生命科学问题, 搞清它们的基本规律和时空联系, 建立“生物学周期表”。广义地说, 生物信息学是一门使用数学和信息学的观点、理论和方法研究生命现象、组织和分析迅速增长的生物信息数据的学科。主要研究的是遗传物质的载体DNA 及其编码的大分子物质量,以计算机为主要工具,对各种学科交叉的生物信息学进行研究, 找其规律性, 发展出适合它的各种软件,对不断增长的DNA (蛋白质)序列和结构进行收集、整理、提取、加工、分析等。狭义地说, 生物信息学是以计算机科学和数学为主要工具,对生物大分子进行信息的获取、加工、存储、分类、检索等。通过分析逐步认识生命的起源、进化、遗传和发育的本质, 逐步破译隐藏在DNA 序列中的遗传语言, 在分子基础上解释人体生理和病理过程, 为人类疾病的诊断、预防和治疗提供最合理、有效的方法或途径是生物信息学研究的最终目的[2]。
2.3生物信息学在昆虫研究中的前景和展望
由于测序技术不断发展,10年左右可能会形成昆虫测序的高峰期。但昆虫种类多、进化关系复杂、基因杂合度大,基因组测序工作和分析工作仍有很大困
难。昆虫基因数据虽然获得了大量积累,并且正成指数级增加,但相对自然界丰富的物种来说,目前已获得的遗传信息仍然是非常渺小的。生物信息学在昆虫学领域的应用仍处于起步阶段,有很大的发展空间。昆虫基因数据的大量积累也将大大推动对资源昆虫的利用和害虫的控制,甚至有可能开辟资源昆虫利用的新领域,催生出全新的害虫控制技术。
3 基于基因组的昆虫学研究应用
随着越来越多的昆虫正在进行或已经完成基因组测序,现阶段昆虫学的研究已经步入了基因组时代。基因组时代的昆虫学研究方法相比于传统的昆虫学研究被赋予新的时代特征,其中基因组时代昆虫学研究应用包括:
1) 在个体生物学研究中的应用
2) 在多物种间及种群研究中的应用
3) 在系统生物学研究中的应用
4 基因组时代昆虫学研究所面临的挑战
4.1 未来发展趋势
预计未来5到10年,随着大量昆虫的基因组测序工作的完成,昆虫学研究将迈向基因组时代的全盛期。昆虫学的研究也将步入后基因组时代,即功能基因组时代。我们要逐渐将个体生物学、进化学向多物种间或种群内不同个体间的基因组学过渡;将传统的生态学、生理学、行为学以及分子生物学相关研究逐步向更为全面的系统生物学过渡。
4.2 存在的问题与对策
综合国内外现有的有关昆虫功能基因组学及结构基因组学的研究,发现虽然基因组时代的昆虫学研究进展迅速,但也存在以下4个主要问题和不足[3]。
1) 全基因组测序的昆虫样本准备起来较为困难, 导致难以获得纯和的DNA 。
2) 测序方法选择错误:选择第2代测序方法产生的片段较短,后续拼接组装难度较大,易造成测序缺口。
3) 多重视结构基因组草图绘制时,遗传图谱构建和遗传标记的筛选被忽视,缺乏绘制精细物理图谱的有效手段。
4) 大多非模式昆虫缺少遗传群体和突变体库。
针对基因组时代昆虫学研究存在的4大困难,提出以下对策:
1) 建立并优化供试昆虫多代自交、继代饲养体系。
2) 采用第2代测序技术与传统建立人工染色体库相结合的方法。
3) 在开展昆虫全基因组测序工作之前,先初步完成对供试昆虫遗传图谱构建和遗传分子标记筛选。
4) 根据昆虫自身特征,对已完成的全基因组测序工作的昆虫基因的编码和非编码基因序列选择合适的功能基因组学研究手段开展研究。
5 进化树
5.1 进化树的定义及形式
利用DNA 序列进行发育分析, 推断并评价在分子水平上物种的进化关系,并用分支图的形式表现出来,这种图就是进化树。进化树有多叉树,但通常情况是二叉树。它分为有根(rooted )树和无根(unrooted )树两种。有根树反映了树上物种的时间顺序,而无根树只反映分类单元之间的距离不涉及谁是祖先的问题。也就是说,有根树的根节点为全部分类单元最近的共同祖先,它反映了分类单元间的进化关系,而无根树仅反映出分类单元间的分类关系[4]。如图5.1,表示了4个物种(A 、B 、C 、D )的2种有根树和1种无根树形式。
图5.1 4个物种(A 、B 、C 、D )的2种有根树和1种无根树形式。
5.2 构建进化树的目的和作用
构建进化树的目的是重塑物种(分类群)之间的进化关系,并以进化树的形
式描述和展现:树的叶子结点代表某个序列的来源物种,数的拓扑结构表现了各物种间亲缘关系的远近,树的分枝长度刻画了进化树距离的大小。根据进化树的拓扑结构可以研究生物蛋白质分子以及病毒、细菌以至大型哺乳动物等各种有机体之间的生物进化关系。通过进化树的支长长度可以近似估计分类群的分化时间。进化树的理论意义是在于它有助于了解物种的进化历史,为生物学中物种的分类提供可靠依据;实际应用价值在于它在预测DNA 分子的高级结构、蛋白质、基因的表达过程、辅助药物设计、结构多序列比对等方面均有重要作用。
5.3 构建进化树的主要过程
构建进化树的主要过程包括:
1) 获取同源序列数据
2) 确定进化模型
3) 进行多序列比对
4) 根据比对结果提取信息
5) 选择建树算法与参数构建进化树
6 同源数据的获取
所谓同源序列简单的说:是指从某一共同祖先经趋异进化而形成的不同序列。NCBI(National Center for Biotechnology Information)美国国立生物技术信息中心,负责管理GeneBank 。GeneBank 是美国国立卫生研究院维护的基因序列数据库,汇集并注释了所有公开的核苷酸序列。GeneBank 是最常用的数据库,同时还有DDBJ (DNA Data Bank of Japan )日本DNA 数据库, 以及EMBL (European Molecular Biology Laboratory )欧洲生物信息研究所的欧洲分子生物学实验室核苷酸数据库。这三个中心都可以独立地接受数据提交,而且3个中心之间逐日交换信息,并制成相同的充分详细的数据库向公众开放。因此,同源序列可以通过网络在NCBI 等数据库中获取。根据NCBI 数据库记载,已经完成拼接或正在绘制的昆虫基因草图已有50多个。
7 构建进化树的方法及特点
系统进化树的构建从方法上可以分为2类:距离法和离散特征法。
距离法首先构造一个距离矩阵,矩阵上的元素代表每两个生物之间的距离,利用不同的聚类方法得到系统进化树。典型的基于距离的方法有UPGMA 、NJ 、和Fitch –Margoliash 。UPGMA 只适用于进化树中各个分支进化速率相同的情况;NJ 法能够较快地生成单一亲缘的进化树,但构建的进化树不够准确;Fitch-Margoliash 法比UPGMA 准确。距离法是基于距离的聚类算法,运算量小,但不能确定进化分支时间。
离散特征法根据最优原则不同分为最大似然法(ML )和最大简约法(MP )。ML 和MP 是基于特征或符号的构建进化树的方法,事先不需要规定距离测度计算距离矩阵,而是直接通过各分类群序列的碱基或氨基酸顺序来构建进化树。其计算量比距离矩阵法大,对于大种系发育分析来说距离矩阵法比较常用。
距离法、最大简约法、最大似然法在分子序列间的分歧度不高且序列较多时往往构建出具有相似拓扑结构的进化树。但是还是存在着存在拓扑差异。因此有必要了解进化树的构建特点,以便对特定的序列选择合适的构建进化树的方法。 UPGMA (Unweighted Pair Group Method with Arithmetic mean)假设存在一个分子钟,即在进化过程中所有的核苷酸或氨基酸有相同的变异率。在不同谱系间进化速率差异较大、有同源序列的平行进化或进化树的状态空间大示,一般不宜采用此方法。
NJ 法经常被使用,它计算速度最快,并且构建的进化树相对准确。但NJ 法建树时序列上的所有位点都被同等对待,而且所分析序列的进化距离不能太大。NJ 法每迭代一次只搜索最近邻居的配对,而不考虑其他可能的配对,最终只生成单一的最优树而遗漏一些拓扑结构更合理的次优树[5]。
FM 法是距离法建树中一种相对较好的方法。它从所有可能的进化树中挑选树长最短的作为最优树。不过,必须注意的是距离法在将原始数据转换成距离矩阵时难免会丢失一些进化信息。
不同的算法有不同的适用目标。一般来说,最大简约性法适用于符合以下条件的序列:1)所要比较的序列碱基差别较小。2)对于序列上的每一个碱基有近似相等的变异率。3)没有过多的转换/颠换倾向。4)所检验的序列碱基数目较多(大于几千个碱基)。最大简约法是一种不依赖任何进化模型的无噪声统计方法[6],能快速分析出大量序列间的系统发生关系,所构建树中的短分支更接近真
实。但所有重建祖先序列中的最小突变数决定了简约树上的数值,而突变是否按照事先约定的核苷酸最少替代途径进行目前还尚未确定。单一的突变图谱可能会得出似是而非的结论。并且,所有分支的突变数不可能完全相同,由于没有考虑核苷酸的突变过程,使得长分支末端的序列由于趋同进化而显示较好的相似性趋同现象,这违背了简约法则,导致对“长枝吸引”的敏感[7]。因此,当序列单位位点上核苷酸替代数较大时,简约法极可能得出拓扑结构错误的进化树[8]。 最大似然法的计算过程极其耗时。如果分析的序列较多的话,有可能要花上几天的时间。最大似然法考虑了所有可能的突变路径,能将数据的系统发生信息完全利用。然而,最大似然法在很大程度上依赖于对核苷酸替代模型的选择[9]。最大似然法没有评估拓扑结构的优劣,而是假定分支长度估计最精确的拓扑结构为最优树。最大似然法的目的是寻找参数空间中的最高点,因而最终结果仅取决于某个点的数值。实际上,系统发生所关心的是树的拓扑结构,分支长度反而成为干扰参数,忽略分支长度应该更合理些[10]。
8 多序列比对的软件
对于一个完整的进化树分析需要:1)要对所分析的多序列目标进行比对(alignment )。2) 建立一个进化树。3) 对进化树进行评估。 多序列比对用来区分一组序列之间的差异,但主要用来描述一组序列的相似关系,以便对一个基因家族的特征有一个简明扼要的了解。做alignment 的软件很多,经常使用的有clustalx 和clustalw ,前者是在window 下的,后者是在dos 下的。
CLUSTAL 比对是将所有序列两两比较,得出两两间差异值(最粗的距离),再根据序列间的差异把差异越小的序列放在一起构建一个分类树。最终操作是以这个分类树作为引导树,从各个相似序列的组作为起点做多重联配,直到所有序列被联配上。下面介绍用CLUSTALX 软件对已知DNA 序列做多序列比对的步骤:
1) 以FASTA 格式准备好DNA 序列*.txt文件(例sequence.txt )如下所示:
2) 进入CLUSTALX 程序,点击File 进入Load Sequence, 打开上述txt 文件, 值得注意的是txt 文件的保存路径中不能出现中文,且序列的名称不宜过长否则可能导致无法打开。
3) 点击Alignment-Alignment parameters,在其子目录下选择Multiple Alignment parameters ,可以设置比对参数。设置好后点OK 确认。点击Alignment –Output Format Option可以选择输出格式,如PHYLIPH 或FASTA 格式。点击Alignment-Do Complete Alignment,进行完全比对。这时输出两个Clustal 默认的两个格式:比对文件sequence.aln 和向导树文件sequence.dnd 。其中aln 格式文件是默认输出,可以转换成各种格式,而且很多软件都支持这种格式。dnd 格式文件是根据两两序列相似值构造的一个指导后面多重连配的启发树,它不能做进化分析。因为进化分析要考虑的所有同源位点的综合效应,所以用aln 格式
本 科 毕 业 论 文
文件专门做进化分析。 第 9 页 共 29 页
9 序列分析软件
通常用DNAMAN 软件进行分析预测结果。DNAMAN 是一种常用的核酸序列分析软件。由于它功能强大,使用方便,已成为一种普遍使用的DNA 序列分析工具。使用DNAMAN 可以画质粒图谱,两个或多个序列的相似性比较,测序后两端序列拼接成一个完整的序列。还有引物设计、序列的限制性内切酶图谱、寻找序列的二级结构等功能。DNAMAN 的操作界面有20个Channel ,点击Channel 工具条上相应的数字,即可击活相应的Channel 。每个Channel 可以装入一个序列。将要分析的序列(DNA 序列或氨基酸序列)放入Channel 中可以节约存取序列时间,加快分析速度。而且DNAMAN 提供自动载入功能,用户只需激活某个Channel ,然后打开一个序列文件,则打开的序列自动载入被激活的Channel 中。下面仅介绍如何用DNAMAN 进行序列分析。
双击打开DNAMAN, 激活某个Channel 后,通过sequence-load sequence选择路径如from sequence file 分别导入到不同的Channel 中。要比较两个序列,可以使用DNAMAN 提供的序列比对工具Dot Matrix Comparison(点矩阵比较)通过执行Sequence-Dot matrix comparison命令,设置比对序列,进行比较。进行序列同源性分析,通过执行Sequence-Alignment 可以选择Two Sequence Alignment 进行两个序列同源性分析,和Multiple Sequence Alignment进行多序列同源分析。下图为三个序列的同源分析结果:
本 科 毕 业 论 文 第 10 页 共 29 页
从结果可以看出序列间的相似度。并且可以点击右上角的Options 按钮,可以从弹出的对话框中选择不同的结果显示特性选项。点击Output, 出现一些可选项,可以通过这些选项,绘制同源关系图(例如Tree-homology tree 命令)。显示蛋白质二级结构(Protein Secondary Structure命令),绘制限制性酶切图(Restriction Analysis命令)等。
10 进化树的构建算法
10.1 UPGMA算法
10.11 UPGMA算法特点
不加权算术平均组对方法(Unweighted Pair Group Method using
Arithmetic mean ,UPGMA )是一种重要的进化树构建方法,源于数值分类学中表征图的构建。应用该方法的基本假设是各分类群的进化速率相同,因此从祖先节点分离出来的两个分类群到该祖先节点的分支长度相等。当多个种群之间距离相等时改变物种的输入顺序可能产生不同拓扑结构的进化树[11]。不加权算术平均组群方法(UMGMA )是对UPGMA 的改进算法,可以产生多叉树,和解决不唯一性问题。 10.12 UPGMA算法方法及基本思想
UPGMA 将距离矩阵中最小距离的两个分类群聚类,并计算新分类群与剩余分类群之间的距离得到新的距离矩阵,至此完成算法的一轮计算,重复此过程,最
终会得到一个以所有原始分类群为叶节点的进化树。其基本思想是反复把两个距离最近的种群合并,直到形成一个包含所有物种(通常用DNA 或蛋白质序列表示) 的种群。
10.13 UPGMA中的距离矩阵
按照一定的遗传模型,从已知的生物序列中推断出各个物种之间的进化历史,把任意两个序列间的进化历史转化成数字,就得到两两之间的进化距离,把所有的距离用矩阵的形式表示出来,就得到距离矩阵。根据距离矩阵将构建出一棵进化树。UPGMA 的初始矩阵可以如表10.1表示:
表10.1 初始序列距离矩阵
从该矩阵开始迭代计算新的距离矩阵,并完成分类群单元的聚合。其中D ij 表示第i 个物种和第j 个物种之间的距离。假设D 12为最小值,则先将物种1、2聚合成新的分类群单元u =(1,2) 。得到新的距离矩阵。如表10.2所示:
表10.2 迭代一次的距离矩阵
∑d [x , k
x ∈U =
U ∑∑x ∈U y ∈V
U . V
]
新的分类群与物种k 之间的距离为;d [U , k ]而两个新的分
d [x , y ]
类群单元U 、V 之间的距离可以用公式 d [U , V ]
=
给出。所有
物种合成一个分类群时,UPGMA 算法结束[12]。
例:用MEGA 软件对11种果蝇的adh DNA 序列数据用UPGMA 方法可以得到下图的进化树。
10.2 Fitch-Margoliash 法
Fitch-Margoliash 法首先估计树中各长度,如果仅有3个分类群,先假定分类群A 和B 之间的进化距离最小,然后与C 聚合。如果多于3个分类群,则将其他分类群合并为C , AB 的距离不变,A (或B )与C 的距离d AB (或d B C )为A (或B )与C 中所有分类群的距离算术平均值。
Fitch-Margoliash 的拓扑结构如图10.1,AB 间的距离为分枝x 和y 之和,即d AB =x +y ,同样有d BC =y +z,d AC =x+z,解出x ,y ,z 。将AB 合并成一个新的分类群,用类似的方法计算出新的x ,y ,z 值。依次可计算出所有分支长
n
n
度的估计值。所有可能的进化树的残差平方和为:R s =
∑∑
i =1
j =1
d ij -d ij
d ij 2
。其中n
为分类群个数d ij 分别为分类群i ,j 间的观察距离和估计距离,取R s 值最小的为最优树。
图10.1 Fitch-Margoliash 法示意图
10.3 NJ算法
本 科 毕 业 论 文
10.31 NJ算法基本思想
第 13 页 共 29 页
NJ 法从星形树出发,假定X 为共同祖先,如图10.2,可计算出星形树的分
n
支长度和:S 0=
∑
i =1
L iX =
1n -1
n
∑
i
d ij (式中i 和j 指不同物种,n 为所分析的
物种数目,下同) 。
图10.2 NJ的星形图
图10.3 分类群1和2相聚的树
定义由单一的内节点连接的分类群为一对邻居,假定分类群1、2为一对邻居,如图10.3,利用最小二乘法估计树的分支长度为S 12。其中
n
S 1
2=
L XY +L
X 1
+L
X 2
+∑L iY =
i =3
(d ∑2(n -2)
I =3
1
n
i +1
d i ) +2
12
d +
1
n
n -2
∑
d ij 。由于事先不
知道哪两个分类群为一对邻居,因而首先应分别计算任意两个分类群为第一对邻居时树的分支长度和,选取值最小时的两个分类群作为邻居,连接它们的内节点,然后在剩余分类群中选取任意两个分类群的最小分支长度和为第二对邻居,类似地找到第二个内节点。循环下去,直至找到所有的内节点。 10.32 NJ算的的步骤
NJ 算法构建进化树的步骤[13]分为: 1)计算任意序列i 和j 之间的距离,d ij =
⎛3⎫
ln ⎪, 其中q 4⎝4q -1⎭3
为两个序列中
相应位置上相同碱基的概率,即q =
两个比对序列中相同位置碱基相同的总数
max (两个序列的总长)
。
N
2)计算第i 个叶子节点(即第i 个序列)的净分歧度:r i =叶子节点的个数,d ik 为叶子节点i 和k 之间的距离。
∑d ik ,其中N
k =1
是
3)计算任意两节点i 和节点j 之间的速率校正距离M ij ,并挑选出最小的速率校正距离M ij :其中M ij =d ij -
r i +r j N -2
。
4) 定义一个新的节点t ,t 的左、右孩子分别是第i 个节点和第j 个节点。t 到
i ,j
的距离为d it =
d ij 2
+
r i -r j 2(N -2)
,d jt =d ij -d it 。
5) 从叶子节点集合L 中删除节点i 和j ,并在集合中添加t 节点,N 的值减去1。
6) 如果叶子节点集合L 中剩余的叶子节点个数大于2,则重复步骤2继续计算,直到叶子节点个数为2,即进化树完成建成。 10.33 NJ算法的问题与改进
NJ 算法的计算效率优于其他算法,但是当剩下最后3个节点时,无论3个节点间的距离是多少,它们的速率校正距离M ij 都是相同的。因此,如果选择错误,可能会产生误差较大的进化树。
NJ 的改进算法借鉴UPGMA 算法的思想即每次挑选具有最小d ij 的节点,把它们作为相邻的节点,同时吸收NJ 算法的精髓,采用校正距离,从而产生一个新的选择标准:R ij =αM ij +βd ij (其中α+β=1)。这样即考虑了M ij ,又解决了最后三个节点M ij 全部相同的问题,改进算法的效率大大提高了。
10.4 ML法建立系统进化树的模型和算法
10.41 ML法的基本思想
ML 法建立进化树,是通过引入一个2m 维的向量,该向量经转化后将似然函数变成一个二次多项式,再利用代数方法求近似解。 10.42 变量说明
设N ={'A ', ' T ', ' C ', ' G '}表示基因集, X ={ X 1,X 2 , ⋯, X m } 表示给定的
m 个生物体(系统进化树中的叶子节点) ,X i ={x i 1, x i 2,..., x in
},x ij ∈N ,
(i =1, 2,..., m ; j =1, 2,..., n ) , X i 表示第i 个生物的基因序列;Y 表示祖先生物体
(系统进化树中的内部节点) 的集合。生物进化的基本单元结构如图10.4表示。
A
、B 、C 的基因序列分别为A =(a 1, a 2,..., a n ), B =(b 1, b 2,... b n )
, C =(c 1, c 2,..., c n );
t 1, t 2为进化时间(图中节点上的字母代表生物体, 边上的字母代表进化时间) 。
p 1, p 2 为C 相对A 、B 的保守概率(基因不发生变异的概率) 。进化时间越长,
[14]
生物发生变异的概率越大
。
图10.4 生物进化示意图
p k
与t k k =(1,2) 之间的关系为:
⎛ C = G T ⎝A
A p (t ) s (t ) s (t ) s (t )
C s (t ) p (t ) s (t ) s (t )
G s (t ) s (t ) p (t ) s (t )
T ⎫s (t ) ⎪
⎪s (t ) ⎪s (t ) ⎪
⎪p (t ) ⎭
p (t )
-4αt
⎧) /4⎪p (t ) =(1+3e
, ⎨ -4αt
) /4⎪⎩s (t ) =(1-e
如果保守概率等于0.25, 表示下一代没有遗传;保守概率等于1,表示下一代没有变异。所以, 0.25
1
n
∑n
i =1
u i , 其中u i =⎨
⎧1, ifa i =b i ⎩0, ifa i ≠b i
。
10.43 确定父亲生物的ML 模型
根据生物进化关系, 在每个时间层面上, 最相似的两个生物一定是最近从某一生物(其父亲) 进化来的。所以, 构建系统进化树的过程就是从现有生物体中找出最相似的两个个体, 通过它们估计出其父亲的基因序列, 并将时间推移到其父亲的时代, 重复这一过程, 直至找出全部生物的共同祖先。对10.4 所示的进化单
元结构, 以保守概率p 1 和p 2 为参数, 通过极大似然估计得到序列C , 即求C , 使条件概率P (A , B /C ) 达到最大, 具体模型为:设似然函数L (C , p 1, p 2) 表示给定序列C 出现A 、B 的概率, 即:
n
L (C , p 1, p 2) =P (A , B |C , p 1, p 2) =
∏P (a , b
i
i =1
i
|c i , p 1, p 2) , 其中
⎧p 1, p 2, ifc i =a i =b i ; ⎪
⎪p 1(1-p 2) 3, ifc i =a i ≠b i ;
P (a i , b i |ci , p 1, p 2) =⎨(i =1, 2,..., n )。记
⎪(1-p 1) p 23, ifc i =b i ≠a i ; ⎪⎩(1-p 1)(1-p 2) 9, else
n
, p 2) =l (c i , p 1, p 2) =P (a i , b i |c i , p 1, p 2) ,则L (C , p 1
∏l (c , p , p
i
1
i =1
) 2
。使L (C , p 1, p 2) 达到
最大,得到:c i (p1, p 2) ={z |l(z,p 1, p 2) =max(l (x,p 1, p 2))},i =1, 2,... n 。由上式可知如果已知p 1, p 2,可通过序列A 、B 就可以确定其父亲C 的基因序列。
10.5 MP算法
最大简约法(MP)来自Occam 剪刀哲学原理:对一现象的解释,应选择最简单(假设最少) 的解释而不是更复杂的。对系统发生分析而言,该法则有两层含义[15]:一是用最少的进化事件(如突变) 去解释观察到的数据。二是在任何进化模型或机制下,假设尽可能少。具体来说最优树就是在解释整个进化过程时,全部位点的最小核苷酸替代数之和最小的树。假定有四条序列i :AGGGTAACTG ;j :ACGATTATTA ;k :ATAATTGTCT ;l :AATGTTGTCG 。可构建出3种无根树。如图10.5。
10.5 3种可能的四分类群无根树
将四条序列相对应的各位点核苷酸依次代人树的相应位置,遍历所有可能的进化路通径,分别计算出其中最小的核苷酸替代数:树A(0 3 2 2 0 1 1 1 l 3);树B(0 3 2 2 0 1 2 1 2 3);树c(0 3 2 1 0 1 2 l 2 3),树A 中全部位点的最小核苷酸替代数之和最小,为最优树。上述的计算仅位点4、7、9对构建MP 树提供了有用的信息,其他位点在3种拓扑结构中都有一致的最小核苷酸替代数,为简
化计算,可以只考虑要求选取的位点。MP 法利用每个可变位点的所有进化路径的分支(包括内部或外部) 的核苷酸平均替代数估算分支长度。当分类群较少时,可穷尽式搜索所有可能的树从而找到最优树;当分类群较多时,可以采用分支限界法和启发式搜索两种搜索方法,只搜索那些可能正确的树。
11 构建进化树的软件及实例
11.1 用MEGA 软构建进化树的步骤
1)序列文本
建树之前先将每个样品的序列都分别保存为txt 文本中,序列只包含序列字母(ATCG 或氨基酸简写字母)。
2)序列导入MEGA 软件中,界面如下:
然后导入需要构建进化树的序列:执行Align-Edit/Build Alignment 点击OK 。如果是DNA 序列,点击DNA ,如果是蛋白质序列,点击Protein 。
执行Edit-Insert Sequence From File 插入要建立进化树的序列。
3)序列比对分析
点击Alignment-Align by ClustalW 或直接点击W, 开始比对,比对完成后删除序列两端不能完全对齐的碱基。然后点击Data-Phylogentic Analysis 进行系统分析。关闭该窗口并保存。
4)进化树的构建
选择Analysis-Phylogeny ,然后选择需要采用的分析方法构建进化树,以NJ 树为例,Bootstrap 选择1000,点Computer ,开始计算,计算完毕,生成系统进化树。
11.2实例
在NCBI 上搜索小胸鳖甲Microdera puctipennis(Genebank 登录号:JF421286.1) 光滑鳖甲Anatolica punctipennis (Genebank 登录号: EF569673.1)、赤拟谷盗Tribolium castaneum (Genebank 登录号:XP974442)、
半目大蚕蛾 Antheraea yamamai (Genebank 登录号:AB179657.1)、二化螟(Chilo suppressalis ,Genebank 登录号:ADE05296.1),粘虫(Mythimna separate,Genebank 登录号:ABY55233.1),谢氏宽漠王(Mantichorula semenowi ,Genebank 登录号:ADB44081.1),果蝇(Drosophila auraria ,Genebank 登录号:CAA55168),拟南芥(Arabidopsis thaliana ,Genebank 登录号:CAA05574.1)等昆虫的HSP70氨基酸序列。通过MEGA 软件,采用Neighbor-Joining 法,重复1000次,构建Bootstrap 验证的系统进化树:
此方法结果表明,HSP70基因的保守性很强,小胸鳖甲与同科的光滑鳖甲亲缘关系最近,其次是赤拟谷盗,而与同科的谢氏宽漠王关系较远,推测结果可能与它们的地理环境不同有关。小胸鳖甲与光滑鳖甲的采集地点相同,生活环境也比较相似,而与谢氏宽漠王的区别较大[16]。此例也说明了环境对生物进化的影响。
11.3 用CLUSTALX 和PHYLIP 软件构建进化树
一. 用CLUSTALX 软件对已知DNA 序列做多序列比对。由于前面已有详细的比对方
法,这里只就比例,简单的说明。步骤如下:
1) 以FASTA 格式准备好DNA 序列*.txt文件(例sequence.txt )如下图所示:
2) 进入CLUSTALX 程序,点击File 进入Load Sequence, 打开上述txt 文件。
3) 点击Alignment-Alignment parameters,在其子目录下选择Multiple Alignment parameters,设置比对参数。设置好后点OK 确认。点击Alignment –Output Format Option 选择PHYLIPH 输出格式。点击Alignment-Do Complete Alignment ,进行完全比对。文件在PHYLIP 软件目录下以*.phy存在。用记事本方式打开文件部分序列如下:
本 科 毕 业 论 文 第 21 页 共 29 页
其中8表示有8个序列,1915表示每个序列有1915个碱基。
二. 用PHYLIP 软件推到进化树
1) 将PHYLIP 软件目录下刚建立的*.phy文件拷贝到exe 文件夹中。依次使用seqboot ,dnadist ,neighbor ,consense 。seqboot 为重复抽样,现已成为分析分子树置信区间最常用的方法。双击进入seqboot ,输入建立的*.phy回车。输入r ,回车,把bootstrap 改为适当的值,如果值太低,树不可靠。这里将bootstrap 改为1000。输入1000,回车。再输入y ,回车。输入random number seed (随机种子数,必须是奇数)。例如5,回车。可以看到计算过程。
关闭seqboot ,在exe 文件中出现一个outfile 文件。
本 科 毕 业 论 文
2) 用dnadist 计算核苷酸距离。
第 22 页 共 29 页 将outfile 文件重命名为infile 。双击打开dnadist ,输入m ,回车。再输入d ,回车。输入刚才设定的republicate 的数目,设定好参数后输入y 确认。程序开始运行,并在exe 中生成一个outfile 。将outfile 文件名改为infile ,为避免与原来infile 文件重复,将原先文件改为infile1。
3) 在exe 下选择通过距离构建进化树的算法。点击neighbor 程序,输入m ,回车。输入设定的republicate 值,回车。输入设定的奇数种子数,回车。输入y 确认后,程序开始运行。并在exe 文件夹中产生outfile 和outtree 两个结果输出。outtree 文件是一个树文件,可以用treeview 等软件打开。outfile 是一个分析结果的输出报告,包括了树和其他一些分析报告,可以用记事本直接打开。打开部分内容如下:
4) 将outtree 文件名改为intree ,点击drawtree 程序,输入font1文件名,作为参数,输Y 回车确认参数。程序开始运行,并出现Tree Preview图。
5) 点击DRAWGRAM 程序,输入font1文件名,作为参数。输Y 确认参数。程序开始运行,并出现Tree Preview图。
6) 将exe 文件夹中的outfile 文件名改为outfile1,以避免被新生成的outfile 文件覆盖。点击consense 程序。输入y 确认设置。exe 文件夹中新生成outfile 和outtree 。Outfile 文件用记事本打开,内容如下:
7) 将exe 文件夹中的intree 文件名改为intree1,将outtree 改intree 。点击drawtree 程序,输入font1文件名,作为参数。输Y 确认参数。程序开始运行,并出现Tree Preview图。
8) 点击DRAWGRAM 程序,输入font1文件名,作为参数。输Y 确认参数。程序开始运行,并出现Tree Preview图。
上述例子中在NCBI 中搜索了8种昆虫的CytP4 450 。并对CYP4家族中的某些基因进行了进化树的构建。可以看出果蝇(Drosophila melanogaster)与致倦库蚊(Culex quinquefasciatus )被聚类在同一进化分支中,表明它们的亲缘关系最近。而与美散白蚁(Reticulitermes flavipes)的亲缘关系最远。从而对研究生物的抗药性有一定的帮助。
12 进化树的评估
常用的评估进化树的方法有两种,分别是Bootstrap 和Jackknife 。所谓Bootstraping 法就是从整个序列的碱基(氨基酸)中任意选取一半,剩下的一半序列随机补齐组成一个新的序列。这样,一个序列就可以变成了许多序列,一个多序列组也就可以变成许多个多序列组。根据某种算法(最大简约性法、最大似然法、邻位相连法)每个多序列组都可以生成一个进化树。将生成的许多进化树进行比较,按照多数规则(majority-rule )我们就会得到一个最“逼真”的进化树。其数值反应了该树枝的可信的百分比。 所谓Jackknife 则是另外一种随机选取序列的方法。它与Bootstrap 法的区别是不将剩下的一半序列补齐,只生成一个缩短了一半的新序列。
结 论
分子生物学和比较基因组学是生物信息学研究的主要内容。昆虫的研究对于我们人类本身具有重大的意义。本文以昆虫基因为载体,介绍了几种常用的构建进化树的软件和方法。从基因的序列搜索、比对、及建立进化树都有说明。而且比较了建立进化树方法的优缺点和适用情况。但在实际中我们要根据所分析的序列特点,合理的选择构建进化树的方法。
在构建进化树的软件中,我们还可以根据需要美化构建的进化树,还可以使进化树呈现不同的形式。这些在本文中没有过多的介绍,可以根据自己的爱好,自己动手尝试。
进化树能很好的表现进化关系。可以直观的看出亲缘关系的远近。对人们了解进化有很大的用途。在昆虫基因研究中,我们可以通过构建进化树了解害虫的耐药性,从而更好的控制害虫。还对经济类昆虫提高产量有很大的意义。总之,昆虫基因的研究与人们的生产、生活方面都有重大的作用。
致 谢
在本毕业论文中,xx 老师的实事求是、刻苦钻研、认真负责的态度和深厚的理论水平都使我受益匪浅。在以后的学习和工作中我们都应该向他学习,始终保持这种精神,不骄不躁、脚踏实地得做好每一件事。所以,我首先要感谢x 老师给我的帮助和精神动力。感谢他不止在本论文中,还有大学四年期间对我的教导和鼓励。
在此还要感谢在这四年中我的所有任课老师,正是您们的辛勤教诲才能让我们在这四年中获得了知识,能够给自己正确的定位。在以后的学习、工作中能充分实现自己的价值。
最后,感谢同组同学的关心,感谢他们在百忙之中抽出时间给予我热心的帮助。
参 考 文 献
1 张赞,刘金定,黄水清,李飞. 生物信息学在昆虫学研究中的应用. 应用昆虫学报,2012,49(1):1-11
2 黄科,曹家树,吴秋云. 生物信息学. 情报学报,第21卷,第四期,2002 3 叶恭银,方琦. 基因组时代的昆虫学研究. 应用昆虫学报,2011,48(6):1531-1538
4 李涛,赖旭龙,钟杨. 利用DNA 序列构建系统树的方法. 遗传,2004,26(2):205~210
5 12 Maier D.The complexity of some problems on sobsequences and supersequences .ACM ,1978,25(7):322~336
6 Sourdis J,Nei M Relative effcieiencies of the maximum parsi- mony and distance matrix methods in obtaining the correct phyiogenetic tree MolBiol Ev01.1988,5(3):298~311
7 Holder M,Lewis P O.Phylogeny estimation:traditional and Bayesian approaches.Nature ,2003,4:275~284
8 Li W H Evolutionary change of restriction cleavage sites and phylogenetic inference Genetics,1986:187~213
9 Mike S.David P.Pasimony ,likelihood ,and the role of models In molecular phylogenetics.MolBioL Evol,2000,17:839~850
10 Goloboff P A Parsimony.1ikelihood ,and simplicity.Cladistics 2003,19:91~103
11 赵国屏.生物信息学.科学出版社,2002
12 徐立业. 多叉进化树构建方法的研究与实现. 硕士学位论文.2007 13 谭严芳,金人超. 一种基于NJ 的高效构建系统进化树算法.2004
14 董安国,高琳,赵建邦,呼加路. 基于DNA 序列的系统进化树构建. 自然科学
版,第36卷,第10期,2008 15 钟扬,张亮,赵琼.简明生物信息学北京:高等教育出版社,2001.116~149
16 马文静,马纪. 小胸鳖甲热激蛋白基因( Mphsp70) 的克隆、序列分析及温度对其表达的影响. 应用昆虫学报,2012,49(2):439~447