2009年  第54卷  第9期: 1250 ~ 1261
www.scichina  csb.scichina
《中国科学》杂志社
SCIENCE IN CHINA PRESS
论 文
美姑脊蛇Achalinus meiguensis 线粒体基因组全序列及系统发育地位研究
王广力①, 何舜平②, 黄松③, 何苗①, 赵尔宓①④*
① 四川大学生命科学学院, 成都 610064; ② 中国科学院水生生物研究所, 武汉 430072; ③ 黄山学院生命与环境科学学院, 黄山 245041; ④ 中国科学院成都生物研究所, 成都 610041 * 联系人, E-mail: zem006@163.
bootstrap 软件com  2008-10-21收稿, 2009-03-13接受
摘要  报道了美姑脊蛇Achalinus meiguensis 线粒体基因组全序列. 美姑脊蛇线粒体全序列长17239 bp, 由22个tRNA, 2个rRNA 和13个蛋白质基因及2个非编码的控制区或D-loop 组成, 存在着基因重排现象. 对已报道蛇类线粒体基因组全序列进行比对分析后, 发现一些蛇类线粒体基因组进化规律: 双控制区现象在爬行动物进化历史中独立地发生, 有不同的演化历史; tRNA 假基因是在真蛇下目(Caenophidia)中进化形成的; T ΨC 臂的相对较短(一般少于  5 bp)和缺失“DHU”臂造成蛇类tRNA 较短. 通过MP 和BI 分析确定美姑脊蛇系统发育位置应该位于瘰鳞蛇Acrochordus granulatus 和其他真蛇下目蛇类之间, 而不应该属于游蛇科Colubridae 或者眼镜蛇科Elapidae. 由于美姑脊蛇与真蛇下目中蝰科Viperidae 、游蛇科以及眼镜蛇科之间的单系性都被统计检验拒绝(P  < 0.01), 由此认为美姑脊蛇所在的闪皮蛇亚科(Subfamily, Xenodermatinae)应该提升到科或者更高的分类阶元. 依据系统发育统计检验结果, 我们选择Bayesian 树来进行分歧时间的估算. 原蛇下目与真蛇下目的分化时间为109.50 Mya; 而瘰鳞蛇与其他真蛇下目种类的分化时间为106.18 Mya; 美姑脊蛇的分化时间在103 Mya; 蝰科与眼镜蛇科和游蛇科构成的单系的分化时间发生在96.06 Mya.
关键词  美姑脊蛇 线粒体全序列系统发育 分化时间
在不同动物类中, 线粒体基因组在基因结构和基因排列方式等方面均显示了极大的多样性, 这种多样性
可能反映了真核细胞不同的进化路线[1]. 就目前的研究而言, 线粒体基因组是一个能够从基因组水平上来分析动物系统发生的分子标记, 可以从线粒体基因组序列信息、基因组成及基因排列方式等进行多方位的分子进化研究. 因而线粒体基因组全序列成为研究动物分子系统发生最有力的证据之一[1]
. 蛇的线粒体基因组大体构造与其他脊椎动物相似, 2个rRNA, 22个左右的tRNA 和13个左右的蛋白编码基因, 至少一个非编码区(control region), 但是蛇类线粒体全序列在其组成、结构与功能等方面又有着一
些与众不同的特点, 如紧凑的基因结构、双控制区、
缩短的tRNA 以及较快的进化速度等[2~5].
脊蛇属Genus Achalinus  Peters, 1869在分类上隶属于游蛇科Family Colubridae, 闪皮蛇亚科Subfamily Xenodermatinae [6,7]. 该属的模式种是黑脊蛇(Acha- linus spinalis ). 脊蛇属蛇类已知有9种, 除1种外, 在中国均有分布, 其中5种是中国特有种. 美姑脊蛇(Achalinus meiguensis  Hu and Zhao, 1966)是中国的特产蛇类, 区别于其他脊蛇的主要形态特征是具鼻间鳞和一枚极小的眶后鳞[8]. 目前仅见于四川的美姑、安县、洪雅、峨眉、宝兴、汶川、卧龙、屏山等以及云南省水富县, 已报道的分布海拔在1500~3000 m 之
间[7,9]. 由于分布的局限性及生境的特殊性, 关于美姑脊蛇研究的文献很少, 而且这些研究大多为形态描述, 对其生态学或者分子生物学研究尚未见报道.
本文首次报道美姑脊蛇的线粒体基因组全序列, 并通过与蛇类已有的线粒体基因组全序列比较, 探讨蛇类线粒体基因组全序列的结构及进化等特点. 闪皮蛇亚科及脊蛇属的系统发生研究中存在基因少、支持率低等问题, 本文使用线粒体基因组全序列来进行系统分析, 试图更好地解决它们的系统关系.
1  材料与方法
(ⅰ) 基因组DNA的提取.  美姑脊蛇采自四川省安县千佛山自然保护区, 实验材料用95%的酒精保存, 标本用6%福尔马林溶液固定, 并保存于四川大学. DNA的提取参照标准酚-氯仿法进行[10], 所提取的DNA置于−20℃保存备用.
(ⅱ) 引物设计.  首先, 使用蛇类编码基因通用引物[11~15]扩增出美姑脊蛇的ND1, ND2, ND4和Cyt b片段. 然后与从GenBank中下载的阿尔卡链蛇(Dinodon semi-carinatus, NC_001945)、冲绳烙铁头(Ovophis oki-navensis, NC_007397)、瘰鳞蛇(Acrochordus granula- tus, NC_007400)和球蟒(Python regius, NC_007399)等线粒体基因组全序列进行排列比对, 选取保守性较高的区域, 设计了可以覆盖美姑脊蛇线粒体基因组全序列的11对PCR 引物(表1).
(ⅲ) PCR扩增和序列测序.  PCR反应体系大概含有约50 ng的模板DNA, 10×缓冲液5 μL, 2.5 mmol/L dNTP各2 μL, 引物(10 pmol/L)各2.5 μL, Taq聚合酶2.0 U, 最后补足灭菌双蒸水至终体积. PCR反应条件为: 95℃预变性3 min; 94℃变性30 s, 48~58℃退火30 s, 72℃延伸90 s, 35个循环; 最后72℃延伸8 min. 为防止污染, 所有PCR反应都设立空白对照组; 扩增产物使用0.8%的琼脂糖凝胶电泳检测, 用BioStar Glass-milk DNA纯化试剂盒(加拿大)回收, 纯化后的线粒体基因使用ABI3730测序仪(美国)测序. 为保证所测序列的可靠性, 每一样品均进行双向测序.
(ⅳ) 序列分析.  获得的序列用Blast软件与GeneBank数据库的序列比对, 以及用Clustal X 1.83软件[16]与已发表的近缘物种的线粒体基因组全序列进行比对, 编码蛋白序列的准确起止位置根据脊椎动物线粒体基因组起止密码子确定, 并使用Sequin (Version5135)软件确定了13个蛋白质基因的位置和DAMBE软件来计算序列的碱基组成和密码子使用频率. 通过tRNA Scan-SE软件[17]定位了21个tRNA基因, 通过Clustal X 1.83软件与4种蛇线粒体基因组全序列进行比对, 定位了剩余的一个tRNA-Ser (AGC)和2个rRNA基因, 并且对其余的基因位置进行了修正和确定. tRNA的二级结构用RNAstructure4.5 (hester.edu/rnastructure.htmLmL)软件推测.
(ⅴ) 系统发育分析.  从GenBank 上下载16种蛇的线粒体基因组全序列(表2), 在13个编码基因中, 由于ND6在轻链编码, 在核苷酸和氨基酸组成上与其他12个编码基因差异明显, 因而在构建系统发生序列矩阵时该基因被排除[20,22,28,29]. 为了考察蛇亚目的单系性, 我们分别选择两栖类、龟类、鳄鱼类、鸟类、
蜥蜴等类的代表种作外类(表2).
系统发育分析分别采用最大简约法(maximum
表1  美姑脊蛇线粒体基因组全序列测定所使用的引物
引物编号上游序列(5′→3′) 下游序列(5′→3′) 退火温度/℃测序长度/bp
E39734/E39735 TAATCTCCGTGCCAGCGAC CTTCACAGGGTCTTCTGGTCTT 58 1776
E45266/E45267 GCCTGCCCAGTGAACATT GCGTTCTAGGAGGGTGAGG 55 779
E46268/E45269 TCCGCTACGACCAACTCAT TGCTCAGGCTATTAGTCAGTG 52 1654
E40978/E40979 CTTCTACAGTAAGGTCAGC TTAGTCAGTTGCCAAAGCC 50 1847
E40980/E40981 TGGACATTGACAGACGAGC TATGCGGTGGTCTACTTC 49 1260
E54599/E54600 CGATCCCTCTACCTAATAGAAGAAGT GGAGTAGGAGGATGATTTCTAGGTCAA 48 2393
E42227/E42228 ACCCACGCCCTAATACTCA ATGTGGCTGATAGAGGAG 52 1919
E42229/E42230 CATCTCTGACTACCAAAAGC GTCCAATGTCTCCGATACG 51 1435
E42046/E42047 CAACGAACAAGACATCCG TCCACCAGGTTGAGATGTT 52 1398
E42231/E45098 CATCCTACGCTCTATCCC ATGGTAGTCAGGTGAAAGG 48
995
E45099/E42232 CCACTGGTTACACTCTCAAG GTGATTAGTTTGGCTTATGG 50 1020
1251
2009年5月第54卷第9期
表2  本研究所用的线粒体基因组序列
物种名称 GenBank登录号来源及出处
Acrochordus granulatus NC_007400 [18]
Agkistrodon piscivorus NC_009768 [5]
Boa constrictor NC_007398 [18]
Cylindrophis ruffus NC_007401 [18]
Deinagkistrodon acutus NC_010223 Yan J, et al. 未发表
Dinodon semicarinatus NC_001945 [3] Enhydris plumbea NC_010200 Yan J, et al. 未发表
Eunectes notaeus AM236347 [19]
Leptotyphlops dulcis NC_005961 [20]
Naja naja NC_010225 Yan J, et al. 未发表
Ovophis okinavensis NC_007397 [18]
Pantherophis guttatus guttatus AM236349 [19]
Pantherophis slowinskii NC_009769 [5]
Python regius NC_007399 [18]
Ramphotyphlops braminus NC_010196 Yan J, et al. 未发表
Xenopeltis unicolor NC_007402 [18]
Pelomedusa subrufa NC_001947 [21]
Dogania subplana NC_002780 Farajallah A, et al.未发表Iguana iguana NC_002793 [22]
Sceloporus occidentalis NC_005960 [20]
Rhineura floridana NC_006282 [23]
Geocalamus acutus NC_006285 [23]
Corvus frugilegus NC_002069 [24]
Vidua chalybeata NC_000880 [25]
Smithornis sharpei NC_000879 [25]
Acanthisitta chloris NC_003128 [26]
Alligator sinensis NC_004448 Wu X, et al.未发表
Gavialis gangeticus NC_008241 已投稿
Xenopus laevis NC_001573 [27]
Achalinus meiguensis 本研究
parsimony, MP)和贝叶斯法(Bayesian inference, BI)进行分析. MP分析在PAUP*4b10软件[30]中完成. 分析选择启发式搜索(heuristic search)、树二等分再连接选项(tree-bisection reconnection, TBR)和100次随机加入序列(addition sequence replicates)等参数进行搜索最大简约树. 所有位点均作无序特征处理, 分析前排除所有无信息位点. 检验MP 树中各节点可信度采用非参数自展法(non-parametric bootstrap)重复检测, 参数设定: nreps = 1000, search = heuristic, conlevel = 50, addseq = random, swap = TBR.
贝叶斯分析采用MrBayes3.1.2软件[31], 以后验概率(posterior probability, PP)来表示各分支的可信度. 替代模型根据MODELTEST v.3.7软件[32]中等级似然率(hLRT)和AIC检验标准来选择替代模型. 起始树设为
随机树, 马尔科夫链的蒙特卡洛方法设置为4条链同时运行3×106代, 3条热链1条冷链. 每100代对系统树进行抽样, 最终得到30001棵系统发育树, 重复一次以确保MCMC的收敛. 将运行过程中所得的冷链对数似然值(log-1ikelihood scores)与相应的代数进行作图, 到对数似然值达到饱和的位置, 到达饱和前的数据被作为老化样本(burnin samples)而舍弃. 在舍弃老化样本后, 根据剩余的样本构建一致树(consensus tree)并计算相关参数.
(ⅵ) 系统发育树检测.  系统发育树构建完毕后, 采用Shimodaira-Hasegawa (SH)检测[33]判断不同构树方法所构成的系统发育树之间是否收敛; 而从一系列相互竞争的拓扑树或系统育假设中选择最优的系统发育树时, Shimodaira-Hasegawa检测和approximately unbi-ased (AU)检测[34]两种统计检验同时被采用, 前者包含在PAUP*4.0b10软件中, 后者则通过CONSEL 0.1f 软件[35]完成, 并且所有统计检验都是基于同一核苷酸替代模型, 即GTR + G + I模型(根据等级似然率和
1252
AIC 选择标准选择最适合的替代模型).
(ⅶ) 美姑脊蛇分化时间.  为了检验是否偏离分子钟, 将有分子钟约束的ML 树与不受分子钟约束的ML 树的似然对数值进行比较[36]; 如果两者差异显著, 则分子钟被拒绝. 对数似然比检验表明, 不同谱系联合线粒体蛋白编码基因数据(χ2 = 1443.035, d f  = 28, P  < 0.001)序列的分子钟假设被拒绝, 不宜采用全局分
子钟(global clock)来估算分支时间.
因此本研究采用松散分子钟约束的Bayesian 方法
[37]
和罚分似然法[38]估算分支之间的分歧时间. 这两种
方法都放宽了对严格分子钟中系统发育类间相同的、稳定的替代速率的约束, 而是允许数个分子钟及分化时间校正. Bayesian 方法在MULTIDISTRIBUTE 软
件[39]中执行, 并将不同基因、
不同进化谱系的进化速率的差异进行整合. 这一分子钟估算方法的运算过程可以分为两步: 第一步, 从数据和基于F84+G 进化
替代模型获得的稳定拓扑结构的系统树中, 运行ESTBRANCHES 来估计分支的枝长. 这一过程允许位点进化速率按照分散的γ 分布而变异, 这一分布依照4种替代形式速率的方差-协方差数据矩阵计算[40,41]
. F84+G 进化替代模型中的参数用PAML 软件的BASEML 程序[42]来估计. 第二步, 用MULTI- DIV-TIME [39]估算分化事件的时间及标准差, 并且通过Markov chain Monte Carlo 方法来估算95%水平的置信区间. Markov 链运行100000代, 在舍弃初始的300000次循环的代数后, 每100代抽样一次. 在运行MCMC 之前, 我们先对根节(tree root)年龄、根节速率以及分支自主相关速率的平均值及标准误等参数进行预设: 以30
0 Mya 前(标准误SD = 50 Mya)作为根节到顶端分化时间的期望值(rttm = 300), 但不对节点分化时间进行限定; 每位点每百万年0.002655的替换(SD = 0.002655)作为tree root 节点的进化速率(rtrate = 0.002655; rtratesd=0.002655); 系统树上沿着衍生枝, 控制每百万年速率与分支自主相关程度尺度的参数设置为.0361 (brownmean = 0.0361, brownsd = 0.0361). 根节到顶端的分化时间最大可能值设置为350 Mya (bigtime = 350). 为检测Markov chain Monte Carlo 的收敛性, 这一过程以不同的起始树至少运行2次.
0MULTIDIVTIME 程序允许对标定点进行最大、最小时间尺度的限定. 本研究采用两个校正点:
(1) 蛇类和蜥蜴类分化时间大约在130~150 Mya [43];
(2) Acanthisitta chloris , Corvus frugilegus 和Vidua cha-lybeata 的分化时间最少82 Mya [44].
罚分似然法在R8S 软件[45]中执行, 该种方法通过最优树的枝长来推算分支之间的分化时间, 采用 Powell 和TN 方法进行优化, 最佳平滑值由交叉验证法(Cross-validation)获得. 同时, 分子钟的分化时间校正也采用以上两个标定点.
2  结果与分析
2.1  美姑脊蛇mtDNA 结构
美姑脊蛇线粒体基因组序列全长为17239 bp, 为一环形的双链DNA (图1). 与大多数蛇类的基因组成及顺序相同[3], 美姑脊蛇线粒体基因组由13个蛋白编码基因, 22个tRNA 基因和2个rRNA 基因(16S rRNA 和12S rRNA)组成. 其中主链(majority-strand), 即H-链编码了29个基因, 而次链(minority-strand), L-链包含了8个基因. 另外将在tRNA-Ile 和tRNA-leu 间以及tRNA-Pro 和tRNA-Phe 间区域的1062和1063 bp 确定为两个非编码区(non-coding region), 即控制区(control region, CR).
(ⅰ) 蛋白编码基因.  在13个蛋白质编码基因
图1  美姑脊蛇(Achalinus meiguensis )线粒体基因组基因
构成
tRNA 基因以相应氨基酸的单个大写字母表示. 其他基因的缩写分别为: ND1~6为NADH 脱氢酶亚基1~6; CO Ⅰ~Ⅲ为细胞素氧
化酶亚基Ⅰ~Ⅲ; ATP6~8为腺苷三磷酸酶亚基6~8; Cyt b 为细胞 素b. O H 和O L 分别表示重链和轻链的复制起点
1253
2009年5月  第54卷  第9期
1254
中, 除ND6定位在轻链外, 其他均定位在重链上. 这些基因都没有内含子, 有些与其他基因有少许重叠(表3). 细胞素氧化酶亚基Ⅱ和Ⅲ(CO Ⅱ, CO Ⅲ), ND4, ND4L, ND5, ATPase 6和Cyt  b 等蛋白编码基因以ATG 为起始密码子, CO Ⅰ和ATPase 8以GTG 为起始密码子, ND1, ND2, ND6和ND3则分别以ATA 和ATT 为起始密码子, 这与其他蛇亚目种类对起始密码子的选择无显著差异. 终止密码子包括4个TAA 和4个AGG, 某些蛋白编码基因如ND1, ND2, CO Ⅲ, ND3, Cyt b 等则使用了非标准终止密码子, 即不完全终止密码子T.
美姑脊蛇线粒体基因组13个蛋白编码基因的密码子使用情况如表  4. 对于氨基酸四重简并的位点, 碱基C 出现的频率最高, 其次是T 或A. 而两重简并的位点, 碱基C 出现的频率稍高于T. 除了氨基酸Arg 密码子外, 碱基G 在其他密码子第三位上都是最
少的, 这种低G 偏倚现象在脊椎动物中很普遍. 在美姑脊蛇蛋白编码基因中, C 的含量最高, 平均为32.8%;
G 的含量最低, 平均仅为13.7%, 蛋白质编码基因密码子的第一、第二和第三位点的A + T 的平均百分含量(50.7%, 57.8%, 52.2%)均高于G + C 含量, 其中第二位点占的比例最高(57.8%). 这一点与线粒体基因组有高AT 偏向性的特征是吻合的.  (ⅱ) RNA 基因.  美姑脊蛇含有的22个tRNA 基因, 可以满足线粒体蛋白质翻译中所有密码子的需要, 其长度为49 bp (tRNA-Cys)~73 bp (tRNA-Leu, tRNA-Asn). 其中tRNA-Glu, Ala, Asn, Cys, Tyr, Ser, Gln, Pro 由L-链编码, 其余由H-链编码. H-链编码的tRNA 基因散布于蛋白质基因和rRNA 基因之间, 相邻基因间隔1~30个碱基或紧密相连, 甚至也发生重叠. 一般说来, 脊椎动物的线粒体tRNA-Met, His, Leu(CUN)等序列保守性最高, 可达80%以上, 而
表3  美姑脊蛇线粒体蛋白编码基因
基因/区域
起始位置
终止位置
间隔(+)/重叠(−)
起始密码子
终止密码子
序列长/bp
氨基酸(aa)
编码链a)
ND1 2545 3493  ATA T++ 949 318 H ND2 4826 5849  ATA T++ 1024 341 H COI 6210 7810 +1 GTG AGG 1601 533 H COII 7936 8622 +4 ATG AGG 687 229 H ATPase 8 8686 8853  GTG TAA 168 56
H
ATPase 6 8844 9518 −10 ATG TAA 675 225 H COIII 9518 10301 −1
ATG T++ 783 261 H ND3 10362 10704  ATT T++ 343 114 H ND4L 10769 11059  ATG TAA 291 97 H ND4 11059 12411 −1
ATG AGG 1353 451 H ND5 12604 14352 +1 ATG TAA 1749 583 H ND6 14348 14851 −5
ATA AGG 504 167 L CYT b  14924 16040 +1 ATG T++ 1117 372
H
a) H, 重链; L, 轻链
表4  美姑脊蛇蛋白质编码基因碱基使用a)
物种
基因 T% C% A% G% (A+T)1% (A+T)2% (A+T)3% ND1 25.4 32.0 29.1 13.5 52.4 55.8 55.3 ND2 22.0 35.6 29.7 12.7 57.1 51.3 46.6 CO1 25.7 30.0 25.3 18.9 49.8 55.3 48.2 CO2 22.2 32.7 28.3 16.7 41.7 59.9 50.2 ATP8 19.6 31.0 35.1 14.3 50.0 53.6 60.7 ATP6 25.2 35.5 27.4 11.9 48.5 59.2 50.4 CO3 25.4 31.8 24.7 18.1 45.8 59.7 44.8 ND3 23.0 35.6 26.2 15.2 42.6 48.2 57.1 ND4L 24.8 32.1 28.3 14.8 52.6 62.9 43.8 ND4 23.1 34.8 29.1 12.9 46.1 51.9 58.7 ND5 24.0 33.4 30.7 11.9 57.9 58.0 48.0 ND6 17.0 28.1 51.3 3.7 61.3 76.5 67.1 美姑脊蛇Achalinus meiguensis
Cyt b
26.4 33.8 26.8 13.1 52.9 58.7 47.9 平均
23.4 32.8 30.2 13.7 50.7 57.8 52.2
a) (A + T) 1%, (A + T) 2%, (A + T) 3%分别指蛋白质编码基因密码子的第一、第二和第三位点的A +T 的百分含量
tRNA-Ser (AGY)的变异最大, 同源性只有54%~64%. Mt-tRNA可形成典型的三叶草形二级结构, 一般含有7个碱基的氨基酸接受臂, 5个碱基的TΨC臂, 反密码子臂及4个碱基的双氢尿嘧啶(DHU)臂, 各臂碱基对存在高比例的错配.
美姑脊蛇线粒体基因组包含了两条rRNA基因, 分别是16S rRNA和12S rRNA. 其位置与其他蛇类相同. 16S rRNA和12S rRNA基因序列分别为1486和933 bp, 这两个基因A + T含量分别为56.6%和51.1%, 均高于各自G + C含量(43.4%和48.9%).
(ⅲ) 非编码区.  至今为止, 蛇亚目中已测的mtDNA全序列, 除Leptotyphlops dulcis, Ramphotyph- lops braminus和Typhlops murius有一个控制区外, 其余种类都有两个完全复制的控制区[5]. 美姑脊蛇也有两个控制区, 占全序列的12.4%, 这个区域富含AT (A + T为62.5%). 由于缺少编码序列所具有的碱基约束性, 控制区被认为是线粒体基因组中变异最大的部分[46]. 但是在所有脊椎动物线粒体控制区也存在若干共同的保守区, 如重链复制起点以及重链和轻链转录的驱动子[47]. 控制区通常分为3个区域: 终止序列区、中央序列区和保守序列区[48]. 这3个区域中, 碱基T 的含量变化不大, 保持在35.4%~37.8%之间; 碱基A 的平均含量在终止序列区(30.2%)和保守序列区(28.9%)稍高于中央保守区(27.2%); 而碱基G的含量最低, 在7.6%~11.8%之间, 这与整个线粒体基因组序列相似.
2.2  系统发育关系
最大简约分析产生了1个最简约的系统发育树, 树长为52912, CI = 0.3203, RI = 0.4361. 图2显示了50%主要合意树. 简约分析支持蛇亚目为一个单系, 支持率为100%. 钩盲蛇和细盲蛇形成姊妹, 支持率为100%, 处于基部. 其他蛇类可以分为两个主要分支: 分支Ⅰ包括美姑脊蛇(Achalinus meiguensis)、瘰鳞蛇(Acrochordus granulatus)、黄水蚺(Eunectes notaeus)、球蟒(Python regius)、闪鳞蛇(Xenopeltis unicolor)、红尾筒蛇(Cylindrophis ruffus)和巨蚺(Boa constrictor), 支持率仅为43%; 分支Ⅱ包括玉米蛇(Pantherophis guttatus)、Pantherophis slowinskii、半棱链蛇(Dinodon semicarinatus)、冲绳烙铁头(Ovophis okinavensis)、铅水蛇(Enhydris plumbea)、印度眼镜蛇(Naja naja)、尖吻蝮(Deinagkistrodon acutus)和食鱼蝮(Agkistrodon piscivorus), 支持率为100%. 在分支Ⅰ中, 黄水蚺和巨蚺有最近共同祖先(BP = 100%),
它们与由闪鳞蛇、球蟒和红尾筒蛇组成的单系(BP
= 80%)构成姊妹, 支持率为98%. 这个类再和由
美姑脊蛇与瘰鳞蛇组成的亚枝互为姊妹关系. 在
分支Ⅱ中, 蝰科形成一单系, 与由眼镜蛇科和游蛇
科的种类构成姊妹. 游蛇科铅水蛇与眼镜蛇科
的印度眼镜蛇有最近的共同祖先, 但支持率仅为
49%.
Bayesian分析方法所揭示的蛇类系统发育关系
与最大简约分析结果相似(图3). 蛇亚目构成一个单
系, 得到高的后验概率支持(PP = 1.00). 不同于最
大简约树, 基于形态学划分的盲蛇下目 (Scoleco-
phidia)、原蛇下目(Henophidia)和真蛇下目(Caeno-
phidia)[49,50]各自形成单系, 瘰鳞蛇处于真蛇下目的
基部, 其次是美姑脊蛇, 它与蝰科、眼镜蛇科和游蛇
科组成的单系互为姊妹, 后验概率为1.00. 无论
是最大简约分析, 还是Bayesian分析, 两种系统发育
分析方法都较好地解决了蛇类下目间的系统发育关
系.
2.3  分子钟估算
依据系统发育统计检验结果(表5), 我们选择
Bayesian树来进行分歧时间的估算. 如材料和方法中
所阐述, 树中分支节点的分歧时间的估算是在松散
分子钟约束下, 采用Bayesian方法分析完成的. 估算
出蛇亚目的分化时间在128.03 Mya, 其95%置信区
间为124.32~131.60 (C2, 表6, 图4). 原蛇下目与真
蛇下目的分化时间为109.50 Mya, 其95%置信区间
为106.07~113.71 (C3, 表6, 图4); 而瘰鳞蛇与其他
真蛇下目种类的分化时间为106.18 Mya, 其95%置
信区间为103.17~110.22 (C4, 表6, 图4); 美姑脊蛇
的分化时间在103 Mya (C5, 表6, 图4); 蝰科与眼镜
蛇科和游蛇科构成的单系的分化时间发生在96.06
Mya (C6, 表6, 图4), 这与蝰科的化石记录一致.
3  讨论
3.1  蛇类线粒体基因组进化特点
在已测出全序的17种蛇中, 除了属于盲蛇下目
Leptotyphlops dulcis, Ramphotyphlops braminus和
Typhlops murius有单一的控制区, 其他种类都具有双
控制区, 并且控制区加倍事件发生约为128 Mya
(95%置信区间: 124.32~131.60). 这与Dong等人[18]的
观点不一致. 由于控制区加倍事件不仅仅在蛇亚目
1255