CVtree构建进化树算法简介
频数(frequency)的定义
我们可以使用的序列包括DNA序列和氨基酸序列，下文以氨基酸序列为例。
将\(K\)-string记为\(a_1 a_2 \cdots a_K\)。一段总长为\(L\)的序列一共可以得到\((L-K+1)\)个\(K\)-string。对于某一个\(K\)-string \(a_1 a_2 \cdots a_K\)而言，它的总数（或称频数，frequency）记为\(f(a_1 a_2 \cdots a_K)\)。于是某一个\(K\)-string出现的概率为
\(f(a_1 a_2 \cdots a_K) = \frac{f(a_1 a_2 \cdots a_K)}{(L-K+1)}\).
减去随机背景(random background)
这一步处理的目的是为了消除对进化不起作用(neutral)的随机变异的影响。
我们在计算取值为\(K\)的情况的这一步时，先计算出string长度为\((K-1)\)和\((K-2)\)时的概率。在马氏性假设下，序列第\(K\)位出现某个\(a_K\)和\(a_1\)无关，可以得到：
\(p(a_2 a_3 \cdots a_K \vert a_2 a_3 \cdots a_{K-1}) = p(a_1 a_2 \cdots a_K \vert a_1 a_2 \cdots a_{K-1})\).
上式展开整理后，我们可以估计\(K\)这一步时的概率：
\(p^0 (a_1 a_2 \cdots a_K) = \frac{p(a_1 a_2 \cdots a_{K-1}) p(a_2 a_3 \cdots a_{K-1})}{p(a_1 a_2 \cdots a_K)}\).
这里的数学处理的背景详见原文中所列的参考文献。
减去随机背景的这一步如下：
\[c^i(a_1 a_2 \cdots a_K) =\begin{cases} \frac{p(a_1 a_2 \cdots a_K)-p^0(a_1 a_2 \cdots a_K)}{p^0(a_1 a_2 \cdots a_K)}, & p^0 \neq 0, \\ 0, & p^0 = 0. \end{cases}\]这里\(i\)的取值从\(1\)到\(N=20^K\)，遍历了所有可能的\(K\)-string。
对于物种\(A\)，我们可以组合上面一步遍历\(i\)后的计算结果，得到一个\(20^K\)维向量
\(A = (a^1,a^2,\cdots,a^N)\),
对于物种\(B\)，同样可以得到向量
\(B = (b^1,b^2,\cdots,b^N)\).
计算两物种间的距离
这一步的处理是相当常规的。
首先我们计算两个向量之间的夹角的余弦值：
\(C(A,B)=\frac{A \cdot B}{\vert\vert A \vert\vert\, \vert\vert B \vert\vert}\).
然后做一个线性变换，将距离映入\([0,1]\)：
\(D(A,B)=\frac{1-C(A,B)}{2}\).
最后我们用邻接法(neighbor-joining method)构建进化树。
参考文献
Zhao Xu, Bailin Hao, “CVTree update: a newly designed phylogenetic study platform using composition vectors and whole genomes”, Nucleic Acids Res. published online on April 26, 2009. doi:10.1093/nar/gkp278.
声明: 本文采用 BY-NC-SA 协议进行授权. 转载请注明转自: CVtree构建进化树算法简介