Restricted DCJ Model: Rearrangement Problems with Chromosome Reincorporation

Abstract

We study three classical problems of genome rearrangement—sorting, halving, and the median problem—in a restricted double cut and join (DCJ) model. In the DCJ model, introduced by Yancopoulos et al., we can represent rearrangement events that happen in multichromosomal genomes, such as inversions, translocations, fusions, and fissions. Two DCJ operations can mimic transpositions or block interchanges by first extracting an appropriate segment of a chromosome, creating a temporary circular chromosome, and then reinserting it in its proper place. In the restricted model, we are concerned with multichromosomal linear genomes and we require that each circular excision is immediately followed by its reincorporation. Existing linear-time DCJ sorting and halving algorithms ignore this reincorporation constraint. In this article, we propose a new algorithm for the restricted sorting problem running in O(n log n) time, thus improving on the known quadratic time algorithm. We solve the restricted halving problem and give an algorithm that computes a multilinear halved genome in linear time. Finally, we show that the restricted median problem is NP-hard as conjectured.

1. Introduction

During evolution, genomes undergo large-scale mutations, for example, a segment of DNA can get reversed or moved to another position. In genome rearrangement problems, we define a genome model and a fixed set of allowed operations, and then we try to minimize the number of rearrangement operations in the evolutionary history. In this article, we study the classical problems of genome rearrangement—the sorting, halving, and the median problems—in a restricted double cut and join model.

Consider two extant species Π and Γ with a common ancestor Λ (Fig. 1). Starting with the common ancestral gene order, after speciation the two species evolve independently and undergo mutations. Nowadays, we can determine the gene order of the extant species Π and Γ, and we can ask how many mutations are needed to transform one genome into the other. This is called the (genomic) distance between Π and Γ. Note that the rearrangement operations are symmetric, so the distance from Π to Γ is the same as the distance from some intermediate genome Λ′ to Π plus the distance from Λ′ to Γ. The genomic distance between Π and Γ is thus a lower bound and a good estimate for the number of mutations throughout the evolution of Π and Γ from their common ancestor Λ. Computing the genomic distance and one particular sequence of rearrangement operations (a rearrangement scenario) is the sorting problem.

FIG. 1.

Evolution of gene orders Π and Γ from the common ancestor Λ.

Note that the gene orders Π and Γ do not provide enough information to reconstruct the ancestral genome Λ. Therefore, let us consider gene order Θ of an outgroup species. The ancestor of Π and Γ is then reconstructed according to the parsimony principle by finding a genome M that minimizes the total distance from Π, Γ, and Θ. Such a genome is called median, and the corresponding problem is called the median problem.

In the halving problem, we imagine a genome that underwent a whole genome duplication and then evolved by large-scale rearrangements (Fig. 2). Thus, even though immediately after the duplication, the genome had two perfect copies of each chromosome; nowadays, we observe the copies scattered over the genome. The goal of the halving problem is, given a present genome with two (paralogous) copies of each marker, to reconstruct the ancestral genome immediately after the duplication. In other words, the goal of genome halving is to find a perfectly duplicated genome with the smallest genomic distance from the given genome.

FIG. 2.

The ancestral genome undergoes a whole genome duplication and subsequently evolves by rearrangements. The goal of genome halving is to reconstruct the genome immediately after the whole genome duplication, given the gene order of the extant species.

In this article, we study these problems in the double cut and join (DCJ) model, introduced by Yancopoulos et al. (2005). A DCJ operation models most of the large-scale mutation events—such as reversals, translocations, fusions, and fissions—in a unified way. Furthermore, transpositions and block interchanges can be simulated by two operations: an appropriate segment of a chromosome is extracted, creating a temporary circular chromosome, which is then reinserted at the proper place in the next step.

The sorting algorithm given by Yancopoulos et al. (2005) running in quadratic time guarantees that each new circular chromosome is immediately reincorporated, thus mimicking transpositions and block interchanges.

Bergeron et al. (2006) restated the model and gave a simple linear-time algorithm for DCJ sorting ignoring the reincorporation constraint. However, the algorithm finds a sequence of DCJ operations without any explicit mention of the underlying operations (e.g., reversals, translocations, block interchanges), and many circular chromosomes may coexist in intermediate stages of the sorting process. Such sorting sequences are not biologically plausible, e.g. in eukaryotic organisms, which typically have only linear chromosomes (Fig. 3).

FIG. 3.

Two optimal DCJ sorting scenarios transforming the first genome into the last one in 8 operations. (a) In the unrestricted model, we can get the scenario on the left. Even though the number of operations is minimal, the scenario is not very biologically plausible. Note that we started and ended with linear chromosomes, but the intermediate genomes contain several circular chromosomes. (b) In contrast, the scenario on the right is restricted—each circular excision is immediately followed by its reincorporation. Thus, we can, for instance, treat the first circular excision plus reintegration as a transposition of the block (h, i, j, k) after marker b.

