Volume 14 Supplement 15

## Proceedings of the Eleventh Annual Research in Computational Molecular Biology (RECOMB) Satellite Workshop on Comparative Genomics

# Phase change for the accuracy of the median value in estimating divergence time

- Arash Jamshidpey
^{1}Email author and - David Sankoff
^{1}

**14(Suppl 15)**:S7

https://doi.org/10.1186/1471-2105-14-S15-S7

© Jamshidpey and Sankoff; licensee BioMed Central Ltd. 2013

**Published: **15 October 2013

## Abstract

We prove that for general models of random gene-order evolution of *k* ≥ 3 genomes, as the number of genes *n* goes to ∞, the median value approximates *k* times the divergence time if the number of rearrangements is less than *cn*/4 for any *c <* 1. For some *c** ≥ 1, if the number of rearrangements is greater than *c***n*/4, this approximation does not hold.

## Introduction

The iterative improvement of approximate solutions to the Steiner tree problem by optimizing one internal vertex at a time has a substantial history in the "small phylogeny" problem for parsimony-based phylogenetics, both at the sequence level [1] and the gene order level [2]. It has been generalized to iterative local subtree optimization methods such as "tree-window-hill" [3] and "disc covering" [4, 5]. Here we focus on the "median problem" for gene order where we estimate the location of a single point (the median) in a metric space given the location of the three or more points connected to the median by an edge of the tree. Given *k* ≥ 3 signed gene orders *G*_{1}, ..., *G*_{
k
} on a single chromosome or several chromosomes, and a metric *d* such as breakpoints [6], inversions [7], inversions and translocations [8], or double-cut-and-join [9], find the gene order *M* such that ${\sum}_{i=1}^{k}d\left({G}_{i},M\right)$ is minimized.

Although it plays a central role in gene order phylogeny, the median suffers from several liabilities. One is that it is hard to calculate in most metric spaces. Not only is it NP-hard [10], but exhaustive methods are costly for most instances, namely unless ${G}_{1}\dots \phantom{\rule{2.77695pt}{0ex}},{G}_{k}$ are all relatively similar to each other, which we will refer to generically as the *similar genomes condition*. Another problem is that heuristics tend to produce inaccurate results unless a suitable similar genomes condition holds [11]. Still another, is the tendency in some metric spaces to degenerate solutions [12] unless the same conditions prevails.

In this paper we add to this litany of difficulties by showing that as *k* genomes evolve over time, as modeled by any one of several biologically-motivated random walks, there is a phase change after *n*/4 steps, where *n* is the number of genes. With *u < n*/4 steps, the sum of the normalized distances ${\sum}_{k}d/n$ from each of the genomes to the starting point - the ancestor - converges to *ku/n* in probability, and this is the median value. When *u > c***n*/4 steps, for a constant *c** ≥ 1, the sum of the normalized distances to the median converges in probability to a value less than *ku/n*, and that the ancestor is no longer the median.

Our proof is inspired by a result of Berestycki and Durrett [13] in showing that the reversal distance between two signed permutations converges in probability to the actual number of steps, after rescaling, if and only if *u < n*/2. The technique is to construct a graph with genes as vertices and edges added between vertices according to how they are affected by transpositions. Properties of the number of components of random Erdös-Renyi graphs can then be invoked to prove the result.

## Definitions

We represent a unichromosomal genome by a signed permutation, where the sign indicates whether the gene is "read" from left to right (tail-to-head) or from right to left (head to tail) on the chromosome. Let ${S}_{n}^{\pm}$ be the signed symmetric group of order *n*, i.e. the space of all signed permutations of length *n*. A *reversal* operation applied to a signed permutation reverses the order, and changes the signs, of one or more adjacent terms in the permutation. A DCJ operation, which can apply not only to signed permutations but to more general genomes containing linear and circular chromosomes, *cuts* the genome in two places and rejoins pairs of the four "loose ends" in one of two possible new ways (one of which may be equivalent to a reversal). We define the reversal and DCJ distances, *d*_{
r
} and *dcj*, to be the minimum number of reversal and DCJ operations, respectively, needed to transform one genome to another.

