基于共同前缀位置的哺乳动物mtDNA序列系统发育分析

Phylogenetic Analysis of Mammalian mtDNA Sequences Based on Common Prefix Positions

摘要: 本文基于共同前缀标识符，通过环形前缀树法构建DNA序列的环形前缀树，得到能够表征序列信息的前缀集，去除生物序列长度对序列间距离计算的影响的同时并考虑序列对的共同前缀位置等因素，创建了一种新的非比对方法——基于共同前缀的位置差法，选取含有80条哺乳动物的mtDNA序列的数据集作为测试对象，运用邻接(Neighbor-Joining)法构建系统发育树，快速且准确地得到符合生物学分类的结果。

Abstract: This paper is based on the common prefix identifier. A circular prefix tree of DNA sequences is constructed by the circular prefix tree method to obtain a prefix set that can characterize sequence information, which removed the influence of biological sequence length on the calculation of the distance between sequences, and considered the common prefix position of sequence pairs. And other factors, created a new non-comparison method—the position difference method based on the common prefix. The data set containing 80 mammalian mtDNA sequences was selected as the test object, and the phylogenetic tree was constructed using the Neighbor-Joining method, which could fast and accurately get the results that meet the biological classification.

1. 引言

20世纪50年代以来，随着分子生物学和分子遗传学的兴起，生物学家根据各种生物间在分子水平上的进化关系，建立分子进化的系统发育树(Phylogenetic Tree)，估测物种间的亲缘关系，直观地阐明生物间的进化历程。

系统发育树重建传统上取决于多个或成对的序列比对。但是，使用基于对对的方法分析大型数据集时会遇到各种限制。高等真核生物的全基因组比对可能超过计算内存 [1]。因为存在基因组重排和重复的组合等因素使得不能进行整个基因组的比对，并且过多的序列信息增加了计算生物学中基因组比较的计算量和时间要求。因此，必须在初始步骤中确定正在研究的基因组的可比同源片段，并且迫切需要更快的序列分析算法。为此，提出了非比对方法克服基于比对方法的局限性 [2] [3] [4]。序列比较的非比对方法可以定义为在算法应用的任何步骤中量化序列相似性/相异性的任何方法，非比对方法不使用或产生比对(残基—残基对应的分配)。它们算法复杂度低且可以抵抗改组和重组等比对方法不能处理的事件。

非比对方法可大致分为两组：基于定义长度子序列频率的方法(基于词的方法)和评估全长序列之间信息内容的方法(基于信息理论的方法)。还有一些方法不能分类在任何一组中，包括那些基于匹配词长度的方法，混沌游戏表示，迭代映射，以及DNA序列的图形表示，它以定量的方式捕获基本组成的实质和序列的分布。本文是通过一种新的非比对方法——基于共同前缀的位置差算法(CPPA算法：Common Prefix Position-difference Algorithm)构建80中胎盘类哺乳动物的系统发育树，能够快速而准确地得到聚类结果，并且优于分布向量法 [5] 的结果。

2. 材料与方法

2.1. 材料

2020年10月22日从NCBI的Genebank数据库中下载的80条有胎盘类(Placentalia)哺乳动物的线粒体完整基因组，基因组的平均长度为1.65 × 10³碱基对。因为它是在大多数多细胞生物中从母亲那里继承而来的，不是高度保守，而且突变率不是很高，所以它对于研究基于线粒体基因组的进化关系非常有效。本文将Bo等人(2011) [5] 中80个物种中含有非ACGT的序列进行了同属或同科替换，保证本方法的有效性，替换后的数据集详细信息见附录。

2.1.1. 前缀标识符

在有限字母 $\Sigma =\left\{{\sigma }_1,{\sigma }_2,\cdots ,{\sigma }_t\right\}$ 上指定的长度为m符号序列 $P=P\langle 1\rangle P\langle 2\rangle \cdots P\langle m\rangle$，其中 $P\langle i\rangle \in \Sigma$ 是序列中的第i个符号，并且被称为位于第i个位置。 $P=P\langle i:j\rangle$ 表示从P的位置i开始到j的位置的字符子串。与 $\Sigma$ 上的字符串S相关联的前缀树是t叉树 $T_P=\left(N,n_0,S_P\right)$， $W_p\left(n\right)$ 是 $S_P$ 的最小长度唯一子串，其最左边的符号位于S的位置j。 $W_p\left(n\right)=S\langle j:k\rangle$ 在S中仅发生一次，而 $S\langle j:k-1\rangle$ 在S中至少发生两次。则称 $I\left(j\right)=S\langle j:k\rangle$ 为与位置j相关的前缀标识符 [6]。