In this work, we revisit the original study of Yancopoulos et al. (2005), restricting the DCJ model to genomes with multiple linear chromosomes. We require that each excision of a circular chromosome is always immediately followed by its reincorporation.

We borrow techniques from other studies on sorting by reversals and block interchanges (Christie, 1996; Feng and Zhu, 2007; Swenson et al., 2010) and propose a new algorithm that sorts multichromosomal linear genomes in the DCJ model, reincorporates circular chromosomes, and runs in O(n log n) time. This improves the original result by Yancopoulos et al. (2005).

Furthermore, we present a new result on the halving problem. If no restriction on the linearity of chromosomes is imposed and no guarantee concerning circular reintegration is required, we can use linear-time algorithms proposed by Warren and Sankoff (2009) and Mixtacki (2008) to reconstruct the ancestor. However, given a multilinear genome, these algorithms may predict some circular chromosomes in the ancestral genome. In the worst case, these algorithms may even produce Ω(n) circular chromosomes given a single linear chromosome of length n. Again, this is not biologically plausible, when organisms with linear genomes are considered.

The restricted halving problem has not been studied previously and is stated as open in Tannier et al. (2009). In this article, we propose a new algorithm that reconstructs a multichromosomal linear perfectly duplicated ancestor in linear time. One particular halving scenario (reincorporating circular chromosomes) can be found in O(n log n) time.

Finally, we show that the median problem is NP-hard in the restricted DCJ model, as conjectured by Tannier et al. (2009).

The article is organized as follows: in Section 2, we review the previous results. We introduce the DCJ model and the restricted versions of sorting, halving, and median problems. New results on these problems are presented in Sections 3, 4, and 5, respectively. We conclude in Section 6.

2. Preliminaries

2.1. The double cut and join model

Genome model

In the DCJ model, genomes Π and Γ consist of the same set of markers (genes, synteny blocks). Every marker g has two ends, called extremities, which we denote g⁻ and g⁺.

Each extremity p is either adjacent to some other extremity q (two consecutive markers on a chromosome) or to a telomere—the end of a linear chromosome. In the first case, we say the extremities form an adjacency pq; in the second case, we have a telomeric adjacency p○. Thus, a genome is a set of adjacencies such that every extremity is in exactly one (possibly telomeric) adjacency.

For genome Π, we can draw a genome graph G_Π where vertices are extremities and edges connect either adjacent extremities or two extremities of the same marker. Every vertex in this graph has degree 1 or 2, so the components of G_Π are paths and cycles. These components represent chromosomes—linear and circular, respectively.

Multilinear genomes

In this article, we will be mainly interested in multilinear genomes, i.e., genomes containing only linear chromosomes. These are more comfortably written as signed permutations: Choose a direction of a linear chromosome, and list the markers from left to right; write \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $${\overrightarrow g}$$\end{document} , if extremity g⁻ precedes g⁺ and \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $${\overleftarrow g}$$\end{document} otherwise. Thus, chromosome \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$({\overrightarrow 1} , {\overrightarrow 3} , {\overleftarrow 2})$$\end{document} (which is the same as \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$( { \overrightarrow 2} , { \overleftarrow 3} , { \overleftarrow 1} )$$\end{document} ) corresponds to the set of adjacencies \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$\{ \circ 1^ - , \ 1^ + 3^ - , \ 3^ + 2^ + , \ 2^ - \circ \} $$\end{document} . We will write −g for the reversed marker g, i.e., \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$- { \overleftarrow g} = { \overrightarrow g}$$\end{document} and \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$- { \overrightarrow g} = { \overleftarrow g}$$\end{document} .

DCJ operation

A double cut and join operation acting on adjacencies pq and rs replaces them by either adjacencies pr, qs, or ps, qr (the adjacencies pq and rs can be telomeric or even an empty chromosome ○○). We say that the operation cuts pq and rs and joins either pr, qs, or ps, qr. By these operations, we can mimic every common rearrangement operation in genomes: To invert a segment, we cut at its ends and join reversed. By cutting and joining adjacencies on different linear chromosomes, we get a translocation (Fig. 4a). By cutting two telomeric adjacencies p○ and q○ and joining pq, we can fuse two chromosomes into one or create a circular chromosome from a linear one (as a byproduct, we get an empty chromosome ○○; see Fig. 4b for the reversed process—fission). By two DCJ operations, we can mimic transpositions and block interchanges: First, we cut out an appropriate segment and by joining its ends create a temporary circular chromosome. In the next step, we reincorporate it into the original chromosome (Fig. 5).