*breakpoint graph BP*(Π, Π') of two genomes represented by Π and Π' contains vertices for the head and tail of each gene,

*black edges*edges defined by the adjoining heads or tails of two adjacent genes in the genome Π and

*grey edges*defined by two adjacent genes in the genome Π'. Let

*id*=

*I*, the identity permutation, and

*BP*(Π) =

*BP*(Π,

*id*). It is well-known that

*BP*(Π). We traverse a cycle

*c*∈

*cBP*(Π) in a counter-clockwise manner if we start at the left-most vertex of

*BP*(Π) (in the usual representation), travel along its unique adjacent black edge and end at the same vertex through its unique adjacent grey edge. Then we say a black edge in

*c*is positively oriented if we move along it from left to right in a counter-clockwise traversal. Otherwise we say it is negatively oriented. Similarly, for the grey edge (

*i*

_{ t }, (

*i*+ 1)

_{ h }) we say it is positively oriented if during a counter-clockwise traversal we move along it from

*i*

_{ t }to (

*i*+ 1)

_{ h }. Otherwise it is negatively oriented. We define the orientation function

*ξ*on the edges of

*BP*(Π) to be:

We say the black (grey) edges *e*, *e*' are parallel, denoted by *e* || *e*' if *ξ*(*e*) = *ξ*(*e*'). Otherwise we say they are crossing. This is just a reformulation of Hannenhalli and Pevzner's original concept of oriented cycles. An oriented cycle in this definition is a cycle including at least one positively and one negatively oriented black edge. The mechanism by which a reversal affects a genome can easily be seen using the *BP* graph. Let *ρ* be a reversal acting on two black edges *e*, *e*' in *BP*(Π). If they are in two different cycles we have a merger of the two to construct a new cycle. But if *e*, *e*' are in a same cycle, that cycle either splits, if *e* ∦ *e*', or does not split if *e* || *e*'.

## Limit Behavior of the Median Value

*d*

^{ n }be a metric on the space of all signed permutations length

*n*. For a set

*A*of these permutations, define

*m*^{
d,n
}(*A*) is called the median value of A under the metric *d*^{
n
}. A signed permutation which makes ${g}_{A}^{d,n}$ minimum is called a median solution of *A*. Denote by *d*_{
r
} and *dcj* the reversal and DCJ distances on ${S}_{n}^{\pm}$.

*X*

_{0}=

*id*, the identity permutation, and let ${X}_{t}^{n}$ be a stochastic process on ${S}_{n}^{\pm}$, where at random Poisson times

*τ*

_{ κ }, with rate 1, we choose two elements of ${X}_{{\tau}_{\kappa}}^{n},$ namely

*i*,

*j*and let

*ρ*(

*i*,

*j*) operate on ${X}_{{\tau}_{\kappa}}^{n},$ that is

*ρ*(

*i, j*) is the reversal acting on

*i*and

*j*. We call ${X}_{t}^{n}$ a reversal random walk (r.w.) on. ${S}_{n}^{\pm}.$ Suppose ${X}_{t}^{1,n},\phantom{\rule{2.77695pt}{0ex}}\dots ,\phantom{\rule{2.77695pt}{0ex}}{X}_{t}^{k,n}$ be

*k*independent reversal r.w. all starting at the identity element,

*id*. Define

*m*

^{ d,n }(

*A*

_{ t }), remains a good estimator for the total divergence time,

*kt*, as well as to the total distance of points in

*A*

_{ t }to

*id*, namely ${g}_{{A}_{t}}^{d,n}\left(id\right)$. To answer this question we use the fact that the speed of escape of the r.w. up to some particular time, is the same from any point of the space and is close to 1, the maximum value. Berestycki and Durrett studied speed of transposition and reversal random walks with the related edit distances while in the latter they used "approximate reversal distance" instead of reversal itself, ignoring the effect of hurdles and fortresses. This turns out to be the same as DCJ distance on single chromosomes. We have

where *h*(*π*) and $\stackrel{~}{f}\left(\pi \right)$ are the number of hurdles and fortresses, respectively.

Although Berestycki and Durrett only proved their theorem for the random transposition r.w. on *S*_{
n
}, they suggested that same method should carry over to reversal r.w. The following proposition is proved in [13] for *approximate* reversal distance (i.e., DCJ distance).

In this result and in the ensuing discussion *a*_{
n
} is an arbitrary sequence such that *a*_{
n
} → ∞ as *n* → 0. When it is unambiguous we drop *n* from ${A}_{t}^{\left(n\right)}$ and ${X}_{t}^{n}$.

**Propostition 1**

*[Berestycki-Durrett] Let c be fixed and let Xt be a reversal r.w. on*${S}_{n}^{\pm}$

*starting at id. Then*

*where*

*and* $\frac{w\left(n\right)}{{a}_{n}\sqrt{n}}\to 0$*in probability*.

**Remark 1**

*The function*1 -

*f is linear for c <*1,

*f*(

*c*) = 1 -

*c*/2,

*and sublinear for c >*1, 1 -

*f*(

*c*)

*< c*/2

*This means that for c*≤ 1

*and r.w. travels on an approximate geodesic (or parsimonious path) asymptotically almost surely. f is the function counting the number of tree components of an Erdös-Renyi random graph with n vertices for which the probability of having each edge is* $\frac{c}{n}$, *denoted by G*(*c, n*). *See Theorem 12 in* [14], *Chapter V*.

We extend the above theorem for the bonafide reversal distance. To do so we need to estimate the number of hurdles of ${X}_{\frac{cn}{2}}.$ Recall that an oriented cycle in a breakpoint graph is a cycle including an orientation edge, that is a grey edge with two black adjacency edges *e*, *e*', where a reversal involving *e* and *e*' splits the cycle [15]. As we discussed this is equivalent to saying *e* ∦ *e*'. It is not difficult to show

**Lemma 1** *Let C* ∈ *cBP*(π), *then C is oriented if and only if there exists exactly two equivalence classes of black edges, that is there exist at least two black edges with different signs*.

Then

**Theorem 1**

*Let c >*0

*be fixed and let X*

_{ t }

*be a reversal r.w.starting at id. Define h*

_{ t }:=

*h*(

*X*

_{ t })

*to be the number of hurdles in BP*(

*X*

_{ t }).

*Then*

**Proof**. Cycles of the BP that have never been involved in a fragmentation event must be oriented, since the two rejoined black edges resulting from an inversion-induced merger of cycles cannot be parallel.

Therefore we need only to count the number of edges that have been involved in a fragmentation event. To do so we apply the method of counting cycles in [13], Theorem 3. Hurdles occur only in those cycles with length more than one that have been involved in a fragmentation up to time $\frac{cn}{2}.$ We call such cycles fragmented cycles. The number of fragmented cycles with length more than $\sqrt{n}$ is always less than $\sqrt{n}$. But to count all fragmented cycles in ${X}_{\frac{cn}{2}}$ with size less than $\sqrt{n}$ we need to find an upper bound for the rate of a fragmentation up to time $\frac{cn}{2}$. Since a fragmentation occurs when two black edges in one cycle are chosen, to fragment a cycle in BP, for any chosen black edge *e* we only can pick another black edge *e*' in the same cycle whose graph distance in the breakpoint graph is less than $2\sqrt{n}$. (The coefficient 2 arises from the fact that the cycles are alternating in BP.)

Thus the rate of fragmentation at an arbitrary time *t* is not more than $\frac{n}{n}\cdot \frac{2\left(\sqrt{n}\right)}{n}\phantom{\rule{2.77695pt}{0ex}}=\frac{2}{\sqrt{n}}$. Integrating up to time *t*, this gives us the expected number of fragmented cycles at time *t* is $\frac{2}{\sqrt{n}}\xb7t$. For $t=\frac{cn}{2}$ this expectation is $c\sqrt{n}.$ Now, dividing by ${a}_{n}\sqrt{n}$, the result follows from Chebyshev's inequality and the fact that hurdles only occurs in fragmented cycles. ■

**Theorem 2**

*let c >*0

*be fixed and let X*

_{ t }

*be a reversal r.w. on*${S}_{n}^{\pm}$

*starting at id and let*${d}_{r}:={d}_{r}^{\left(n\right)}$

*denote the reversal distance on*${S}_{n}^{\pm}$

*. Then*

*where f is the same function as in the statement of Proposition 1 and w*' (*n*) *is a function with* $\frac{{w}^{\prime}\left(n\right)}{{a}_{n}\sqrt{n}}\to 0$*in probability*.

**Proof**. Since

*d*

_{ r }(Π) =

*dcj*(Π) +

*h*(Π) +

*f*˜(Π) by the proposition we have

*d*

_{ r }(

*X*

_{ cn/2 }) = (1

*− f*(

*c*))

*n*+

*w*(

*n*) +

*h*

_{ cn/2 }+

*f*˜(

*X*

_{ cn/2 }). But

in probability, by the convergence of $\frac{w\left(n\right)}{{a}_{n}\sqrt{n}}$ and $\frac{{h}_{cn/2}}{{a}_{n}\sqrt{n}}$ in Proposition 1 and Theorem 2 and $\stackrel{~}{f}\left({X}_{cn/2}\right)\phantom{\rule{0.3em}{0ex}}\le 1.\u25a0$

**Theorem 3**
*Let*
${X}_{t}^{1,n},\phantom{\rule{2.77695pt}{0ex}}\dots ,\phantom{\rule{2.77695pt}{0ex}}{X}_{t}^{k,n}$
*be k independent reversal r.w in*
${S}_{n}^{\pm}$
*starting at id. Suppose either*

*a) d* := *dcj dcj distance*

*or*

*b)* $d:={d}_{r}^{\left(n\right)}$ *reversal distance*.