本文的CPPA算法的前缀树模型是以Weiner Peter的 [6] 前缀标识符理论为基础提出的，并通过将张欣 [7] 的相异度量模型中线性前缀树算法改进得到的。

2.1.2. 环形生物序列前缀集

哺乳动物的线粒体基因组是长度约为1.65 × 10³碱基对的双链DNA环 [8]。但是从NCBI (National Center for Biotechnology Information Search database)网站中获取的序列是线性的，因此将序列首尾进行拼接，还原为环形序列，再基于前缀标识符，将线性前缀树算法 [7] 改进为环形前缀树算法构建前缀树，提取前缀集。

对任一序列 $S=S_1,S_2,S_3,\cdots ,S_{n-1},S_n,S_i\in \left\{A,C,G,T\right\},1\le i\le n,n=length\left(S\right)$，将序列S足够长的前m个碱基拼接在S末端，构建新序列 $S^{\prime }=S_1,S_2,S_3,\cdots ,S_{n-1},S_n,S_{n+1},S_{n+2},\cdots ,S_{n+m},m\ll n,n+m=length\left(S^{\prime }\right)$，对 $\forall i,1\le i\le n$， $\exists$ 唯一前缀 $S_i\_prefix=S_i,S_{i+1},\cdots ,S_{i+l_i},l_i=length\left(S_i\_prefix\right)$，s.t.当 $l_j=l_i\left(j\ne i\right)$ 时，满足 $S_j\_prefix\ne S_i\_prefix$，且至少存在一个 $S_k\_prefix$，满足 $S_k\_prefix=S_i,S_{i+1},\cdots ,S_{i+l_i-1},l_k=l_i-1$。由 $S_i\_prefix,1\le i\le n$ 构成的集合 $I\left(S\right)$ 称为序列S的前缀集。

2.1.3. 序列间基于共同前缀的位置法

DNA序列 $S_i,S_j\left(i\ne j\right)$ 的长度分别为 $n_1,n_2$。首先应用环形前缀树法构建两条序列的前缀树，得到前缀集 $I\left(S_i\right),I\left(S_j\right)$。再求出两条序列共同前缀集 $I\left(S_i\right)\cap I\left(S_j\right)$ 的元素在 $S_i,S_j$ 中的标准化位置，去除序列长度对计算序列间距离的影响。第t个共同前缀在 $S_i$ 中的标准化位置定义为：

$P_i^t=P{r}_i^t/n_i$ (1)

$P{r}_i^t$ 为第t个共同前缀在序列 $S_i$ 中的实际位置。

假设 $S_i,S_j$ 共有m个共同前缀，所有共同前缀在 $S_i,S_j$ 中标准化位置构成如下位置向量： $P_i=\left(P_i^1,P_i^2,\cdots ,P_i^m\right),P_j=\left(P_j^1,P_j^2,\cdots ,P_j^m\right)$。

两条序列的共同前缀的位置做差，得到位置差向量。考虑到哺乳动物的线粒体基因组是环形的，第t个共同前缀都的位置在序列中的位置如图1所示。

Figure 1. The relative position of the common prefix in the sequence pair

图1. 共同前缀在序列对中的相对位置

$S_i,S_j$ 的第t个共同前缀的位置一同放在一个长度为1的环形DNA中，则距离差会出现两种结果红色弧线段的长度或绿色弧线段的长度，因为共同前缀位置之间的距离越小，则两序列越相似，因此取两种位置差中较小的数值作为所求的位置差。因此，两条序列 $S_i,S_j$ 的位置差向量 $position\_diff$ 定义为：

$position\_diff=\left(d_1,d_2,\cdots d_m\right)$，其中 $d_t=\min \left\{\left|P_i^t-P_j^t\right|,1-\left|P_i^t-P_j^t\right|\right\},t=1,2,\cdots ,m$ (2)

本文采用基于共同前缀的位置差算法来计算进化距离。因为序列间共同前缀越多，序列越相似，所以将序列对 $S_i,S_j$ 间的相似性距离定义为：

$dis\left(S_i,S_j\right)=\frac{1}{m}\left|position\_diff\right|=\frac{1}{m}\sqrt{{\displaystyle {\sum }_{t=1}^m{\left(d_t\right)}^2}}=\frac{1}{m}\sqrt{{\displaystyle {\sum }_{t=1}^m{\left(\min \left\{\left|P_i^t-P_j^t\right|,1-\left|P_i^t-P_j^t\right|\right\}\right)}^2}}$ (3)