FIG. 4.

Examples of operations in the DCJ model (a, b).

FIG. 5.

To interchange blocks a and c (left) in the DCJ model, we cut before a and after b and create a temporary circular chromosome (center). The next operation cuts between a and b and after c and reincorporates the blocks in the correct order (right).

2.2. Sorting

DCJ distance and scenarios

A sequence of k DCJ operations transforming a given genome Π into Γ is called a DCJ scenario of length k. A scenario of minimum length is called optimal and its length is the DCJ distance between Π and Γ, denoted d(Π, Γ). A sequence of k DCJ operations transforming Π into Π′ is optimal (with respect to Γ), if d(Π, Γ) = d(Π′, Γ) + k.

The distance and a sorting scenario can be calculated using an adjacency graph AG(Π, Γ) (Bergeron et al., 2006). This is a bipartite multigraph where vertices are adjacencies of Π and Γ; an adjacency in Π is connected with an adjacency in Γ by one edge for each extremity that they share. Since every adjacency is connected with one (telomeric) or two other adjacencies, this graph consists of paths and cycles only. If Π and Γ share a common adjacency, this corresponds to a cycle of length 2 or path of length 1 (common telomere) in the adjacency graph. Note that, when Π and Γ are equal, their adjacency graph consists of 2-cycles and 1-paths only.

The following theorem gives the DCJ distance between two genomes:

Theorem 1 (Bergeron et al., 2006)

Given two genomes Π and Γ on n markers, let c be the number of cycles and p_o the number of paths of odd length in the adjacency graph AG(Π, Γ). Then the distance between Π and Γ is \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} \begin{align*}d ( \Pi , \Gamma ) = n - ( c + p_o / 2 ).\end{align*}\end{document}

Given multilinear genomes, we call a sorting DCJ scenario restricted if every DCJ operation that creates a circular chromosome is immediately followed by another operation that reintegrates it into the original chromosome (Fig. 3b). Such scenarios can be viewed as sequences of reversals, translocations, fusions, fissions (with weight 1), and block interchanges, which have weight 2, i.e., count as two operations. In the restricted sorting problem, we are searching for a restricted sorting scenario of minimal length.

2.3. Median problem

Given three genomes Π₁,Π₂,Π₃, the median problem is to find genome Π_M, called median, which minimizes the median score \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} \begin{align*}score ( \Pi_M ) = d ( \Pi_M , \Pi_1 ) + d ( \Pi_M , \Pi_2 ) + d (\Pi_M , \Pi_3 ).\end{align*}\end{document}

In the restricted median problem, we are given multilinear genomes Π₁, Π₂, Π₃, and we require the median Π_M to be multilinear too.

2.4. Halving problem

In the halving problem, we imagine a genome that underwent a whole genome duplication and then evolved into genome Γ. We are given Γ, and our goal is to reconstruct the genome right after the whole genome duplication. More formally: In a duplicated genome, every marker g has two copies − \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$g_1^ - g_1^ +$$\end{document} and \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$g_2^ - g_2^ +$$\end{document} . If p is an extremity, we will denote by \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $${ \bar p}$$\end{document} the other copy—the paralogous extremity. Similarly, if x = pq is an adjacency (possibly telomeric), then \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $${ \bar x} = { \bar p}{ \bar q}$$\end{document} is the paralogous adjacency, and if C is a chromosome (set of adjacencies), then \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $${ \overline C}$$\end{document} is the set of paralogous adjacencies.

We say that genome Θ is a doubled (or a perfectly duplicated) genome, if for each adjacency pq in Θ, adjacency \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $${ \bar p}{ \bar q}$$\end{document} is also in Θ and \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$p \ne { \bar q}$$\end{document} . This is the same as saying that, if we ignore the subscripts (1's and 2's), every linear chromosome has an identical copy and every circular chromosome has either an identical copy, or is itself composed of two successive identical copies.

The genome halving problem can be stated as follows: Given a duplicated genome Γ, find a doubled genome Θ such that d(Θ, Γ) is minimal.

The halving distance and scenario can be calculated using a concept similar to an adjacency graph—a natural graph NG(Γ) introduced by El-Mabrouk and Sankoff (2003). Vertices of this multigraph are adjacencies of Γ, and two adjacencies are connected by an edge, if they share a paralogous extremity. The natural graph consists of paths and cycles only, and Θ is doubled if and only if NG(Θ) consists of 2-cycles and 1-paths only.

The following theorem gives the DCJ halving distance:

Theorem 2 (Mixtacki, 2008)

Let Γ be a duplicated genome with 2n markers. The minimal distance between Γ and any doubled genome Θ is \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} \begin{align*}d ( \Gamma , \Theta ) = n - ( c_e + \lfloor\, p_o / 2 \rfloor ) ,\end{align*}\end{document}