Then for c < $\frac{1}{4}$*we have* $\frac{{\epsilon}_{{}^{{c}_{n}}}^{d,n}}{{a}_{n}\sqrt{n}}\to 0$*in probability*.

**Proof**. We prove the theorem only for

*d*

_{ r }. The proof of the DCJ case is similar. For all

*i*,

*j*∈ {1, ...,

*k*} and for a median solution

*x*of ${A}_{t}^{\left(n\right)}$

*i, j i*≠

*j*

*k*'. Also $|{g}_{{A}_{cn}^{\left(n\right)}}\left(id\right)-kcn|\phantom{\rule{2.77695pt}{0ex}}\le kw\left(n\right)$. Therefore, there exists a constant

*k**such that

This proves the theorem. ■

**Remark 2** *The statement of the theorem suggests ignoring the error of order* $o\left({a}_{n}\sqrt{n}\right)$*for a*_{
n
} → ∞. *id remains as the median of leaves of k independent stochastic processes* ${X}_{t}^{1,n},\phantom{\rule{2.77695pt}{0ex}}\dots ,\phantom{\rule{2.77695pt}{0ex}}{X}_{t}^{k,n}$*up to time* $\frac{n}{4}$*asymptotically almost surely*.

**Theorem 4**

*Let*$c\phantom{\rule{2.77695pt}{0ex}}\le \frac{1}{4}$