获得了用于系统发育分析的所有遗传序列之间的距离构成的距离矩阵后，使用Matlab2013a通过邻接法(NJ)绘制系统发育树 [9]。

3. 结果

3.1. 系统发育分析

Figure 2. The Phylogenetic tree of 80 mammals by CPPA

图2. 由CPPA法得到的80种哺乳动物的系统发育树

80条哺乳动物线粒体全基因组mtDNA序列数据集中，全部是有胎盘类，其中非洲兽总目(Afrotheria) 7条序列，北方兽类(Boreoeutheria) 73条序列。北方兽类包含灵长总目(Euarchontoglires) 19条序列，劳亚兽总目(Laurasiatheria) 54条序列 [10] [11] [12]。

从整体上看，图2中CPPA法得到的系统发育树中，非洲兽总目正确聚类，灵长总目未聚为一个分支，劳亚兽总目的序列来自6个目，5个目的52条序列聚在一起；其它目的2条序列与灵长总目的啮齿目聚到一起。

数据集中，劳亚兽总目的序列具体是：食肉目(carnivora) 27条序列、鲸偶蹄目14条序列、奇蹄目(Perissodactyla) 6条序列、翼手目(Chiroptera) 4条序列、鼩形目(Soricomorpha) 1条序列和猬形目(Erinaceomorpha) 2条序列，6个目都分别聚成一个分支。食肉目的犬形亚目(Caniformia)与猫型亚目(Feliformia)分别聚为一个分支 [13]，鲸偶蹄目中的3 条序列鲸目(cetacea)与偶蹄目(Artiodactyla)的6条反刍亚目(Ruminantia)和1条猪次目(Suina)正确聚类，并与交错在一起，支持鲸偶蹄目的分类 [14]。

灵长总目的序列来自3个目，灵长目(Primates)、兔形目(Lagomorpha)和啮齿目(Rodentia)。灵长总目虽然没有形成一个分支，灵长目和兔形目各自形成自己的分支，且聚类正确 [15]。啮齿目的6条序列来自3个亚目，鼠形亚目(Myomorpha) 4条序列，松鼠亚目(Sciuromorpha) 3条序列和豪猪亚目(Hystricomorpha) 1条序列 [16]。每个亚目各自形成单独的分支。

在科的层级上，含2条及以上的科共16个，均正确聚类。

3.2. CPPA法与分布向量法系统发育树比较

Bo Zhao和Rong L. He等人的分布向量法构建的系统发育树 [5]，仅考虑了生物之间的聚类关系，因此本文在聚类的角度上进行比较。

一、目的层级上，分布向量法的聚类结果中27条食肉目序列被分为两部分，且猫型亚目与犬形亚目均未聚类正确，14条鲸偶蹄目序列被分为鲸目和偶蹄目两部分，不支持鲸偶蹄目的分类，CPPA法的食肉目与鲸偶蹄目均聚类正确，6条啮齿目序列的聚类CPPA法略差于分布向量法。

二、科的层级上，分布向量法的聚类结果中4条犀牛科(Rhinocerotidae)序列、3条海狮科(Otariidae)序列、7条犬科(Canidae)、3条须鲸科(Balaenopteridae)序列、5条人科(Hominidae)序列、4条猕猴科(Cercopithecidae)序列和2条兔科(Leporidae)序列7个科未聚类正确，CPPA法均聚类正确。

以上结果说明由80条哺乳动物线粒体全基因组mtDNA序列数据集构建的系统发育树中，CPPA法的聚类结果远好于分布向量法的聚类结果。

4. 结论

本文基于共同前缀标识符提出的CPPA法，首先通过环形前缀树算法获得每条序列的前缀集，挖掘到能表征序列的定量代表——共同前缀在序列对中所处的标准化位置，去除了序列长度对序列间距离的影响，并充分考虑了共同序列个数在计算序列间相似性距离的重要性，应用Matlab运行算法速度快，且聚类结果准确，说明CPPA法是一种好的非比对方法。通过CPPA法得到的系统发育树中物种间的进化关系还不够好，可以考虑影响序列间距离的其它因素进行改进。

致谢

感谢Bo等人为本研究提供了80条哺乳动物线粒体序列的NCBI检索号。同时，也要对本文中引用文章的作者表示衷心的感谢。

附录

Table 1. Description of mitochondrial genome of 80 mammalians

表1. 80条哺乳动物线粒体序列