where c_e is the number of even cycles and p_o the number of odd paths in the natural graph NG(Γ).

In the restricted halving problem, we are given a multilinear duplicated genome Γ, and the goal is to find a multilinear doubled genome Θ with a restricted halving scenario of minimal length.

2.5. Data structure for handling permutations

Our sorting algorithm uses a data structure for handling permutations by Kaplan and Verbin (2005). It can be traced back to Chrobak et al. (1990), where it was used to improve heuristics for the traveling salesman problem. It supports the following three operations in logarithmic time: find the i^th marker in a linear chromosome, return the position of marker g, and perform a reversal operation.

Linear chromosomes can be represented by a balanced tree supporting operations split and merge (e.g., red-black tree or splay tree). The order is the same as the left-to-right order of markers on the chromosome. In each node of the tree, we store one marker, its orientation, number of descendants, and a reverse flag. A reverse flag being “on” signifies that the whole subtree is reversed. The reverse flag of node v can be cleared (“pushed down”) by changing v's orientation, swapping its children and flipping their reverse flags.

Reversing a segment from i to j can be implemented as follows:

1. Find the i^th and j^th marker (using the information about sizes of subtrees and reverse flags).

2. Split the tree into three parts: T₁ with markers before i, T₃ with markers after j, and T₂ with the segment from i to j.

3. Flip the reverse flag in the root of T₂, and

4. Merge T₁, T₂ and T₃.

We store a lookup table with a pointer to the corresponding node of a tree for every marker. In this way, we can find the position of any marker in logarithmic time.

To support multilinear genomes, we simply concatenate the chromosomes with a delimiter between each pair, and in each node we store the number of delimiters in its subtree. This way, given a marker g, we can tell on which chromosome it is by counting the number of delimiters before g. To support different rearrangement operations, we can express them as a sequence of reversals. For example, block interchange can be mimicked by four reversals; if we add sufficiently many delimiters at the end of the sequence (representing empty chromosomes), we can also mimic fusions and fissions.

3. Restricted DCJ Sorting

Bergeron et al. (2006) gave a linear-time algorithm for DCJ sorting disregarding the constraint of reincorporating circular chromosomes immediately. The solution can be easily adapted to a quadratic-time algorithm for the restricted version: after each step, check whether a circular chromosome was created and if so, find the appropriate DCJ operation acting on adjacencies in the circular and the original linear chromosome that reintegrates the circular chromosome. It is not obvious how to do this efficiently (say in polylogarithmic time).

Yancopoulos et al. (2005) proposed to transform Π into Γ by restricted sorting in four stages: (0) Add caps to the ends of linear chromosomes. (1) By translocations, fusions and fissions transform Π into Π′ such that chromosomes in Π′ and Γ have the same marker contents. (2) Perform oriented reversals to get Π″ with all markers in the same direction as in Γ. (3) Finally, use block interchanges to transform Π″ into Γ.

Stages 2 and 3 can be implemented in O(n log n) time using the data structure described in Section 3 (Swenson et al., 2010; Feng and Zhu, 2007). Thus, a unichromosomal restricted DCJ sorting can be solved in O(n log n) time. However, it is not obvious how to implement stage 1 efficiently.

Our algorithm is based on the following observation:

Observation 1

Let g, h be two markers that are adjacent in Γ, but not in Π. If g and h are on different chromosomes in Π, there is a translocation that puts them together. If g and h are on the same chromosome and have a different orientation, there is a reversal that puts them together. These operations are optimal in the DCJ model. Transposition and block interchange take two DCJ operations. These operations are optimal if they create two new non-telomeric common adjacencies and destroy none.

This is simply because, even more generally, k operations, that create k new non-telomeric adjacencies and destroy none, create k new cycles in the adjacency graph, and thus decrease the distance by k.

Theorem 3

A restricted optimal DCJ scenario transforming multilinear genome Π into Γ can be found in O(n log n) time.

Proof

The ends of linear chromosomes, telomeres, produce some difficulties and nasty special cases. Capping is an elegant technique to deal with them: we adjoin new markers (caps) to the ends so that we do not change the distance and we do not have to worry about telomeres any more.

We find all the paths in the adjacency graph AG (Π, Γ). Paths of odd length have one end in Π and one in Γ—simply adjoin a new marker (properly oriented) to the two telomeres. This increases the number of markers by one, but instead of an odd path, we have a cycle and a 1-path, so the distance does not change. For paths starting and ending in Π, add two new markers to the ends of Π and add a new chromosome consisting of just these two markers (properly oriented) to Γ. The case with a path starting and ending in Γ is symmetric. The number of markers increases by 2, but instead of an even path, we have a cycle and two odd paths, so the distance does not change. Capping of all chromosomes can be done in linear time.