*be fixed. Suppose d is either DCJ or reversal distance. Then by the hypothesis of Theorem 3*

**Proof**. This follows directly from the fact that

in probability. ■

Now, it is natural to ask whether the statement of Theorem 4 also holds for some time after $\frac{n}{4}$. In other words, is the median value *kcn* a fair estimator for the total time of divergence? We conjecture not, that the property is lost after time $\frac{n}{4}$, but for now can only prove a weaker upper bound for this time.

**Theorem 5**

*Let*$c>\frac{1}{2}$

*be fixed. Suppose d is either DCJ or reversal distance. Then by the same hypothesis as in Theorem 3*

*where*

*is strictly positive for*
$c>\frac{1}{2}$

**Remark 3** *This theorem shows after time* $\frac{n}{2}$*the error is of order n and so the median value is not a good estimate of k times the divergence time*.

**Proof.**

where $\frac{w\left(n\right)}{{a}_{n}\sqrt{n}}\to 0$ in probability. Dividing by *n*, the result follows. ■

In fact, since *f* (*c*), *c >* 0 is decreasing and for *c <* 1, $f\left(c\right)=1-\frac{c}{2}$, it is easy to see that in the case *k* = 3, for *c >* 0.75, ${\epsilon}_{cn}^{d,n}$ is of order ${\beta}_{c}^{d}n$ for some ${\beta}_{c}^{d}\ge 0$.

**Theorem 6**

