实验七 微生物分子生态学常用软件使用方法
微生物生态学研究中测序已经成为一项常规的必不可少的分析手段，实验后常常会得到大量的核酸序列，有的是细菌基因组上随机的序列片断，有的是16S rRNA基因的克隆文库，有的是功能基因序列等等，如此海量的序列数据，需要进行正确、快速和有效的分析，熟练掌握各种生物学软件的使用方法就显得尤为重要。这里我们主要介绍如何进行序列同源性分析，如何构建系统进化树，如何对克隆文库进行分析，如何对DNA 指纹图谱进行比较分析，介绍相关软件的使用方法。
一、实验原理
这里简要介绍序列数据分析过程中用到的软件：
BLAST 是NCBI(the National Center for Biotechnology Information)的一项服务。BLAST 在网络上可以直接使用，我们可以提交序列，并与NCBI 数据库（GenBank+EMBL+DDBJ+PDB sequences ）进行比对，之后会将一系列的结果返回给用户。
GeneTool 可以进行核酸分析，本文中主要用于去除载体序列。
ClustalX 1.8：广泛使用的多序列比对程序，在ClustalW 多序列比对程序的基础上增加了图形用户界面。输入为多序列的Fasta 格式文件，进行多序列全局比对生成结果的同时，在指定文件夹生成“.dnd”和“.aln”格式文件。
PhyloDraw 0.8：构建进化树的绘图工具，它支持多种多序列比对软件的Multiple Alignment结果。本实验采用ClustalX 进行多序列比对，生成“.dnd”格式的比对文件，最后用PhyloDraw 画出序列进化树。它支持Unrooted tree(无根树) 、Rooted tree(有根树) 、Radial tree(放射状树) 、Rectangle cladogram(矩形进化分支树) 、Slated cladogram和Phylogram(序列进化树) 。这些都是不同的树型，结果是一致的。
下面简要说明Blast 、Fasta 、Cluastx 、PhyloDraw 等进行序列比对以及构建进化树的算法等，作为深入研究的理论基础。
DNA 序列的比对是生物信息学的基础之一，寻找序列相似性的过程称为序列比对。 系统进化推断是通过生物间可观测的性质来建立物种之间进化关系假说的方法。我们的目的是构建系统进化树，它已成为相似性比对为基础表示进化关系的很直观的方法。系统进化树是严格的二叉树，二叉分支假设极大的简化了建树算法。在系统进化树中，序列之间的进化距离可以作为树枝长度的度量。构建系统进化树的方法很多，主要有以下四种方法：
（1）基于成对距离比对的系统进化树：这种方法能够生成有根的树，这种方法首先通过定义每一对序列之间的距离矩阵初始化，然后按照距离分组，最后建立从树枝到树根的树。
（2）基于相邻连接的系统进化树：这种方法不仅根据距离矩阵搜索最小的成对距离，而且会搜索使整个树高最小的相邻集，最适合于进化距离较短的情况。
（3）基于最大节约法的系统进化树：这种方法是从一系列可能的树中找到一个需要最少的核苷酸替换就可以解释所看到差异的树。
（4）基于最大似然法的系统进化树：这是一种概率法，它通过在信息位点的每种可能的进化改变的概率排列并使树的总概率最大化来寻找最佳选择。
最常见的UPGMA ，它的全称是使用算数平均数的未加权对群法(Unweighted Pair Group Method with Arithmetic)，该算法属于基于成对距离比对生成系统进化树算法。
z UPGMA 算法描述：
初始化过程 ：
(1) 为每个物种建立一个群(Cluster)；
(2) 每个群的大小赋初值n i =1，即只包括一个物种；
(3) 计算任意两个群之间的距离D ij ，采用二维数组存储该距离矩阵；
(4) 输出树为T ，为每个物种分配一个叶结点。
循环过程：每一个循环都能将将其中两个群合成为一个群
(1) 遍历所有距离值，找出具有最近距离D ij 的两个群i 和j ；
(2) 创建一个新的群(ij)，它共有n (ij)=n i +n j 个物种；
(3) 连接树上的i 和j 到一个新的结点，该结点对应于新的群(ij)，连接i 和j 的树枝长度为
D i , j
2；
(4) 按照下面公式计算从新的群到其它每个群的距离（不包括i 和j ）
n j n i D (ij ), k =(D i , k +(D j , k n i +n j n i +n j
(5) 删除距离矩阵中i 和j 的记录，添加新的记录D (ij),k
(6) 返回1直到只剩下一个群；
综上所述，该算法主要思想是首先将每个序列被分配到自己的群中，从树的零高度开始这个序列的分支，找出距离最近的两个群合并为一个群，直到剩下一个群为之。树枝的长度
反应两者之间的距离，即进化时间的长短，构造的顺序是从树枝到树根逐渐构造。
二、实验目的
分析实验数据，大量的序列数据信息分析整理，进行同源性比较、构建进化树、分析指纹图谱的相似性等。
三、实验材料
(1) 以本室的序列数据为例，介绍相关软件的使用方法。
(2) 生物学软件:如上列举的软件。
四、操作步骤
4.1 序列分析及进化树构建
4.1.1. 去除载体序列，目标序列经克隆（以克隆到 Promega 公司生产的pGEM-T 载体为例）通用引物（T7/SP6）测序后，测序结果中带有部分载体序列，在进行序列分析以前，要首先去除载体序列，可以使用DNAMAN 和GeneTool 等等，这里以GeneTool 为例进行说明。 下图是GeneTool 软件的主界面：
将需要去载体序列的文件（Raw Sequence）打开，搜索EcoR I 酶切位点 GAATTC （pGEM-T easy 载体两端均有该位点）, 以下图为例，上游和下游的序列均用蓝色标出
从上游位点向后第5个碱基开始，下游位点向前第11个碱基开始为我们的目标片断, 选中之后输出，保存为FASTA 格式的”TXT”文件，这样就完成了我们去载体序列的过程。
4.1.2. 到GenBank 数据库中进行Blastn 分析, 找到其Closest Relatives： 打开http://www.ncbi.nlm.nih.gov/BLAST/如下图所示：
选择做Nucleotide-nucleotide BLAST (Blastn)弹出界面：
将要比对的序列填入Search 中，下面以GenBank 中公开的NC_003045序列（Bovine coronavirus, complete genome）为例进行分析。
参数可以采用默认值，之后点击BLAST ，以可以得到下图的结果：其中Query = (31,028
letters) 表示我们序列全长为31028bp ，查询的ID 为
1089336352-4777-[1**********]0.BLASTQ4，
点击Format 开始搜索GenBank 数据库，采用Blast 的方法，并将结果返回给用户,
下表显示了GenBank 中和序列NC_003045做Blast 之后得到的序列，从上到下同源性降低。
下图表示AF391541.1序列和NC_003045序列的配对情况。
我们可以将排在最前面的（1－2条）序列下载下来，这就是与测许序列亲缘关系最近的序列信息（Closest Relatives）。值得注意的是，有时候序列之间是反向匹配的，我们需要将测许序列顺过来，使其方向从小到大，这样就完成了第二步——在数据库中寻找同源序列。
4.1.3. 通过Clustalx 软件和PhyloDraw 构建系统进化树。
第二步中，我们可以将一个文库中测序得到的每个序列都到GenBank 中找到与其相似的同源性最高的序列，用这些序列进行构建进化树时，进化树中将包含大量已知序列，可以作为目标序列进化地位的参考。在下图的范例中，我们L-46c 、L-33c 、L-18、L-67、L-92、L-ASa 、L-ASb 、L-46d 、L-33d 均为我们实验室测序得到的序列，有了从GenBank 中得到的亲缘关系较近的序列，生成的进化树中能够一目了然地看到这些未知序列的分类地位。
以下介绍如何使用ClustalX 软件和PhyloDraw 构建系统进化树
首先，将所有的要构建系统进化树的序列存储到一个文件中，采用FASTA 格式，FASTA 格式又称为Pearson 格式，这是比较简单而使用最多的序列格式。序列文件的第一行是由大于符号开头的任意文字说明，主要为标记序列用，从第二行开始为序列本身。碱基名称大小写均可，如下所示：
> sequence1
acttaaaaagattttctatctacggatagttagctctttttctagaccttgtctactcaa
ttcaactaaacagaaattttgtccttccttccggccgcatgttcatgctgctggaagctg
> sequence2
acttaaaaagattttctatctacggatagttagctctttttctagaccttgtctactcaa
ttcaactaaacagaaattttgtccttccttccggccgcatgttcatgctgctggaagctg
> sequence3
acttaaaaagattttctatctacggatagttagctctttttctagaccttgtctactcaa
ttcaactaaacagaaattttgtccttccttccggccgcatgttcatgctgctggaagctg
其次，使用ClustalX 打开这个存储多个序列的文件
之后做Alignment, 将输出“.dnd”格式的文件作为PhyloDraw 的输入。
最后，使用PhyloDraw 打开上一步生成的“.dnd”文件，
可以有不同的树型供选择：
Unrooted tree
Rooted tree
Radial tree
Phylogram
通过上面的去除载体序列，到GenBank 下载Closest relatives序列，以及最后的构建进化树三个步骤，就完成了对序列进行分析的一个主要工作。
4.2 使用RDP 在线分析数据及构建进化树
RDP（Ribosomal Database Project）主要是提供关于核糖体相关的序列数据，它可以实现在线的构建进化树，序列比对等。实验中得到的16S rRNA序列信息可以在线使用RDP 进行分析。
下面简要介绍如何在线使用RDP 分析序列并构建系统进化树。
打开http://rdp.cme.msu.edu/cgis/phylip.cgi这是一个通过web 接口提供Phylip 和weighbor 服务的程序。我们可以使用这个程序来创建距离矩阵（相似性矩阵）和构建系统进化树。
我们可以使用自己的序列和RDP 数据库中序列结合在一起构建系统进化树，步骤如下： 1. 编辑要构建进化树的数据集合Edit Data Set:
首先将序列上传到RDP 中，点击Edit Data set,之后，可以通过browse 打开本地序列，见下面的左图。另外还可以选择和数据库里的0－10条最接近的序列进行比对，这里我们选择Include 2 neighbors，见下面右图
将自己要构建系统进化树的序列一个个上传到数据库中
2. 计算距离矩阵Distance Matrix
点击Distance Matrix，弹出新的对话框，选择Calculate Matrix,一共有4中计算距离矩阵的方法：Kimura 2-parameter； Jin/Nei, coeff；
Maximum Likelihood； Jukes-Cantor。
根据需要选择其中一种即可。
3. 构建系统进化树Phylogenetic Tree:
第二步完成之后点击Phylogenetic tree，弹出新的对话框，可以选择Neighbor-joining 或者Weighbor-joining 之后点击Calculate tree可以得到系统进化树了。
通过以上步骤，我们实现在在线使用RDP 构建系统进化树的过程。
4.3 文库分析（richness evenness rank abundance coverage ）
4.3.1 生态学家根据两个参数定义了物种多样性species richness：（1）种群当中物种的个数，生态学家通常把它叫做物种丰度species abundance，和（2）物种的相对丰富程
度，或者称作物种均度 species evenness 。物种丰度对于种群的多样性的影响是十分明显的。一个拥有20个树种森林的种群显然比一个拥有80个树种森林的种群的多样性要差。
4.3.2 一种物种多样性的定量指标
生态学家建立了许多关于物种多样性的指标，它的数据值取决于物种丰度和物种均度的级别。一种最常使用的测量物种多样性的指标是Shannon －Wiener （香侬）指标：
H ' =−∑p i log e p i
i =1
其中：
H ' 表示Shannon －Wiener 多样性指标的值
p i 表示第i 个物种所占的百分比
log e 表示p i 的自然对数
s s 表示种群当中含有多少种不同的物种数目
4.3.3 丰富程度等级曲线(Rank aboundance)
我们还可以描述一个种群的各个物种之间的相对丰富程度和物种的多样性，方法是根据它们丰富程度的等级来划分物种的相对丰富程度。得到的丰富程度等级曲线(Rank aboundance) 为我们提供了关于一个种群的重要信息。
此外，有关文库分析的详细说明，参见 “实验四”部分。
4.4 指纹图谱分析（UVI 凝胶成像系统）
通过PCR －TGGE ，LP-RAPD, T-RFLP等得到的微生物群落结构的指纹图谱，常常需要分析图谱中条带间的相似性，Di Giovanni 在分析ERIC-PCR 指纹图谱时，用的是SYSTAT V. 7.0 (SPSS, Inc., Chicago, IL) (1) 基于 1 和 0 的算法；Franklin 对RAPD 指纹图谱的分析是基于Jaccard coefficient算法，计算的是2个图谱之间共有条带的情况 (2)；Zoetendal 对TGGE 指纹图谱的分析是通过arithmetic averages (UPGMA)软件分析 (3)。我们采用UVI (UVItec, Cambridge, U K)凝胶成像系统自带的分析软件，进行聚类分析。该软件提供了2种算法，① Nei and Li coefficient (也就是Dice Coefficient), 基于Coeff: a=2nxy/(nx+ny),其中nx 和ny 分别代表泳道x 和y 的条带数，nxy 为两个泳道共有的条带；② Jaccard Coefficient, 基于Coeff: b=nxy/(nx+ny-nxy)来分析算法。
4.4.1 打开 UVI BandMap 软件界面，打开要分析指纹图谱，如下图，
4.4.2 根据泳道数设定Lanes ，然后Detect Bands，需要手动选定条带，
4.4.3 点击Detect Bands,生成图谱，选择算法：Jaccard Coefficient或Dice
Coefficient 算法, 生成树状图，
4.4.4 拷贝出树状图，即得到指纹图谱的聚类分析图, 两两相似性就可以直接读
出。
4.5 计算机模拟PCR
4.5.1 PCR简介
聚合酶链式反应(PCR)是体外扩增DNA 序列的技术：DNA 复制不是像分子克隆技术中那样在细胞中进行，而是用纯化的酶进行。基本的PCR 是利用一对根据所扩增的目的片断所设计的引物进行的。引物与对应的DNA 链退火，DNA 合成向中心区域进行。反应涉及三个不同的温度阶段：双链DNA 的变性（90℃以上进行），引物与单链模板退火（50℃左右进行），引物延伸，合成新的DNA 越过靶区域（70℃左右进行）。每组三个反应称为一个PCR 循环。我们在研究中往往会对基因组中部分序列片断感兴趣，有了这项技术，就可以针对研究目的为保守区域设计引物，用它从基因组中“扩增”出目标片断，进行分析研究。
4.5.2 计算机模拟PCR 过程
生物学的研究工作过去是完全依赖于实验的，现在依靠计算科学的发展和生物信息学的进步，越来越多的实验工作可以在计算机的指导下进行了。分子生物学的研究方法将转变为理论分析的指导下，将实验与理论相结合的研究过程。
以PCR 实验为例，随着计算机模拟PCR 算法的研究和实现，预测PCR 产物这部分的工作完全可以由计算机来做。通过计算机模拟可以大大减少我们实际实验的盲目性。到目前为止，已经出现了若干个用于PCR 产物预测的算法和程序。
SPCR 软件是我们实验室自主开发的用于预测PCR 产物的软件。
它是基于信息系数的PCR 产物预测的方法。我们利用信息系数(Information Coefficient)的理论设计了一种比较DNA 序列之间相似性的算法，并把这一算法用于PCR 产物的预测。信息系数是一种基于信息论的系数，通常用于计算两个信息源之间共有信息量的多少，目前多用于数量分类。
DNA 分子本身是一个遗传信息的载体，而DNA 的复制和转录可以看作是信息的传输，而PCR 的过程其实就是特定DNA 分子大量被复制的过程，因此也可以把它看作是一种信息传输的过程。我们可以把模板和引物看作是两个信息源，用信息系数来计算两个信息源之间共有的信息量，以此来描述两个DNA 的相似性，并通过这种相似性来选择引物在模板上的退火位点。我们在碱基数理论的基础上利用信息系数设计了PCR 产物预测的算法，并开发了用于PCR 产物预测的程序——SPCR。
该程序在Windows 平台开发完成，用C++编写，编译之后生成一个可执行文件，直接运行，无需安装。SPCR 可以识别IUPAC 的核酸代码，因此可以用IUPAC 来表示简并碱基。SPCR 的主要参数有上游引物退火系数(I up ) 、下游引物退火系数(I dn ) 、产物生成概率(P a ) 、最大和最小产物长度(L max , Lmin ) 以及引物序列文件、基因组序列文件和结果数据文件。计算完成之后SPCR 会把计算结果数据(包括所有产物序列) 以及计算过程中使用的所有参数以文本格式保存到指定的结果数据文件中。另外SPCR 可以对模拟扩增结果模拟Agarose 凝胶电泳图谱。
SPCR 的时间复杂度和基因组序列的长度成正比，并且可以模拟多基因组的扩增。 下面是该软件的主界面：
Primer 为引物信息：分别为上游引物和下游引物
Template 中填入模板信息
Up Threshold和Down Threshold为上游引物和下游引物阈值
Max Product用来设定PCR 产物长度的上限
Min Product用来设定PCR 产物长度的下限
点击Sim Gel可以模拟凝胶电泳图
下面我们结合一个实例来讨论SPCR 软件如何进行模拟PCR 实验的。
实验材料：
上游引物：IN-2(+)：5'-GGGTTGGGACTATCCTAAGTGTGA-3'
下游引物：IN-4(-) ：5'-TAACACACAACNCCATCATCA-3'
模板序列为：NC_004718
实验步骤：
将上游引物填入Primer1, 下游引物填入Primer2。点击Add Template选择模板序列路径，比如F:\NC_004718.seq,如果要删除该模板可以选择Remove 。将SPCR 的实验结果保存在Result Data 中，我们可以新建一个result.txt 文件，然后在Result Data一栏中点击Browse 将result.txt 打开。这样就可以将结果就存储到result.txt 文件中了。Up Threshold, Down Threshold, Pa Threshold, Max Product和Min Product都有初始值, 我们可以根据需要更改在Sim Gel一栏中点击Browse 打开result.txt 可以得到如下右图所示的模拟凝胶电泳图
使用SPCR 软件进行PCR 产物预测准确快速，在进行实验前，使用SPCR 模拟PCR ，模拟结果可以作为实验的重要参考。
附录：网上生物信息资源
主要数据库的地址
数据库
MEDLINE
GenBank 组织 National Library of Medicine 地址 www.nlm.nih.gov National Center for Biotechnology
www.ncbi.nlm.nih.gov
EMBL
RDP
DDBJ
SWISS-PROT European Bioinformatics Institute Michigan State Univesity National Institute of Genetics, Japan Swiss Institute of Bioinformatics www.ebi.ac.uk http://rdp.cme.msu.edu/html/ www.ddbj.nig.ac.jp www.expasy.ch
PIR National Biomedical Research www.nbrf.georgetown.edu
Foundation PRF
PDB Protein Research Foundation, Japan www.prf.or.jp Research Collaboratory for Structural www.rcsb.org
Bioinformatics CSD Cambridge Crystallographic Data www.ccdc.cam.ac.uk
Center
新一代分子生物学数据库
信息 数据库 地址
www.genome.ad.jp/dbget/ligand.html www.genome.ad.jp/dbget/aaindex.html
www.expasy.ch/sprot/prosite.html
www.blocks.fhcrc.org/
www.biochem.ucl.ac.uk/bsm/dbbrowser/PRINTS/
www.sanger.ac.uk/Pfam/pfam/wustl.edu protein.toulouse.inra.fr/prodom.html
scop.mrc-lmb.cam.ac.uk/bsm/cath www.biochem.ucl.ac.uk/bsm/cath
www.ncbi.nlm.nih.gov/COG/ www.genome.ad.jp/kegg
www.genome.ad.jp/kegg
www.mcs.anl.gov/WIT2/ ecocyc.PangeaSystems.com/ecocyc/ 化合物和反应 LIGAND AAindex 蛋白质家族和Motif 序列 PROSITE Blocks PRINTS Pfam ProDom 三维折叠子分类 SCOP CATH 直系同源基因 COG KEGG 生化途径 KEGG WIT EcoCyc
基因组多样性 NCBI
Taxonomy
OMIM
www.ncbi.nlm.nih.gov/ www.ncbi.nlm.nih.gov/Taxonomy/ .www.ncbi.nlm.nih.gov/Omim/
五、参考文献
1. Di Giovanni, G. D., L. S. Watrud, R. J. Seidler, and F. Widmer. 1999. Comparison of Parental
and Transgenic Alfalfa Rhizosphere Bacterial Communities Using Biolog GN Metabolic
Fingerprinting and Enterobacterial Repetitive Intergenic Consensus Sequence-PCR
(ERIC-PCR). Microb Ecol 37:129-139.
2. Franklin, R. B., D. R. Taylor, and A. L. Mills. 1999. Characterization of microbial communities
using randomly amplified polymorphic DNA (RAPD). J Microbiol Methods 35:225-35.
3. Zoetendal, E. G., A. von Wright, T. Vilpponen-Salmela, K. Ben-Amor, A. D. Akkermans, and
W. M. de Vos. 2002. Mucosa-associated bacteria in the human gastrointestinal tract are
uniformly distributed along the colon and differ from the community recovered from feces.
Appl Environ Microbiol 68:3401-7.
4. R.Duibin, S.Eddy, A.Krogh, G.Mitchison, Biological sequence analysis, Cambridge University
Press, 1998
5. Minoru Kanehisa, Post-genome Informatics, Oxford University Press, 2001
6. Paul A.Rota et al,Characterization of a Novel Coronavirus Associated with Severe Acute
Respiratory Syndrome,Sciencexpress, May 2003
7. Ron Shamir, Algorithms for Molecular Biology, 2001
8. 郝柏林，张淑誉编著 生物信息学手册 上海科学技术出版社。
（本章作者：张宇镭）