Without loss of generality, we may assume that the markers in chromosomes of Γ are consecutive numbers \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$(\overrightarrow {k_0} , \ldots , \overrightarrow {k_1 - 1} ) , (\overrightarrow {k_1} , \ldots , \overrightarrow {k_2 - 1} ) , \ldots , ( \overrightarrow {k_{s - 1}} , \ldots , \overrightarrow {k_s - 1} )$$\end{document} where \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$0 = k_0 < k_1 < k_2 < \cdots < k_s = n$$\end{document} (otherwise renumber the markers).

We will be transforming Π into Γ gradually “from left to right”: once we have transformed the beginning of a chromosome in Π to \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$\overrightarrow {k_i} , \overrightarrow {k_i + 1} , \ldots , \overrightarrow {j}$$\end{document} , we extend it by moving j + 1 next to \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $${ \overrightarrow j}$$\end{document} .

There are several cases we need to consider:
1. If \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$\overrightarrow {j + 1}$$\end{document} is already next to \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $${ \overrightarrow j}$$\end{document} , we are done.

2. If j + 1 is on a different chromosome than \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $${ \overrightarrow j}$$\end{document} , we can always use a translocation. In the rest of the proof, we assume that j + 1 is on the same chromosome, to the right of \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $${ \overrightarrow j}$$\end{document} .

3. If \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $${ \overrightarrow j}$$\end{document} and \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$\overleftarrow {j + 1}$$\end{document} have different orientation, we can use a reversal.

Otherwise, following Christie (1996), find the marker m with the highest number between \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $${ \overrightarrow j}$$\end{document} and \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$\overrightarrow {j + 1}$$\end{document} and find m + 1.
4. If m + 1 is on a different chromosome, we can use a translocation to move it next to m; this operation also moves \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$\overrightarrow {j + 1}$$\end{document} to another chromosome, so we can use another translocation to move it next to \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $${ \overrightarrow j}$$\end{document} .

Otherwise the situation is \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $${ \overrightarrow j} , \ldots , m , \ldots , \overrightarrow {j + 1} , \ldots , m + 1$$\end{document} (since m is the highest number between \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $${ \overrightarrow j}$$\end{document} and \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$\overrightarrow {j + 1}$$\end{document} and the part of the chromosome to the left of \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $${ \overrightarrow j}$$\end{document} is already sorted, m + 1 must be to the right of \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$\overrightarrow {j + 1}$$\end{document} ).
5. If m and m + 1 have different orientation, we can use a reversal to move m + 1 next to m; this will also change the orientation of \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$\overrightarrow {j + 1}$$\end{document} , so in the next step, we can use another reversal to move \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$\overleftarrow {j + 1}$$\end{document} next to \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $${ \overrightarrow j}$$\end{document} .

6. Finally, if m and m + 1 have the same orientation, we interchange blocks
\documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} \begin{align}{ \overrightarrow j} , \left [ \ldots , { \overrightarrow m} \right ] , \ldots , \left [ \overrightarrow {j + 1} , \ldots \right ] , \overrightarrow {m + 1} \quad \rightsquigarrow \quad { \overrightarrow j} , \overrightarrow {j + 1} , \ldots , { \overrightarrow m} , \overrightarrow {m + 1} ,\end{align}\end{document}

if both \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $${ \overrightarrow m}$$\end{document} and \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$\overrightarrow {m + 1}$$\end{document} have positive direction and \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} \begin{align}{ \overrightarrow j} , \left [ \ldots \right ] , { \overleftarrow m} , \ldots , \left [ \overrightarrow {j + 1} , \ldots , \overleftarrow {m + 1} \right ] \quad \rightsquigarrow \quad { \overrightarrow j} , \overrightarrow {j + 1} , \ldots , \overleftarrow {m + 1} , { \overleftarrow m} ,\end{align}\end{document}

if \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $${ \overleftarrow m}$$\end{document} and \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$\overleftarrow {m + 1}$$\end{document} have both negative direction. By two operations we move \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$\overrightarrow {j + 1}$$\end{document} to \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $${ \overrightarrow j}$$\end{document} and \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $${ \overrightarrow m}$$\end{document} to \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$\overrightarrow {m + 1}$$\end{document} .

Every step can be implemented in O(log n) time using an extended version of the data structure from Section 3. We need the data structure to support the following operations: (1) Given a marker, find the chromosome that contains it. (2) Perform a DCJ operation. (3) Given interval \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$i , \ldots , j$$\end{document} , find the marker with the highest number on the chromosome between i and j. To support this query, we store the highest number in the subtree in each node. ▪

3.1. Perfect DCJ scenarios