*Let k*= 3

*and d be either dcj or dr. Consider the same hypothesis in Theorem 3. Assume c* be solution of*

*Then for all c > c**

*there exists*${\beta}_{c}^{d}$

*such that*

**Proof.**

*i*= 2, 3 is the same as $d\left(id,\phantom{\rule{2.77695pt}{0ex}}{X}_{\frac{cn}{2}}^{1,\phantom{\rule{2.77695pt}{0ex}}n}\right)$. This is true since the Cayley graph of ${S}_{n}^{\pm}$ w.r.t. reversals is symmetric and regular and so $P\left({X}_{0}=id,\phantom{\rule{2.77695pt}{0ex}}{X}_{\frac{cn}{4}}=\mathrm{\Pi}\right)=P\left({X}_{0}=\mathrm{\Pi},\phantom{\rule{2.77695pt}{0ex}}{X}_{\frac{cn}{4}}=id\right)$. But therefore by symmetry of the Cayley graph we can just consider $d\left(id,\phantom{\rule{2.77695pt}{0ex}}{X}_{\frac{cn}{2}}^{1,\phantom{\rule{2.77695pt}{0ex}}n}\right)$. Hence,

*x >*0 be so that

*x*=

*c*>

*c*

^{∗}. Since

*f*(

*x*)

*>*0 for all

*x >*0

*c*≥

*c**

This proves the statement. ■

*n*goes to

*∞*, showing that the median is only useful in a highly restricted region of the space.. The following theorem is entailed by a theorem in [16]. Let

*c*

_{ n }=

*c*

_{ n }(Π) be the number of cycles in the BP graph of a random $\mathrm{\Pi}\in {S}_{n}^{\pm}$. Let

*d*

_{ n }be a distance (metric) on ${S}_{n}^{\pm}$. Define

to be the ball of radius *cn* in ${S}_{n}^{\pm}$.

**Theorem 7**

*Let*0

*< c <*1

*be fixed. Then*

**Proof.**

- a)$For\phantom{\rule{2.77695pt}{0ex}}all\phantom{\rule{2.77695pt}{0ex}}\mathrm{\Pi}\in {B}_{cn}^{dcj}\phantom{\rule{2.77695pt}{0ex}},\phantom{\rule{2.77695pt}{0ex}}\left|cBP\left(\mathrm{\Pi}\right)\right|\phantom{\rule{2.77695pt}{0ex}}\ge \left(1-c\right)n.$(42)

_{ n }does not converge to 0. Therefore there exists a subsequence ${\left\{{n}_{i}\right\}}_{i\in \mathbb{N}}$ such that ${\gamma}_{{n}_{i}}\ge \epsilon $ for a constant

*ε >*0. This implies

- b)For the second part it suffices to observe that for all $\mathrm{\Pi}\in {S}_{n}^{\pm}$ we have${d}_{r}\left(\mathrm{\Pi}\right)\ge dcj\left(\mathrm{\Pi}\right).$(46)

■

## Conclusion

We have shown that the median value for DCJ and for reversal distance for a reversal r.w.has good limiting properties if the number of steps remains below *cn*/4, for any *c <* 1, but for some value *c* > 1, more than this number of steps destroys these limiting properties. The critical value may indeed be *c* = 1, but for now we can only show that for *c* > 3 (and *c* > 2) the median value is no longer a good estimator of the distance between the *id* and the current position of the r.w. (and *k* times the divergence time, respectively).

Note that a simulation strategy to estimate *c* is not available because of the hardness of calculating the median. As *n* increases even to moderate values all exact methods require prohibitive computing time.

These results imply that the steinerization strategy for the small phylogeny problem may lead to poor estimates of the interior nodes of a phylogeny unless the taxon sampling is sufficient to assure that a "similar genomes condition" holds for every k-tuple of genomes used in the course of of the iterative optimization search. This can be monitored prior to each step in the iterative optimization of the phylogeny through successive application of the median method.

## Declarations

### Acknowledgements

Research supported in part by grants from the Natural Sciences and Engineering Research Council of Canada. DS holds the Canada Research Chair in Mathematical Genomics. Thanks to Armin Jamshidpey and Leili Rafiee Sevyeri for help in preparation of the manuscript.

