1. 引言
鸟类起源于侏罗纪时代的兽脚亚目，是恐龙的唯一现存后裔。同时鸟类也是四足类脊椎动物中物种最丰富的一个纲，目前大约有10,000多个鸟类物种存在于地球的不同环境中。在过去数十年中，人们对鸟类之间物种进化关系的研究有了巨大的进展和重大的突破 [1] 。研究认为，大约在6600万年前小行星撞击地球或是火山爆发而发生了著名的白垩纪–第三纪灭绝事件(Artistic depiction of asteroidal impact)，这次大范围的物种大灭绝事件造成了恐龙的灭亡与哺乳动物的兴起，同时也为鸟类新物种的爆发创造了条件且提供了新的生存环境。在白垩纪时期，所有现代鸟的祖先发生了第一次分化，出现了一些飞行能力不太强的鸟类祖先，如鸵鸟的祖先。这次分化之后形成的一个姐妹支又出现了一次分化，而产生了现在的鸡、鸭等的祖先 [2] 。但是，到目前为止仍有许多关于鸟类分类和进化的问题未得到解决。这是因为鸟类在早期进化过程中经历了那次辐射性的物种大爆发事件，而爆发时间很短使得现存鸟类之间的进化关系难以区分。虽然长期以来研究者众多，但是在2010年前2/3以上的目、科分类级别关系较为混乱。同时，在大爆炸事件前后1000~1500万年内形成的新鸟小纲(Neoaves)是占据现存鸟类物种中绝大比例的一支 [3] [4] [5] [6] ，但目前为止仍然无法确定新鸟小纲中一些早期分支的分类情况和具体的进化规律。因此，鸟类进化树的构建和分析有助于了解鸟类的进化历史，对研究鸟类的进化过程和进化规律具有重要的生物学意义，研究中所用的方法对其它类群进化问题的研究同样具有重要的参考价值。
近年来，随着测序技术水平的快速提高，人类获得了更加丰富的鸟类基因组数据，鸟类分子系统发生学的研究也因此进入了快车道。2014年，Jarvis等人 [3] 基于48种鸟类全基因组数据，利用系统发育基因组学的方法研究了鸟类的起源与进化。他们重建了一个可信度极高的鸟类进化树(可信度90%以上), 证实了此前备受争议的姐妹群，重新确定了鸟类的进化关系。在这个树中，鸟类物种被分成古颚总目(Palaeognathae)和今颚总目(Neognathae)，古颚总目是其余所有现存鸟类的姐妹群，今颚总目由鸡雁小纲(Galloanseres)和新鸟小纲组成。Jarvis等发现，占现存鸟类95%物种的新鸟小纲可以分成两大分支，分别被命名为雀类(Passerea)和鸽类(Columbea)，这两大分支分别独立趋同进化出了各自的陆生和水生鸟类物种。2015年，Prum等人 [7] 利用靶定下一代DNA测序技术，完成了对198种鸟类物种的测序，构建了Zenodo数据库，从而建立了更加全面的系统发育树。在这个系统发育关系中，鸟类物种同样被分成古颚总目和今颚总目，鸡雁小纲和新鸟小纲组成今颚总目。但是在单源的新鸟小纲内，Prum等确定了五个主要分支，分别命名为Strisores, Columbaves, Gruiformes, Aequorlitornithes和Inopinaves，且对它们进行了更为详细的分类。建树过程中运用的时间参数以得到的地质化石较准时间为准。
尽管鸟类的基因组信息越来越全面，但在构建鸟类进化树中也遇到了诸多困难。其中一个主要原因是鸟类物种大爆发的时间很短，导致鸟类基因组之间存在大量的基因不完全谱系分选(Incomplete Lineage Sorting of Genes) [8] [9] 现象，从而给进化树构建带来不确定性：由单个基因构建出来的进化树很难与由全基因组构建的物种树一致。另外，全基因组上某些区域会有较多的噪音，这些噪音会影响物种进化树的推断 [3] 。因此，研究者需要从全基因组序列中仔细挑选出适合构树的高质量区域。
本文选取两侧为可变区域的保守序列数据，基于k联体的非联配方法，构建了一个由47个鸟类物种组成的鸟类基因组进化树。所谓k联体的非联配方法就是生物序列数据频次统计算法，也称为k-mer频次统计算法。遗传语言-DNA序列的信息很大程度上反映在其k-mer频数中，k-mer频次统计信息可以用来揭示生物序列中各种子序列的分布规律，它是一种衡量序列相似性的重要工具，这为研究鸟类物种亲缘关系提供了新的角度 [10] 。我们分别统计了9联体到14联体的频率，计算了不同基因组之间的距离，构成一个多维的距离矩阵，利用邻接法构建出能反映47个鸟类物种之间进化关系的进化树，我们称之为IMU(Inner Mongolia University的简称)进化树。在计算过程中发现，随着k的增加，IMU进化树越来越稳定。当k增大到12时，进化树已经不再随k值的增加而变化，由此说明12联体对于鸟类基因组具有物种特征性，这也意味着12联体频数是能够较好描述进化的动力学变量。同时通过观察IMU进化树的拓扑结构，发现5种古老鸟类处于树的根部，其余42种现代鸟类形成了Neoaves1(陆鸟)，Neoaves2 (陆鸟和水鸟的混合)，Neoaves3(包含猛禽的陆鸟)，Neoaves4(水鸟)等4支。最后，将IMU进化树与Prum和Jarvis等人构建的进化树进行了比较，分析了鸟类的进化和分类情况。
2. 方法
2.1. 数据来源
本文的数据来源于Zenodo(https://zenodo.org/record/28343#.VzCY5PmEBFn)，由Prum等 [7] 人于2015年所构建，是由靶定下一代DNA测序方法得到的序列数据，其中包含可以代表鸟类所有主要分支的198个鸟类和2个鳄目类(crocodilians)共200个物种的数据。本文通过比较由Prum等人研究的198种鸟类与Jarvis等人研究的48种鸟类属的一致性，从而将两种不同的数据名称进行了对应，最终得到了能够代表不同属的47种鸟类，如表1所示。例如，Prum等人研究的198种鸟中只有一种企鹅，为洪堡企鹅; 而Jarvis等人研究的48种鸟中同为企鹅属的有两种企鹅，分别为帝企鹅和洪堡企鹅; 故本文应用洪堡企鹅的数据。因此，由英文常用名和Prum等人命名的鸟类一定对应相同的属，但对应的具体鸟类可能不同。为方便与不同进化树进行比较及分析，本文均用英文常用名表述鸟类。具体地，每一种鸟类的DNA序列包含由靶定基因组中约259个基因所产生的394个序列片段，这些片段多数是两侧为可变区域的保守序列，其中每个基因平均约有1524个碱基，最后每种鸟类的DNA序列约有390,000个碱基。
在建树过程中，时间参数依据Prum等得到的地质化石较准时间，即文献 [7] 附录中图7所示的时间，由此确定了每一物种的产生时间。
2.2. 本文研究的47种鸟类名称和记号
为方便与Prum和Jarvis构建的进化树进行比较，本文对47种鸟类的英文常用名、记号、对应于英文常用名的中文常用名 [2] 、以及Prum对鸟类的命名进行了一一对应(如表1所示)。字母记号的规则为：前面的大写字母表示多层次分类，之后的两位数字对应物种在其分类中的编号，最后一位小写英文字母代表编号后更加细致的分类，若没有则不写。
因此，198种鸟类的字母记号具体如下：P01-P09表示Palaeognathae, G01-G16表示Galloanserae,
Table 1. The name and mark of 47 avians
表1. 47种鸟类的名称和记号
Continued
NS01-NS13表示Neoaves-Strisores，NC01-NC16表示Neoaves-Columbaves，NG01-NG09表示Neoaves-Gruiforms，NA01-NA43表示Neoaves-Aequorlitornithes，NIM01-NIM10表示与Neoaves- In-opinaves-Accipitriformes相关的混合鸟类, NIC01-NIC26表示Neoaves-Inopinaves-Coraciimorphae, NIA01-NIA56表示Neoaves-Inopinaves-Ausralaves, NIA13p-NIA56p中的字母p表示subclass passeriformes。通过分析新命名的字母记号可以直接了解物种分类的从属问题，方便快速判断所构进化树的准确性。
2.3. 数据预处理
k联体是指在核苷酸序列中连续排列的k个碱基。由于下载得到的测序数据中含有碱基不确定的位点(占序列总长的0.56%~5.76%)，如强配对碱基(S)，次黄嘌呤碱基(I)，含氨基碱基(M)等，为了保证序列长度不变，同时避免引入虚假的k联体，将这些位点用“−”表示。这样，经过预处理后的序列数据则是由腺嘌呤(A)、胸腺嘧啶(T)、胞嘧啶(C)、鸟嘌呤(G)四种碱基及“−”组成。
本文不考虑含有“−”的k联体，故仅包含字母的k联体最多只有4k种。我们从经预处理后的47种鸟类的基因组序列资料出发，统计k联体(k = 9 − 14)的频数及其相应频率。据此构建进化树。
2.4. 核苷酸k联体的非联配算法
我们采用k联体的非联配算法 [11] (KNA算法: k-mer Non-aligned Algorithm)来计算进化距离。对于任意两保守序列
和
设其核苷酸
型k联体的频率分别为
和
，则两序列之间的距离定义为:
(1)
根据(1)式，利用k联体频率可计算出47种鸟类两两之间的进化距离值，从而得到一个47 × 47维的距离矩阵。
2.5. 用NJ算法构建鸟类进化树
邻接法(Neighbor-Joining,简称为NJ)是由Saitou和Nei [12] 于1987年提出的，是一种基于距离矩阵的快速聚类方法。该算法不需要关于分子钟的假设，不考虑任何优化标准。其基本思想是进行类的合并时不仅要求待合并的类是相近的，而且要求待合并的类远离其他的类，从而通过对完全没有解析出的星状进化树进行分解，来不断改善星状树。在递推过程中NJ算法不仅考虑一对对的距离，而且考虑了整个树中距离之和的极小化。
特别是，在NJ算法的递推过程中要对一对物种(如物种1和物种2)的分支长度求和，记为S12:
(2)
其中，N表示物种总数，d1k和d2k分别表示物种1和2与其它物种之间的距离值，d12表示物种1和2之间的距离值，dij表示除了物种1和2之外的其他物种两两之间的距离值。当S12为求得的所有Sij中最小值时，表示1和2即为邻近物种。很显然，式(2)中物种1和物种2是对称的，即两物种替换后S12值不变。由此说明一个姐妹支上的两个邻近物种是可以上下互换的。
3. 结果
3.1. 频数统计结果
对47种鸟类基因组中的靶定测序数据，我们计算了从9-mer到14-mer的k-mer的频数分布，如表2所示。结果表明，对于多联体而言，频数为2 − 4和频数大于4的k-mer占k-mer总体的比例随着k的增大而降低，而频次为1的k-mer随着k的增大而增大。
例如，对于12联体，我们计算了暴风鹱(Northern fulmar)、企鹅(penguin)、红头美洲鹫(Turkey vulture)、鸽子(Pigeon)、麝雉(Hoatzin)和非洲鸵鸟(Common ostrich)共6种鸟类频数大于4的12-mer的分布，如图1所示。其中，横坐标表示12-mer的频数，范围从5到20，纵坐标表示给定频数下12-mer的种类数。分析发现，基本满足随着频数的增加对应种类数减小的趋势。通过统计较高频次12-mer的分布可深入探究IMU树中各姐妹支形成的原因。
3.2. 利用距离矩阵构建IMU进化树
本文基于KNA算法，我们构建了47种鸟类的进化树，称为IMU(Inner Mongolia University的简称)树.对所得的6个k-mer (k = 9 − 14)进化树进行比较，可得出较合理的IMU树。经分析发现从9-mer进化树开始,5种古鸟就位于树的根部，与已知的进化关系相符; 紧邻其上有一支NI，树的最上面有一支NA，这些分枝已经达到稳定。当k增加到12时，得到12联体时的NJ进化树，就和Jarvis和Prum等构建的进化树基本一致。当k的值由12增加到13、14时，发现进化树已不再随k的增加而发生任何变化。由此说明，12联体对于鸟类基因组具有物种特征性。我们得到的IMU进化树如图2所示，以“P/G”和“N”
Table 2. The proportion of k-mer with different frequency
表2. 不同频数k-mer的占比
Figure 1. The distribution of 12-mer frequency insix kinds of avian DNA sequences
图1. 6种鸟类DNA序列的12联体频数分布
开头分别代表了古鸟和新鸟，可分为五部分,分别为古老鸟，新鸟小纲的1、2、3、4共五支，各自代表了古老鸟、陆鸟、水陆混合鸟类、陆鸟和水鸟，图中分别以紫、黑、蓝、红、绿色的线表示。其中，古老鸟类(Palae/Gallo)、猛禽(Neoaves3中的NIM02a、NIA04、NIM07a、NIM08a、NIA01、NIM09)、水鸟(Neoaves4)的分类非常明确。此外，物种形成时间取自文献 [7] 中的补充图7。
4. 讨论
4.1. Prum和Jarvis等构建进化树的比较
由Prum等人构建的进化树 [7] 和Jarvis等人构建的进化树 [3] 是目前最为可信的两棵鸟类生命之树。通过对这两棵进化树的拓扑结构进行比较，我们分析了它们之间的异同。
两树一致把鸟类分为古颚总目和今颚总目，且鸡雁小纲和新鸟小纲共同组成今颚总目。两树也同样认为蜂鸟科，雨燕科和夜鹰科形成姐妹分支，其拓扑结构为(蜂鸟科 + 雨燕科) + 夜鹰科) (三者都属于夜鹰目)。沙鸡目，拟鹑目和鸽形目形成姐妹分支，其拓扑结构为(沙鸡目 + 拟鹑目) + 鸽形目)。美洲火烈鸟(火烈鸟科)和凤头鷉(鸊鷉科)形成姐妹分支。日鳽目，鹲蒙目和潜鸟目形成姐妹分支，其拓扑结构为(日鳽目 + 鹲蒙目) + 潜鸟目)。对于鹃形目，蕉鹃目和鸨形目，两树都认为它们是近亲,然而物种的拓扑结构有些差异，Jarvis等人认为(鹃形目 + (蕉鹃目 + 鸨形目)，Prum等人认为(蕉鹃目 + (鹃形目 + 鸨形目)。
对于水生鸟类，两树一致认为鹱形目和企鹅目有密切关系。鸬鹚科，鹮科，鹭科和鹈鹕科形成连续的姐妹分支，其拓扑结构为(鸬鹚科 + 鹮科) + 鹭科) + 鹈鹕科)(四者都属于鹈形目)，且鹈形目与鹱形目-企鹅目进化支形成姐妹群。对于陆生鸟类，两树一致的认为美洲鹫科与鹰科是近亲，形成姐妹分支(两者同属于鹰形目)，在此我们将其称为陆鸟群1。鸮形目，鼠鸟目，鹃鴗目，咬鹃目，犀鸟目，䴕形目和佛法僧目形成连续的姐妹分支，其拓扑结构为(䴕形目 + 佛法僧目) + 犀鸟目) + 咬鹃目) + 鹃鴗目) + 鼠鸟目) + 鸮形目)，在此我们将其称为陆鸟群2。雀形目，鹦形目，隼形目和叫鹤目形成连续的姐妹分支，其拓扑结构为(雀形目 + 鹦形目) + 隼形目) + 叫鹤目)，在此我们将其称为陆鸟群3。在新鸟小纲中，两树一致的认为陆生鸟类的共同祖先是顶端食肉动物。
Figure 2. The Phylogenetic tree of 47 avians by KNA algorithm (IMU tree)
图2. 由KNA算法得到的47种鸟类的进化发生树(IMU树)
对于两树拓扑结构的差异，可归纳为以下四点：其一，两树虽一致的认为麝雉应该自成一目，但是Jarvis等人认为麝雉与鸻形目，鹤形目形成姐妹群；而Prum等人构建的进化树将其放置在陆鸟的基部，作为其余所有陆鸟的姐妹群，而且认为鹤形目应该自成一大分支，与所有的陆鸟和水鸟形成姐妹分支。其二，Jarvis等人认为火烈鸟科-鸊鷉科与鸽形目是近亲，且紧邻鸡雁小纲; Prum等人认为火烈鸟科-鸊鷉科与岸鸟(鸻形目等)是近亲，紧邻鹤形目，紧邻鸡雁小纲的是鹰形目。其三，Jarvis等人认为(蜂鸟科 + 雨燕科) + 夜鹰科)进化支与(鹃形目 + (蕉鹃目 + 鸨形目)进化支形成姐妹群;Prum等人认为(沙鸡目 + 拟鹑目) + 鸽形目)进化支与(蕉鹃目 + (鹃形目 + 鸨形目)进化支形成姐妹群。其四，在整个陆生鸟类这个大进化支内，Jarvis等人认为拓扑结构为(陆鸟群1 + 陆鸟群2) + 陆鸟群3)，而Prum等人认为拓扑结构为(陆鸟群1 + (陆鸟群2 + 陆鸟群3)。
4.2. IMU进化树与Prum和Jarvis等构建进化树的比较分析
Jarvis和Prum等人构建的进化树都是基于序列联配方法构造的联配进化树，而IMU进化树是基于非联配方法构造的。因此，IMU进化树为解析Jarvis和Prum等构建进化树之争提供了一个新的观察角度。通过三树比较，它们的拓扑结构基本一致，但存在一些差异，分析如下：
一，从进化树整体上看，三树一致把鸟类分为古颚总目和今颚总目，且鸡雁小纲和新鸟小纲共同组成今颚总目，鸡雁小纲由鸡形目和雁形目构成，古颚总目和鸡雁小纲位于进化树的基部。
二，在IMU进化树中，鹃形目与蜂鸟科-雨燕科进化支形成姐妹群，这与Jarvis等人构建的进化树和Prum等人构建的进化树都不同。在IMU进化树中，鸻形目和鹤形目形成姐妹群，这与Jarvis等人构建的进化树相同。然而，在IMU树中鸻形目-鹤形目进化支与火烈鸟科-鸊鷉科进化支形成姐妹分支，但在Prum等人构建的进化树，鸻形目与火烈鸟科-鸊鷉科进化支形成姐妹分支，鹤形目自成一大支。麝雉是一个较为特殊的物种，它在三树中的位置各不相同, IMU进化树认为麝雉与红冠蕉鹃(蕉鹃目)是近亲，Jarvis等人构建的进化树中麝雉与鸻形目，鹤形目形成姐妹分支，Prum等人构建的进化树将其放置在陆鸟的基部，与所有陆鸟形成姐妹群。
三，在IMU进化树中，仓鸮，红腿叫鹤，白尾海雕，白头海雕和游隼这些食肉类猛禽形成了连续的姐妹群，而另外两个进化树中分别位于陆鸟基部。同时，䴕形目，佛法僧目，雀形目，鼠鸟目和咬鹃目形成连续的陆鸟姐妹群，这个陆鸟姐妹群不同于在另外两个进化树中处于进化树的中下端，且紧邻鸡雁小纲，这个结果与Slack [13] 得到的相近。对于水生鸟类，三树的拓扑结构是一致的，新奇的是在IMU进化树中水鸟处于进化树的最顶端，其另外两树不同。
四，在IMU进化树中，鹦形目与夜鹰科形成了姐妹分支，日鳽目与拟鹑目-沙鸡目进化支形成姐妹群，这些关系都是独特且新奇的。特别需要指出的是鸠鸽科和美洲鹫科，在IMU树中，它们形成了姐妹分支，而在另外两树中，鸠鸽科与沙鸡科-拟鹑科进化支形成了姐妹群。
以上结果说明部分鸟类在IMU进化树中与在另外两树中的分类有差异，为了探究其原因，本文以鸠鸽科分类为例，从序列数据出发研究了其保守序列中频数大于4的12-mer类型(如表3所示)。在IMU树中鸽子与红头美洲鹫形成姐妹分支，而在另外两树中鸽子与黄喉沙鸡、褐拟鹑形成姐妹分支。对比分析发现，鸽子与红头美洲鹫的12-mer高频次类型具有高度相似性(如表3所示，具有相同颜色、下划线等标记的12-mer相同)，这应该是导致在IMU树中鸠鸽科和美洲鹫科形成姐妹分支的主要原因。我们猜测这些重复性极高的12联体同时出现在两类鸟的基因组中可能说明它们有某些相似的功能。这是一个鸟的分类研究中值得注意的问题。
5. 结论
本文利用靶定下一代DNA测序技术得到的394个保守片段的DNA序列数据，基于k联体的非联配算法构建了47种鸟类的k(k = 9 − 14)联体NJ进化树，我们称之为IMU进化树。结果显示随着k增加到
Table 3. The 12-mers with a frequency greater than 4 in 4 avians
表3. 4种鸟类中频数大于的4的12联体
注：表中用相同颜色和下划线等标记的两个12-mer是相同的，通过比较鸽子与红头美洲鹫的12-mer高频次类型，发现它们具有高度相似性。
12时进化树达到稳定，由此表明12联体对于鸟类基因组具有物种特征性，这也意味着12联体频数是能够较好描述进化的动力学变量。同时通过观察IMU进化树的拓扑结构，发现其明确的分为5支。最后，比较了IMU进化树与Prum和Jarvis等构建进化树的异同，分析了鸟类的进化和分类情况，发现三树的各个姐妹分支基本一致。但处于三个树根部的新鸟分类还有明显不同。这可能暗示发生物种大爆发的过程中尚存在某些有待探索的规律。如果核苷酸频数是好的动力学变量，我们就可以写出关于这些变量的时间进化方程，进而把基因组进化方程量子化 [14] 。这样就可以用量子随机性来解释新鸟的大爆发式的产生。这是一个有趣的有待进一步研究的问题。
致谢
感谢基金项目；国家自然科学基金资助项目(批准号：61462068和31106188)和内蒙古自治区自然科学基金资助项目(批准号：2014MS0103)对本论文的大力支持。感谢Prum等人为本研究提供了鸟类及其对应化石时间数据。同时，也要对本文中引用文章的作者表示衷心的感谢。