Bérard et al. (2009) studied the problem of finding a scenario transforming genome Π into Γ that does not break a given set of common intervals. An interval in genome Π is a set of markers such that the subgraph of G_Π induced by their extremities is connected. Intervals of Π have zero or two borders—adjacencies such that one extremity is inside and one outside. Let I be any set of markers with zero or two borders. A DCJ operation preserves I, if I still has zero or two borders in the resulting genome. Bérard et al. (2009) showed that, for nested sets of common intervals (when the intervals do not overlap), the shortest scenario can be found in polynomial time and, for weakly separable sets, the problem is NP-hard but fixed parameter tractable.

Since their algorithms use algorithms for DCJ distance and sorting as a black box, one can use them in conjunction with our algorithm to get perfect restricted DCJ scenarios.

4. Restricted DCJ Halving

The halving problem was first studied in the reversal model. Solutions were given by El-Mabrouk and Sankoff (2003) and Alekseyev and Pevzner (2007). In the DCJ model, the problem was studied by Warren and Sankoff (2009) and then corrected and simplified by Mixtacki (2008). However, the restricted version of the halving problem has not been studied and is stated as open by Tannier et al. (2009).

The simple approach that works for sorting—perform a DCJ operation, test whether a circular chromosome was created, and reincorporate it—does not work for halving: In some cases, the circular chromosome cannot be reincorporated. For example, take chromosome (1₁, 1₂, 2₁, 2₂)—after excision of circular chromosome [1₁, 1₂] it, is not possible to reincorporate it, and the algorithm of Mixtacki (2008) ends with two (perfectly duplicated) circular chromosomes [1₁, 1₂] and [2₁, 2₂]. On the other hand, by fission and translocation we can get \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} \begin{align}( 1_1 \ \mid \ 1_2 , 2_1 , 2_2 ) \quad \rightsquigarrow \quad ( 1_1 \ \mid \ ) , ( 1_2 , 2_1 \ \mid \ 2_2 ) \quad \rightsquigarrow \quad ( 1_1 , 2_2 ) , ( 1_2 , 2_1 ) .\end{align}\end{document}

Kováč et al. (2010) gave the first algorithm for the restricted halving problem running in time \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$O ( n \sqrt n )$$\end{document} . Here, we improve on this result by giving a new algorithm finding a multilinear doubled ancestor in linear time. A particular restricted halving scenario can then be found in O(n log n) time. This also proves that the halving distance is always the same in the restricted and the unrestricted case, and can be computed in linear time from the natural graph (Mixtacki, 2008).

Our new algorithm works in three steps:
1. First, we use the halving algorithm by Mixtacki that computes a doubled genome in linear time.

2. In general, the result may contain circular chromosomes. Thus, in the next step, we transform it into a multilinear doubled genome. We prove that this can be done in linear time while preserving optimality of the result.

3. Finally, if a particular halving scenario is needed, we can use the sorting algorithm from Section 3 to obtain a scenario transforming the given genome into a doubled genome.

Theorem 4

Given a duplicated genome Π and a result of its halving, a doubled genome Θ, we can find a multilinear doubled genome Θ′ within the same DCJ distance from Π in linear time.

Proof

For convenience, assume that Π and Θ are properly capped. The capping procedure is similar to the one described in Section 3. However, we will treat the corresponding caps at the beginning of two paralogous chromosomes in Θ as paralogs, and for each added chromosome consisting of two caps, we add a paralogous copy to both Π and Θ. This way we ensure that
1. both Π and Θ contain exactly two copies of each marker,

2. Θ is a doubled genome (in particular, every linear chromosome has an exact copy in Θ),

3. the distance between Π and Θ is preserved, and

4. for each (linear) chromosome in Π, there is a linear chromosome in Θ beginning with the same marker.

The main idea of our “linearization” algorithm is to search for adjacencies {x, y} in Π such that in Θ, x belongs to a linear chromosome, while y belongs to a circular chromosome. We will call such an adjacency {x, y} an integrating adjacency, since if {x, p} and {y, q} are the adjacencies in Θ, and if we replace \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} \begin{align} & \{x , p\}, \{y, q \} \quad {\rm by} \quad \{x, y\}, \{p, q\} \ {\rm and} & (1) \\ & \{{\overline x} , { \overline p} \} , \{ { \overline y} , { \overline q} \} \quad {\rm by} \quad \{ { \overline x} , { \overline y} \} , \{ { \overline p} , { \overline q} \} , & (2) \end{align} \end{document}

we incorporate the circular chromosome (and its copy) into the linear chromosome (and its copy). Thus, our algorithm searches for integrating adjacencies in Π and incorporates circular chromosomes.