**Declarations**

Publication of this article was supported by the Canada Research Chair in Mathematical Genomics.

This article has been published as part of *BMC Bioinformatics* Volume 14 Supplement 15, 2013: Proceedings from the Eleventh Annual Research in Computational Molecular Biology (RECOMB) Satellite Workshop on Comparative Genomics. The full contents of the supplement are available online at http://www.biomedcentral.com/bmcbioinformatics/supplements/14/S15.

## Authors’ Affiliations

## References

- Sankoff D, Cedergren RJ, Lapalme G: Frequency of insertion/deletion, transversion and transition in the evolution of 5S ribosomal RNA. Journal of Molecular Evolution. 1976, 7: 133-149. 10.1007/BF01732471.View ArticlePubMedGoogle Scholar
- Blanchette M, Bourque G, Sankoff D: Breakpoint phylogenies. Genome Informatics. Edited by: S. Miyano & T. Takagi. 1997, Tokyo: Universal Academy Press, 25-34.Google Scholar
- Sankoff D, Abel Y, Hein J: A Tree - A Window - A Hill; Generalization of nearest-neighbour inter-change in phylogenetic optimisation. Journal of Classification. 1994, 11: 209-232. 10.1007/BF01195680.View ArticleGoogle Scholar
- Huson D, Nettles S, Warnow T: Disk-covering, a fast-converging method for phylogenetic tree reconstruction. Journal of Computational Biology. 1999, 6: 369-386. 10.1089/106652799318337.View ArticlePubMedGoogle Scholar
- Tang J, Moret B: Scaling up accurate phylogenetic reconstruction from gene-order data. Bioinformatics. 2003, 19: i305-i312. 10.1093/bioinformatics/btg1042.View ArticlePubMedGoogle Scholar
- Sankoff D, Blanchette M: The median problem for breakpoints in comparative genomics. Proceedings of Computing and Combinatorics (COCOON). Edited by: T. Jiang and D.T. Lee. 1997, 1276: 251-263. 10.1007/BFb0045092. Lecture Notes in Computer ScienceView ArticleGoogle Scholar
- Sankoff D, Sundaram G, Kececioglu J: Steiner points in the space of genome rearrangements. International Journal of the Foundations of Computer Science. 1996, 7: 1-9. 10.1142/S0129054196000026.View ArticleGoogle Scholar
- Bourque G, Pevzner PA: Genome-scale evolution: Reconstructing gene orders in the ancestral species. Genome Research. 2002, 12: 26-36.PubMed CentralPubMedGoogle Scholar
- Zhang M, Arndt W, Tang J: An exact solver for the DCJ median problem. Pacific Symposium on Biocomputing. 2009, 138-149.Google Scholar
- Tannier E, Zheng C, Sankoff D: Multichromosomal median and halving problems under different genomic distances. BMC Bioinformatics. 2009, 10: 120-10.1186/1471-2105-10-120.PubMed CentralView ArticlePubMedGoogle Scholar
- Zheng C, Sankoff D: On the Pathgroups approach to rapid small phylogeny. BMC Bioinformatics. 2011, 12: S4-PubMed CentralPubMedGoogle Scholar
- Haghighi M, Sankoff D: Medians seek the corners, and other conjectures. BMC Bioinformatics. 2012, 13 (S19): S5-PubMed CentralPubMedGoogle Scholar
- Berestycki N, Durrett R: A phase transition in the random transposition random walk. Probability Theory and Related Fields. 2006, 136: 203-233. 10.1007/s00440-005-0479-7.View ArticleGoogle Scholar
- Bollobás B: Random Graphs. 2001, Cambridge University Press, 2View ArticleGoogle Scholar
- Hannenhalli S, Pevzner PA: Transforming cabbage into turnip: Polynomial algorithm for sorting signed permutations by reversals. Journal of the ACM. 1999, 46: 1-27. 10.1145/300515.300516.View ArticleGoogle Scholar
- Székely LA, Yang Y: On the expectation and variance of the reversal distance. Acta Univ. Sapientiae, Mathematica. 2009, 1: 5-20.Google Scholar

## Copyright

This article is published under license to BioMed Central Ltd. This is an open access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.