It remains to be proven that:
1. Replacing the adjacencies results in a doubled genome (i.e., the result is still perfectly duplicated) with fewer circular chromosomes.

2. Replacing the adjacencies preserves the DCJ distance.

3. If there is no integrating adjacency in Π, there are no circular chromosomes in Θ.

4. The algorithm can be implemented in linear time.

Figure 6 shows a picture proof of the first claim.

FIG. 6.
We transform a doubled genome by incorporating circular chromosomes (a, b). We cut adjacencies {x, p}, {y, q} and replace them by {x, y}, {p, q} (bold). We do the same with the paralogous copies. There are two cases, depending on whether {y, q} and \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$\{{\overline y}, {\overline q}\}$$\end{document} are in two paralogous chromosomes (left), or in one doubled chromosome (right). Note that in both cases, the resulting genome is a doubled genome with fewer circular chromosomes.

To see the second claim, let Θ be the genome before, and Θ′ the genome after the two replacements. Note that the replacements (1) and (2) are actually two DCJ operations, so the difference between d(Π, Θ′) and d(Π, Θ) cannot be more than 2. Furthermore, note that the first operation decreases the distance between Π and Θ, so d(Π, Θ′) ≤ d (Π, Θ). However, since Θ is optimal, we also have d(Π, Θ) ≤ d(Π, Θ′).

Let m be a marker belonging to a circular chromosome in Θ and let L be the chromosome in Π containing m. We can write L as \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$L = (m_1, m_2, m_3, \ldots, m, \ldots)$$\end{document} . Thanks to capping, m₁ belongs to a linear chromosome in Θ. Let m_k₊₁ be the first marker of L belonging to a circular chromosome in Θ, then there is an adjacency between m_k and m_k₊₁ in Π such that the first extremity belongs to a linear chromosome and the second one belongs to a circular chromosome in Θ. In other words, if Θ has a circular chromosome, there is always an integrating adjacency in Π. This proves our third claim.

We implement the algorithm as follows: First, we traverse all chromosomes in Θ and populate a lookup table that stores the type of chromosome (linear or circular) for each extremity. Let ℓ be the list of all extremities in linear chromosomes of Θ (this list will grow as we incorporate circular chromosomes). For each extremity x in ℓ, let {x, y} be an adjacency in Π and let C be the chromosome in Θ containing y. If {x, y} is integrating, i.e., if C is a circular chromosome, we set the chromosome type of each extremity of C to linear in the lookup table, we append these extremities to the list ℓ and incorporate the circular chromosome(s) by replacing adjacencies as in (1) and (2).

The algorithm obviously runs in linear time, since for each extremity, we determine only once whether it belongs to a linear or a circular chromosome, we change this information at most once (when a circular chromosome is incorporated in a linear one), and once the extremity belongs to a linear chromosome, we check at most once whether the extremity adjacent in Π belongs to a circular chromosome. ▪

Corollary 1

The halving distance is the same in the DCJ and in the restricted DCJ model. The distance and one multilinear doubled genome can be computed in linear time. A restricted halving scenario can be computed in time O(n log n).

5. Restricted DCJ Median

In sorting and halving problems, the reincorporation restriction did not change the distance, and we were just searching for restricted scenarios that better explain the evolutionary history. This is not the case for the median problem: Consider three linear genomes: (1, 2, 3), (2, 1, 3), and (2, 3, 1). Their median in the unrestricted case consists of linear chromosome (2, 3) and circular chromosome [1]. Since we can obtain any of the linear genomes by a single operation (reincorporation of chromosome [1]), the median score is 3. This score, however, cannot be achieved in the restricted model. Generalizing this pattern, we can get genomes of length 3n with unrestricted median score 3n and restricted median score 4n.

Caprara (2003) proved that finding a median in a DCJ model restricted to unichromosomal linear genomes1 is NP-hard. Using the method of Caprara (2003), Tannier et al. (2009) were able to prove an NP-hardness result also for the general model. The median problem in the restricted DCJ model was not studied; it was conjectured NP-hard, but stated as open in Tannier et al. (2009). We confirm this conjecture and show that the NP-hardness follows easily from the result of Caprara (2003):

Theorem 5

The median problem is NP-hard in the restricted DCJ model.

Proof

We show that if we could solve the multichromosomal linear DCJ median efficiently, we could also solve the unichromosomal linear DCJ median efficiently: Let Π₁, Π₂, Π₃ be the given (unichromosomal) linear genomes and let Π_M be their multilinear median. We show that we can join telomeres of Π_M and fuse chromosomes into a single chromosome \documentclass{aastex} \usepackage{amsbsy} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{bm} \usepackage{mathrsfs} \usepackage{pifont} \usepackage{stmaryrd} \usepackage{textcomp} \usepackage{portland, xspace} \usepackage{amsmath, amsxtra} \pagestyle{empty} \DeclareMathSizes {10} {9} {7} {6} \begin{document} $$\Pi_M^{\prime}$$\end{document} while preserving the median score.

Consider two linear chromosomes in Π_M. Each has two telomeres, so we can join them in four different ways. What happens to the distance d(Π_M, Π_i), if we join two telomeres? Recall that the DCJ distance between two genomes is calculated as n − (c + p_o/2), where c is the number of cycles and p_o is the number of odd length paths in their adjacency graph.
1. If the telomeres belong to two different paths of even length, the distance remains unchanged.

2. If the telomeres belong to the same path of even length, by joining them, we create a cycle and the distance decreases.

3. Finally, if the telomeres belong to two (different) odd length paths, by joining the telomeres, we get a path of even length and the distance increases. This is the bad case.

However, since each of the given genomes Π₁, Π₂, Π₃ has only two telomeres, the adjacency graph AG(Π_i, Π_M) may contain at most two paths of odd length. That is, if there is a bad case, all the other three possibilities of joining the telomeres are good (do not change the distance or even decrease it). Since there are four ways of joining the telomeres, but for each genome Π₁, Π₂, Π₃, only at most one way is bad, we can always fuse the linear chromosomes. ▪

6. Conclusion

In this work, we revisited the restricted DCJ model for multichromosomal linear genomes, where a temporary circular chromosome is immediately reincorporated after its excision. We improved on the quadratic-time algorithm by Yancopoulos et al. (2005) and proposed an algorithm that runs in O(n log n) time.

Furthermore, we have solved an open problem from Tannier et al. (2009) by giving an algorithm for the restricted halving problem. The algorithm shows that the halving distance for the restricted version is the same as the distance for the unrestricted version and given a multilinear duplicated genome, an optimal multilinear perfectly duplicated genome can always be found in linear time.

Finally, we confirmed the conjecture from Tannier et al. (2009) and proved that the median problem is NP-hard in the restricted DCJ model.

Footnotes

Acknowledgments

We would like to thank Broňa Brejová for many constructive comments. The research of Jakub Kováč is supported by Marie Curie Fellowship IRG-224885 to Dr. Tomáš Vinař, Comenius University grant UK/121/2011, and a Partnership Stipend from Bielefeld University. A preliminary version of this work appeared in Kováč et al. ()

Disclosure Statement

No competing financial interests exist.

1

Caprara called it the cycle median problem—the DCJ model did not exist yet

References

Alekseyev

, Pevzner

2007. Colored de Bruijn graphs and the genome halving problem. IEEE Trans. Comput. Biol. Bioinform, 98–107.

Bérard

, Chateau

, Chauve

et al. 2009. Computation of perfect DCJ rearrangement scenarios with linear and circular chromosomes. J. Comput. Biol, 16:1287–1309.

Bergeron

, Mixtacki

, Stoye

2006. A unifying view of genome rearrangements. Proc. WABI, 163–173.

Caprara

2003. The reversal median problem. INFORMS J. Comp., 15:93.

Christie

D.A.

1996. Sorting permutations by block-interchanges. Inf. Process. Lett., 60:165–169.

Chrobak

, Szymacha

, Krawczyk

1990. A data structure useful for finding hamiltonian cycles. Theor. Comput. Sci., 71:419–424.

El-Mabrouk

, Sankoff

2003. The reconstruction of doubled genomes. SIAM J. Comput., 32:754–792.

Feng

, Zhu

2007. Faster algorithms for sorting by transpositions and sorting by block interchanges. ACM Trans.Algorithms 3.

Kaplan

, Verbin

2005. Sorting signed permutations by reversals, revisited. J. Comput. Syst. Sci., 70:321–341.

10.

Kováč

, Braga

M.D.V.

, Stoye

2010. The problem of chromosome reincorporation in DCJ sorting and halving. Lect. Notes Comput. Sci., 6398:13–24.

11.

Mixtacki

2008. Genome halving under DCJ revisited. Lect. Notes Comput. Sci., 5092:276–286.

12.

Swenson

K.M.

, Rajan

, Lin

et al. 2010. Sorting signed permutations by inversions in O(n log n) time. J. Comput. Biol., 17:489–501.

13.

Tannier

, Zheng

, Sankoff

2009. Multichromosomal median and halving problems under different genomic distances. BMC Bioinform., 10:120.

14.

Warren

, Sankoff

2009. Genome halving with double cut and join. J. Bioinform. Comput. Biol., 7:357–371.

15.

Yancopoulos

, Attie

, Friedberg

2005. Efficient sorting of genomic permutations by translocation, inversion and block interchange. Bioinformatics, 21:3340–3346.