 Research article
 Open Access
 Published:
Analysis of circular genome rearrangement by fusions, fissions and blockinterchanges
BMC Bioinformatics volume 7, Article number: 295 (2006)
Abstract
Background
Analysis of genomes evolving via blockinterchange events leads to a combinatorial problem of sorting by blockinterchanges, which has been studied recently to evaluate the evolutionary relationship in distance between two biological species since blockinterchange can be considered as a generalization of transposition. However, for genomes consisting of multiple chromosomes, their evolutionary history should also include events of chromosome fusions and fissions, where fusion merges two chromosomes into one and fission splits a chromosome into two.
Results
In this paper, we study the problem of genome rearrangement between two genomes of circular and multiple chromosomes by considering fusion, fission and blockinterchange events altogether. By use of permutation groups in algebra, we propose an $\mathcal{O}$(n^{2}) time algorithm to efficiently compute and obtain a minimum series of fusions, fissions and blockinterchanges required to transform one circular multichromosomal genome into another, where n is the number of genes shared by the two studied genomes. In addition, we have implemented this algorithm as a web server, called FFBI, and have also applied it to analyzing by gene orders the whole genomes of three human Vibrio pathogens, each with multiple and circular chromosomes, to infer their evolutionary relationships. Consequently, our experimental results coincide well with our previous results obtained using the chromosomebychromosome comparisons by landmark orders between any two Vibrio chromosomal sequences as well as using the traditional comparative analysis of 16S rRNA sequences.
Conclusion
FFBI is a useful tool for the bioinformatics analysis of circular and multiple genome rearrangement by fusions, fissions and blockinterchanges.
Background
For the past two decades, genome rearrangements have been studied and can be modelled to learn more about the evolution of mitochondrial, chloroplast, viral, bacterial and mammalian genomes [1]. To evaluate the evolutionary distance between two related genomes in gene order, various rearrangement events acting on genes within or among chromosomes have been proposed, such as reversals (also known as inversions) [1–10], transpositions [11, 12], blockinterchanges [13–15], translocations [16, 17], and fusions and fissions [18, 19]. Most genome rearrangement studies in computation involve the issue of solving the combinatorial problem to find an optimal series of rearrangements required to transform one genome into another.
Recently, the study on the genome rearrangement by blockinterchanges has increasingly drawn great attention, since the blockinterchange event is a generalization of transposition and, currently, its computational models measuring the genetic distance are more tractable than those modeled by transposition. Christie [13] first introduced the concept of blockinterchange, affecting a chromosome by swapping two nonintersecting blocks containing any number of consecutive genes. Blockinterchange can be considered as a generalization of transposition, since any exchanged blocks via transposition must be contiguous in a chromosome, whereas those via blockinterchange need not be. As a matter of fact, the occurrence of an exchange of two noncontiguous blocks has been suggested in the previous studies related to the biological processes of bacterial replication [[20], and references therein]. Christie also proposed an $\mathcal{O}$(n^{2}) time algorithm, where n is the number of genes, to solve the socalled blockinterchange distance problem that is to find a minimum series of blockinterchanges for transforming one linear chromosome into another. Later, we [14] designed a simpler algorithm for solving the blockinterchange problem on linear or circular chromosomes with timecomplexity of $\mathcal{O}$(δ n), where δ is the the minimum number of blockinterchanges required for the transformation and can be calculated in $\mathcal{O}$(n) time in advance. We also demonstrated that blockinterchange events play a significant role in the genetic evolution of bacterial (Vibrio) species. Very recently, based on this algorithm, we have further implemented a tool, called ROBIN, for analyzing the rearrangements of gene orders via blockinterchanges between two linear/circular chromosomal genomes [15]. Not only geneorder data but also sequence data are allowed to be input into the ROBIN system. If the input is the sequence data, ROBIN can automatically search for the common homologous/conserved regions shared by all input sequences.
It should be noted that the above blockinterchange studies were dedicated to genomes containing only one chromosome (i.e., unichromosomal genomes) for evaluating their evolutionary relationships. However, for biological species with different numbers of chromosomes, the evolutionary history must also consider events of chromosomal fusions and fissions. A fusion occurs when two chromosomes merge into one and a fission takes place when a chromosome splits into two. The reason is that different chromosomes may as well exchange their genetic material with each other and, moreover, this exchange can only be achieved via interchromosomal operations such as fusions and fissions, instead of intrachromosomal operations like blockinterchanges. Hence, it is worthwhile to study genome rearrangements considering fusions, fissions and blockinterchanges altogether. In this paper, we solve such a genome rearrangement problem by designing an efficient algorithm to compute and obtain a minimum series of all the events involving fusions, fissions and blockinterchanges that are required to transform one circular multichromosomal genome into another, when both have the same set of genes without repeats. Although most eukaryotic genomes are linear, most prokaryotic (e.g., bacterial) genomes are circular and some of them consist of multiple circular chromosomes and large plasmids^{1}. For example, some important bacterial pathogens like Brucella, Burkholderia, Leptospira and Vibrio species fall into this category. Notably, our approach is based on permutation group in algebra, instead of breakpoint graph, a commonly used approach in the study of genome rearrangement.
Recently, Yancopoulos et al. [21] used breakpoint graph to design an algorithm to solve a genome rearrangement problem in which the considered reversals, translocations (including fusions and fissions) and blockinterchanges were given different weights. Unfortunately, their algorithm cannot be applied to solving our problem in which the events we considered are unweighted, because a series of weighted events with minimum weights in total may not be a minimum series of unweighted events, provided the events are given different weights.
Results and discussion
Based on Algorithm Sortingbyffbi developed in this study, we have implemented a web server, called FFBI [22], in which biologists or scientists in genomics can conduct comprehensive analyses of circular genome rearrangements by fusions, fissions and blockinterchanges for their scientific interests and needs. Furthermore, we used this web server to conduct the rearrangement analyses on the whole genomes of three pathogenic Vibrio species, including V. vulnificus, V. parahaemolyticus and V. cholerae, to infer their evolutionary relationships.
Each of these three Vibrio pathogens consists of two circular chromosomes, and all their genomic sequences have recently been reported in GenBank with proteincoding genes annotated (see Table 1 for their sequence information). As annotated in GenBank (as of April 2006), the genomes of V. vulnificus, V. parahaemolyticus and V. cholerae contain 5098, 4992 and 4008 genes, respectively. From these proteincoding genes, we identified a total of 2393 (onetoone) orthologous genes that are physically located in different positions on the chromosomes (see the Method section for construction of orthologous genes). Inevitably, there can be a high possibility that some genes with misannotated or uncertain protein functions are included in the genome annotation data. We therefore used only those authentic genes whose protein functions are not annotated as hypothetical or putative proteins, or are conserved and not poorly characterized (e.g., not those genes with only general function prediction or unknown function) in the NCBI COGs [23] database of orthologous genes. As a result, there are 1274 authentic orthologous genes in total remained for the further study of genome rearrangement. The relative orders of these orthologous genes along chromosomes, as well as the annotated COGs of their coding proteins, are detailed in the web site of our server.
For each pair of these pathogenic Vibrio species, the variation in their gene orders has suggested that the genome rearrangement events have occurred and their genomes are closely related in evolution. To evaluate the contribution of fusions, fissions and blockinterchanges to these observed rearrangements, we used the server developed in this study to compute the rearrangement distance between the gene orders of any two Vibrio genomes. Consequently, as shown in Table 2, the calculated rearrangement distance between V. vulnificus and V. parahaemolyticus is smaller than that between V. vulnificus and V. cholerae and that between V. parahaemolyticus and V. cholerae, suggesting that V. vulnificus is closer to V. parahaemolyticus than to V. cholerae in evolutionary relationship. Intriguingly, this result of genomewide experiment well coincides with those we obtained in the previous chromosomewide experiments [14, 15]. Recall that our previous experiments were conducted in a chromosomebychromosome style of computing the blockinterchange distance between the landmark orders of any two large/small Vibrio chromosomes, where the used landmarks are the maximal unique matches (MUMs) or the locally collinear blocks (LCBs) that are commonly shared by three large/small Vibrio chromosomes.
In fact, the evolutionary relationships of the three pathogenic Vibrio species revealed in our experiment of analyzing their genome rearrangements also confirms that obtained by the biological community on the basis of the traditional comparative analysis of 16S rRNA gene sequences [24–26]. For confirmation, we here repeated this comparative analysis as follows. The 16S rRNA gene sequences of three Vibrios were first aligned using the Clustal W program [27], from which the distance matrix (as shown in Table 3) was then estimated by the algorithm of Kimura's twoparameter model in PHYLIP package [28].
Conclusion
In this paper, we studied the genome rearrangement problem between circular genomes with multiple chromosomes by simultaneously considering fusion, fission and blockinterchange events. We have shown in the Method section that an optimal series of events required to transform one genome into another can be obtained in a canonical order such that all fusions come before all blockinterchanges, which come before all fissions. Based on this property as well as the concept of permutation groups in algebra, we have successfully designed an $\mathcal{O}$(n^{2}) time algorithm to obtain the minimum number of fusion, fission and blockinterchange events for the transformation and also to generate an optimal scenario of the required rearrangement events. In addition, we have practically implemented this algorithm as a web server and applied it to analyzing by gene orders the whole genomes of three human Vibrio pathogens to infer their evolutionary relationships. As a consequence, our experimental results well coincide with the previous results obtained using the chromosomebychromosome comparisons by landmark orders between any two Vibrio chromosomal sequences as well as using the traditional comparative analysis of 16S rRNA sequences. The algorithm, however, should not be applied on linear multichromosomal genomes, because as mentioned in the Method section, it is not always possible to have an optimal scenario in a canonical order for linear genomes. Further studies in genome rearrangement can still be pursued to solve the problem for linear multichromosomal genomes.
Methods
Permutations versus genome rearrangements
In group theory, a permutation is defined to be a onetoone mapping from a set E = {1, 2, ..., n} into itself, where n is some positive integer. For example, we may define a permutation α of the set {1, 2, 3, 4, 5, 6, 7} by specifying α(1) = 4, α(2) = 3, α(3) = 1, α(4) = 2, α(5) = 7, α(7) = 6 and α(6) = 5. The above mapping can be expressed using a cycle notation as illustrated in Figure 1 and simply denoted by α = (1, 4, 2, 3) (5, 7, 6). A cycle of length k, say (a_{1}, a_{2}, ..., a_{ k }), is simply called kcycle and can be rewritten as (a_{ i }, a_{i+1}, ..., a_{ k }, a_{1}, ..., a_{i1}), where 2 ≤ i <k, or (a_{ k }, a_{1}, a_{2}, ..., a_{k1}). Any two cycles are said to be disjoint if they have no element in common. In fact, any permutation, say α, can be written in a unique way as the product of disjoint cycles, which is called the cycle decomposition of α, if we ignore the order of the cycles in the product [29]. Usually, a cycle of length one in α is not explicitly written and its element, say x, is said to be fixed by α since α(x) = x. Especially, the permutation whose elements are all fixed is called an identity permutation and is denoted by 1 (i.e., 1 = (1) (2) ... (n)).
Given two permutations α and β of E, the composition (or product) of α and β, denoted by α β, is defined to be a permutation of E with α β(x) = α(β(x)) for all x ∈ E. For instance, if we let E = {1, 2, 3, 4, 5, 6}, α = (2, 3) and β = (2, 1, 5, 3, 6, 4), then α β = (2, 1, 5) (3, 6, 4). If α and β are disjoint cycles, then α β = β α. The inverse of α is defined to be a permutation, denoted by α^{1}, such that α α^{1} = α^{1}α = 1. If a permutation is expressed by the product of disjoint cycles, then its inverse can be obtained by just reversing the order of the elements in each cycle. For example, if α = (2, 1, 5) (3, 6, 4), then α^{1} = (5, 1, 2) (4, 6, 3). Clearly, α^{1} = α if α is a 2cycle.
Meidanis and Dias [19, 30] first noted that each cycle of a permutation may represent a circular chromosome of a genome with each element of the cycle corresponding to a gene, and the order of the cycle corresponding to the gene order of the chromosome. Figure 1, for example, shows a genome with two circular chromosomes, one represented by (1, 4, 2, 3) and the other by (5, 7, 6). Moreover, they observed that global evolutionary events, such as fusions and fissions (respectively, transpositions), correspond to the composition of a 2cycle (respectively, 3cycles) and the permutation representing a genome. For instance, let α be any permutation whose cycle decomposition is α_{1}α_{2} ... α_{ r }. If ρ = (x, y) is a 2cycle and x and y are in the different cycles of α, say α_{ p }= (a_{1} ≡ x, a_{2}, ..., a_{ i }) and α_{ q }= (b_{1} ≡ y, b_{2}, ..., b_{ j }) where 1 ≤ p, q ≤ r, then in the composition ρ α, α_{ p }and α_{ q }are joined into a cycle (a_{1}, a_{2}, ..., a_{ i }, b_{1}, b_{2}, ..., b_{ j }), i.e., ρ is a fusion event affecting on α (and also called a join operation of α here). If ρ = (x, y) is a 2cycle and x and y are in the same cycle of α, say α_{ p }= (a_{1} ≡ x, a_{2}, ..., a_{ i }≡ y, a_{i+1}, ..., a_{ j }) where 1 ≤ p ≤ r, then in the composition ρ α, this cycle α_{ p }is broken into two disjoint cycles (a_{1}, a_{2}, ..., a_{i1}) and (a_{ i }, a_{i+1}, ..., a_{ j }), i.e., ρ is a fission event affecting on α (and also called a split operation of α here). If ρ = (x, y, z) is a 3cycle and x, y and z are in the same cycle of α, say α_{ p }= (a_{1} ≡ x, a_{2}, ..., a_{ i }, b_{1} ≡ y, b_{2}, ..., b_{ j }, c_{1} ≡ z, c_{2}, ..., c_{ k }) where 1 ≤ p ≤ r, then in the composition ρ α, the cycle α_{ p }becomes (a_{1}, a_{2}, ..., a_{ i }, c_{1}, c_{2}, ..., c_{ k }, b_{1}, b_{2}, ..., b_{ j }), i.e., ρ is a transposition event affecting on α.
Recently, Lin et al. [14] further observed that a blockinterchange event affecting on α corresponds to the composition of two 2cycles, say ρ_{1} and ρ_{2}, and α under the condition that ρ_{1} is a split operation of α and ρ_{2} is a join operation of ρ_{1}α. More clearly, let α_{ p }= (a_{1}, a_{2},..., a_{ k }) be a cycle of α, ρ_{1} = (a_{1}, a_{ i }) and ρ_{2} = (a_{ h }, a_{ j }), where 1 <i ≤ k,1 ≤ h ≤ i  1 and i ≤ j ≤ k. Then ρ_{2}ρ_{1}α is the resulting permutation by exchanging the blocks [a_{ h }, a_{i1}] and [a_{ j }, a_{ k }] of α_{ p }. That is, ρ_{2}ρ_{1} is a blockinterchange event affecting on α by swapping [a_{ h }, a_{i1}] and [a_{ j }, a_{ k }], two nonintersecting blocks in α.
As discussed above, any series of fusions, fissions and blockinterchanges required to transform one circular multichromosomal genome α into another I can be expressed by a product of 2cycles, say ρ_{ k }ρ_{k1}...ρ_{1}, such that ρ_{ k }ρ_{k1}...ρ_{1}α = I (hence, ρ_{ k }ρ_{k1}...ρ_{1} = I α^{1}). This property implies that I α^{1} contains all information that can be utilized to derive ρ_{1}, ρ_{2},...,ρ_{ k }for transforming α into I.
It is well known that every permutation can be written as a product of 2cycles. For example, (1, 2, 3, 4) = (1, 4) (1, 3) (1, 2). However, there are many ways of expressing a permutation α as a product of 2cycles [29]. Given a permutation α, let f(α) denote the number of the disjoint cycles in the cycle decomposition of α. Notice that f(α) counts also the nonexpressed cycles of length one. For example, if α = (1, 5) (2, 4) is a permutation of E = {1, 2,..., 5}, then f(α) = 3, instead of f(α) = 2, since α = (1, 5) (2, 4) (3). Then the following lemma shows the lower bound of the number of 2cycles in any product of 2cycles of a permutation.
Lemma 1 Let α be an arbitrary permutation of E = {1, 2,..., n}. If α can be expressed as a product of m 2cycles, say α = α_{1}α_{2}...α_{ m }with each α_{ i } being 2cycle, then m ≥ n  f(α).
Proof. We prove this lemma by induction on m. The lemma is true if m = 0, since α = 1 then, meaning that f(α) = n, and hence m = n  f(α) = 0. Suppose now that the lemma holds for any permutation that can expressed as a product of less than m 2cycles, where m > 0. Let α' = α_{2}α_{3}, ...,α_{ m }. Then by the induction hypothesis, we have m  1 ≥ n  f(α'). Since α = α_{1}α' and α_{1} is a 2cycle, α_{1} operates on α' either as a fusion by joining two cycles of α' into one cycle (i.e., f(α) = f(α')  1) or as a fission by splitting one cycle of α' into two cycles (i.e., f(α) = f(α') + 1). Whichever α_{1} operates on α', we have f(α) ≥ f(α')  1. As a result, m = (m  1) + 1 ≥ n  f(α') + 1 = n  (f(α')  1) ≥ n  f(α). □
Optimal scenario in canonical order
As mentioned previously, each circular multichromosomal genome with n genes can be expressed by a permutation of E = {1, 2,..., n}. Given two such genomes G_{1} and G_{2} over the same gene set E, the genome rearrangement distance between G_{1} and G_{2}, denoted by d(G_{1},G_{2}), is defined to be the minimum number of events needed to transform G_{1} into G_{2}, where the events allowed to take place are fusions, fissions and blockinterchanges. In this section, we shall show that there is an optimal series of events required to transform G_{1} into G_{2} such that all fusions come prior to all blockinterchanges, which come before all fissions. Here, such an optimal scenario of genome rearrangements is referred as in canonical order.
Lemma 2 d(G_{1}, G_{2}) = d(G_{2}, G_{1}).
Proof. Let Φ = <σ_{1}, σ_{2},...,σ_{ δ }> be an optimal series of events required to transform G_{1} into G_{2}. Clearly, Φ' = <σ_{ δ }, σ_{δ1},...,σ_{1}> is an optimal series of events for transforming G_{2} into G_{1} by reversing the role of every event σ_{ i }, where 1 ≤ i ≤ δ, such that σ_{ i } is a fission (respectively, fusion) in Φ' if σ_{ i }is a fusion (respectively, fission) in Φ.
Lemma 3 There is an optimal series of events required to transform G_{1} into G_{2} such that every fission occurs after every fusion and blockinterchange.
Proof. Let Φ = <σ_{1}, σ_{2},...,σ_{ δ }> be an optimal series of events needed to transform G_{1} into G_{2}. Of course, if every fission occurs after every fusion and blockinterchange in Φ, then the proof is done. Now, suppose that not every fission occurs after every fusion or blockinterchange in Φ. Then let i be the largest index in Φ such that σ_{ i }is a fission preceding σ_{i+1}that is either a fusion or a blockinterchange. We can then obtain a new optimal series Φ' = <σ_{1},...,σ_{i1}, ${{\sigma}^{\prime}}_{i}$, ${{\sigma}^{\prime}}_{i+1}$, σ_{i+2},...,σ_{ δ }> to transform G_{1} into G_{2} such that ${{\sigma}^{\prime}}_{i}$ is a fusion or a blockinterchange and ${{\sigma}^{\prime}}_{i+1}$ is a fission, as discussed below. Suppose that σ_{ i } splits a chromosome α into α_{1} and α_{2}. If σ_{i+1}is a fusion, then we assume that it joins two chromosomes β_{1} and β_{2}into β; otherwise, if σ_{i+1}is a blockinterchange, then assume that it affects β_{1} such that β_{1} becomes β through a blockinterchange. Clearly, if neither β_{1} nor β_{2} is created by σ_{ i }, then the desired series Φ' is obtained by swapping σ_{ i } and σ_{i+1}in Φ (i.e., ${{\sigma}^{\prime}}_{i}$ = σ_{i+1}and ${{\sigma}^{\prime}}_{i+1}$ = σ_{ i }). If both β_{1} and β_{2} are created by σ_{ i }, then the net rearrangement of σ_{ i } (a split operation) followed by σ_{i+1}(a joint operation) either has no effect on α or becomes a blockinterchange affecting α. By removing σ_{ i }and σ_{i+1}from Φ or replacing them with an extra blockinterchange, we thus obtain a new optimal series of the events transforming G_{1} into G_{2} with strictly less than δ events, a contradiction. Hence, we assume that only one of β_{1} and β_{2} is created by σ_{ i } and without loss of generality, let β_{1} = α_{1}. Now, we consider the following two cases.
Case 1: σ_{i+1}is a fusion. For simplicity of discussion, we let α = (1, 2,..., x  1, x,..., y  1), β_{2} = (y, y + 1,..., z), σ_{ i }= (1, x) and σ_{i+1}= (1, y), where 1 <x <y  1 and y <z. Then the net rearrangement caused by σ_{ i } and σ_{i+1}is to transform α and β into α_{2} = (x,x+1,...,y  1) and β = (1, 2,..., x  1, y, y + 1,..., z). In fact, this rearrangement can also be done by first joining α and β_{2} into (1, 2,..., z) via σ_{i+1}and then splitting (1, 2,..., z) into α_{2} and β via (x, y). In other words, we can obtain Φ' by letting ${{\sigma}^{\prime}}_{i}$ = σ_{i+1}and ${{\sigma}^{\prime}}_{i+1}$ = (x, y).
Case 2: σ_{i+1}is a blockinterchange. Clearly, the net rearrangement caused by σ_{ i } and σ_{i+1}is to transform α into β and α_{2}, which is equivalent to the rearrangement by first applying σ_{i+1}to α and then further splitting it into β and α_{2} via σ_{ i }. Then Φ' is obtained by swapping σ_{ i } and σ_{i+1}in Φ.
In other words, we can always obtain Φ' from Φ according to the method described above. Repeating this process on the resulting Φ', we can finally obtain an optimal series of events that are required to transform G_{1} into G_{2} such that all fissions come after all fusions and blockinterchanges.
Lemma 4 There is an optimal series of events required to transform G_{1} into G_{2} in a canonical order such that all fusions come before all blockinterchanges, which come before all fissions.
Proof. Let Φ = <σ_{1}, σ_{2},..., σ_{ δ }> be an optimal series of events required to transform G_{1} into G_{2}. If there are no fusions or blockinterchanges, then the proof is completed. If not, according to Lemma 3, we may assume that all fusions and blockinterchanges occur earlier than all fissions. Let i be the index of the last nonfission in Φ and also let G' be the resulting genome after all σ_{1}, σ_{2},...σ_{ i }have affected G_{1}. Since Φ is optimal, it is straightforward to see that Φ_{ i }= <σ_{1}, σ_{2},...,σ_{ i }> is an optimal series of fusions and blockinterchanges needed to transform G_{1} into G'. As discussed in the proof of Lemma 2, ${{\Phi}^{\prime}}_{i}$ = <σ_{ i }, σ_{i1},...,σ_{1}> is an optimal series of fissions and blockinterchanges for transforming G' into G_{1}. Moreover, by Lemma 3, we can obtain ${{{\Phi}^{\prime}}^{\prime}}_{i}$ = <${{\sigma}^{\prime}}_{i}$, ${{\sigma}^{\prime}}_{i1}$,...,${{\sigma}^{\prime}}_{1}$> from ${{\Phi}^{\prime}}_{i}$ for transforming G' into G_{1} such that all block interchanges in ${{{\Phi}^{\prime}}^{\prime}}_{i}$ occur prior to all fissions. Consequently, <${{\sigma}^{\prime}}_{1}$, ${{\sigma}^{\prime}}_{2}$,...,${{\sigma}^{\prime}}_{i}$> is an optimal series of fusions and blockinterchanges needed to transform G_{1} into G', and all its fusions occur before all its blockinterchanges. Therefore, there is an optimal series of events needed to transform G_{1} into G_{2} such that all fusions come earlier than all blockinterchanges, which come before all fissions. □
It is worth mentioning that an optimal scenario in a canonical order does not necessarily exist for linear multichromosomal genomes. For example, suppose that G_{1} and G_{2} are two given linear multichromosomal genomes, where G_{1} = (1, 4, 5) (2, 3) and G_{2} = (1, 2, 3) (4, 5). Then the optimal scenario between them is a fission, splitting (1, 4, 5) into (1) (4, 5), followed by a fusion, joining (1) and (2, 3) to (1, 2, 3). However, this optimal scenario can not be transformed into another in the canonical order according to the steps as described in Lemmas 3 and 4. Actually, there is no an optimal scenario between such two linear genomes using any two rearrangement events that begin with a fusion.
Algorithm
Let α and I be two given circular multichromosomal genomes over the same gene set E = {1, 2,...,n}. Here, we assume that the genes in I are sorted in the order of increasing and consecutive numbers, and that gene i + 1 is on the right side of gene i within the same chromosome, where 1 ≤ i ≤ n  1. For example, I = (1, 2) (3, 4, 5) (6, 7, 8, 9) if I has three circular chromosomes with two, three and four genes, respectively. In this case, the computation of d(α,I) and its corresponding optimal scenario can be considered as a problem of sorting α using the minimum set of operations, including fusions, fissions and blockinterchanges.
Suppose that ρ_{ λ }ρ_{λ  1}... ρ_{1} is a product of 2cycles that corresponds to an optimal series Φ of fusions, fissions and blockinterchanges for transforming α into I. Then Φ = ρ_{ λ }ρ_{λ  1}...ρ_{1} = I α^{1} and by Lemma 1, λ ≥ n  f(Iα^{1}). If Φ' = ${{\rho}^{\prime}}_{\lambda}{{\rho}^{\prime}}_{\lambda 1}\dots {{\rho}^{\prime}}_{1}$ is another product of λ 2cycles with Φ' = Iα^{1}, then the number of 2cycles in Φ' that function as fusions or fissions must be greater than or equal to that in Φ; otherwise, the total number of fusions, fissions and blockinterchanges in Φ' for transforming α into I must be less than that in Φ, a contradiction. The reason is that a fusion or fission requires only one 2cycle for rearrangement, whereas a blockinterchange requires two 2cycles. In other words, the number of 2cycles serving as the fusions and fissions is minimum in any optimal series of events.
Based on the above observation as well as Lemma 4, below we design an efficient algorithm for computing d(α,I) and its optimal scenario of rearrangement events in a canonical order. Let χ(α) and χ(I) denote the numbers of chromosomes in α and I, respectively, and let α = α_{1}α_{2}...α_{χ(α)}and I = I_{1}I_{2}...I_{χ(I)}. Then, an undirected graph $\mathcal{G}$(α,I) = ($\mathcal{V}$_{ α }, $\mathcal{V}$_{ I }, $\mathcal{E}$) is constructed from α and I as follows.

$$\mathcal{V}$$(1)
_{ α }= {α_{1}, α_{2},...,α_{χ(α)}}.

$$\mathcal{V}$$(2)
_{ I }= {I_{1}, I_{2},...,I_{χ(I)}}.

$$\mathcal{E}$$(3)
= {(α_{ i }, I_{ j })  1 ≤ i ≤ χ (α), 1 ≤ j ≤ χ (I), and α_{ i }and I_{ j }have at least a common gene}.
For instance, suppose that α = α_{1}α_{2}...α_{5} = (1, 2, 10) (11, 8, 9, 3, 6) (7, 4, 5,12) (13, 15) (14, 16) and I = I_{1}I_{2} ... I_{6} = (1, 2, 3) (4, 5) (6, 7, 8) (9, 10, 11, 12) (13, 14) (15, 16). Then the induced graph $\mathcal{G}$(α,I) is shown in Figure 2, which is a bipartite graph since $\mathcal{V}$_{ α }and $\mathcal{V}$_{ I }are independent sets in $\mathcal{G}$(α, I) (i.e., no edge between any two vertices in $\mathcal{V}$_{ α }or $\mathcal{V}$_{ I }). A connected component of $\mathcal{G}$(α,I) is defined to be a maximal subgraph of $\mathcal{G}$(α,I) such that there exists a path between any pair of vertices in this subgraph. For example, the induced $\mathcal{G}$(α,I) as shown in Figure 2 has two connected components. Notice that if in a chromosome I_{ k }of I there are two genes that appear in two different chromosomes α_{ i } and α_{ j }of α, then (α_{ i }, I_{ k }) ∈ $\mathcal{E}$ and (α_{ j }, I_{ k }) ∈ $\mathcal{E}$, and hence both α_{ i } and α_{ j } belong to the same connected component in $\mathcal{G}$(α,I).
Let {$\mathcal{C}$_{1}, $\mathcal{C}$_{2},...,$\mathcal{C}$_{ ω }} denote the collection of all connected components in $\mathcal{G}$(α,I). For each 1 ≤ i ≤ ω, let β_{ i }and J_{ i } denote the chromosomes in α and I, respectively, whose corresponding vertices belong to $\mathcal{C}$_{ i }in $\mathcal{G}$(α,I). Let gene(β_{ i }) and gene(J_{ i }) be the collections of the genes in all chromosomes of β_{ i }and J_{ i }, respectively. Then gene(β_{ i }) = gene(J_{ i }) and gene(β_{ i }) ∩ gene(β_{ j }) = ∅ for any 1 ≤ j ≠ i ≤ ω. Let n_{ i } be the number of genes in gene(β_{ i }). Clearly, n = n_{1} + n_{2} + ... + n_{ ω }. In addition, it can be verified that Iα^{1} = (J_{1}${\beta}_{1}^{1}$) (J_{2}${\beta}_{2}^{1}$) ... (J_{ ω }${\beta}_{\omega}^{1}$) and f(Iα^{1}) = f(J_{1}${\beta}_{1}^{1}$) + f(J_{2}${\beta}_{2}^{1}$) + ... + f(J_{ ω }${\beta}_{\omega}^{1}$).
According to the properties above, we then find a product Φ of 2cycles so that Φα = I as follows. We first find a product Φ_{ i }of 2cycles that corresponds to an optimal series of the rearrangement events required to transform β_{ i }into J_{ i }(i.e., Φ_{ i }β_{ i }= J_{ i } and hence Φ_{ i }= J_{ i }${\beta}_{j}^{1}$) and then let Φ = Φ_{ ω }Φ_{ω1}... Φ_{1}. Clearly, Φ = Iα^{1} and hence Φ corresponds to a feasible series of events for transforming α into I. Actually, we shall show later that the number of 2cycles in each Φ_{ i }is n_{ i } f(J_{ i }${\beta}_{i}^{1}$), and in Φ_{ i } the number of 2cycles functioning as the fusions and fissions is minimum. This causes that the number of 2cycles in Φ equals to ${\sum}_{i=1}^{\omega}{n}_{i}}f({J}_{i}{\beta}_{i}^{1})=nf(I{\alpha}^{1})$, in which the number of 2cycles serving as the fusions and fissions is minimum. As a result, Φ is an optimal series of events for transforming α into I. The above description indicates that the original problem can be conquered by independently solving the same problem on the smaller instance whose induced bipartite graph is a connected component of $\mathcal{G}$(α, I).
To simplify our discussion, throughout the rest of this section we assume that the induced $\mathcal{G}$(α, I) of a given instance α and I has exactly one connected component. Let Φ = <σ_{1}, σ_{2}, ...., σ_{ δ }> be an optimal series of events for transforming α into I in which all fusions precede all blockinterchanges that further precede all fissions. Let n_{ fu }, n_{ bi }and n_{ fi }denote the numbers of fusions, blockinterchanges and fissions, respectively, in Φ. Then δ = n_{ fu }+ n_{ bi }+ n_{ fi }. In the following, we shall show that Φ can be expressed by a product of n  f(Iα^{1}) 2cycles in which the number of 2cycles functioning as the fusions and fissions is minimum.
It should be noticed that the chromosomes considered here are disjoint (i.e., without gene duplication). Hence, for any two chromosomes α_{ i } and α_{ j }in α with (α_{ i }I_{ k }) ∈ $\mathcal{E}$ and (α_{ j }, I_{ k }) ∈ $\mathcal{E}$, there must exist a fusion in Φ that joins α_{ i } and α_{ j }to one chromosome; otherwise, I_{ k }can not be formed from α by a fission later. Since all needed fusions come together in the beginning of Φ, n_{ fu }= χ(α)  1, which is the lower bound of the number of fusions required in any optimal series of events for transforming α into I. After these n_{ fu }fusions, the resulting α becomes only one chromosome. Since the next n_{ bi }blockinterchanges are intrachromosomal mutations, we have n_{ fi }= χ(I)  1. Actually, χ(I)  1 is the minimum number of the required fissions in any optimal series of events for transforming α into I, since it is the minimum number of the fusions used in the corresponding optimal series of events to transform I into α.
Given any cycle ρ, we use x ∈ ρ to denote that x is a number in ρ. For any two x ∈ ρ and y ∈ ρ, they are said to be adjacent in ρ if ρ(x) = y or ρ(y) = x. Next, we show a way to derive n_{ fu }2cycles from Iα^{1} such that these 2cycles function as the fusions that join all chromosomes of α to a single one, if α has multiple chromosomes, where n_{ fu }= χ(α)  1. For simplicity, later in the text we use "cycle in Iα^{1}" to represent " cycle in the cycle decomposition of Iα^{1}" in meaning, unless a possible confusion may arise.
Lemma 5 Let α_{ i }and α_{ j }be any two disjoint cycles in α. Then there must exist a cycle in Iα^{1} that contains two numbers x and y such that x ∈ α_{ i }and y ∈ α_{ j }.
Proof. Since we assume that the induced $\mathcal{G}$(α, I) contains exactly and only one connected component, and α_{ i }and α_{ j }contain some numbers u and v, respectively, such that both u and v are in a cycle I_{ k }of I. Notice that u ∉ α_{ j }and v ∉ α_{ i }. Suppose that there is no cycle in Iα^{1} that contains two numbers x and y such that they are in these two different cycles of α, say x ∈ α_{ i } and y ∈ α_{ j }. Then all numbers in any cycle of Iα^{1} are contained in some cycle of α. Without loss of generality, let I_{ k }= (u ≡ a_{1}, a_{2}, ..., a_{ p }≡ v, ..., a_{ q }) and let p <q for simplifying the discussion. For each 1 ≤ x ≤ p, let α(a_{ x }) = b_{ x }. Then we have Iα^{1}(b_{ x }) = a_{x+1}(since Iα^{1}α = I), which means that both b_{ x }and a_{x+1}are in the same cycle of Iα^{1} and hence they are also in the same cycle of α. If a_{ x }is in α_{ i }, then b_{ x }is also in α_{ i }, which further leads to a_{x+1}∈ α_{ i }. Since u = a_{1} is in α_{ i }, all a_{2}, a_{3}, ..., a_{ p }are in α_{ i }. As a result, both of u and v are in α_{ i }, a contradiction. Hence, there exists a cycle in Iα^{1} that contains x and y such that x ∈ α_{ i }and y ∈ α_{ j }.
The following lemma can be easily verified.
Lemma 6 (a_{1}, a_{2},..., a_{ i },..., a_{ j }) = (a_{1}, a_{2},..., a_{ i }) (a_{i+1},..., a_{ j }) (a_{ i }, a_{ j }), where 1 ≤ i <j.
According to Lemma 5, for any two cycles α_{ i } and α_{ j }of α, we can find two numbers x and y in a cycle of Iα^{1}, say β, such that x ∈ α_{ i } and y ∈ α_{ j }. Let β = (a_{1}, a_{2},..., a_{ q }), where q ≥ 2. Then we consider the following two cases. Case 1: x and y are adjacent in β. For simplicity, let x = a_{q1}and y = a_{ q }. Then by Lemma 6, β = (a_{1}, a_{2},..., a_{q1}) (a_{ q }) (x, y). Case 2: x and y are not adjacent in β. Let x = a_{ p }and y = a_{ q }, where 1 ≤ p <q  1. Then β = (a_{1}, a_{2},..., a_{ p }) (a_{p+1},..., a_{ q }) (x,y) according to Lemma 6. In other words, we can derive a 2cycle (x, y) from β such that it can join α_{ i } and α_{ j }to one cycle. After α_{ i }, and α_{ j }are joined together via (x, y), the number of the cycles (including 1cycles) in the resulting Iα^{1} increases by one. Repeatedly based on the procedure above, we can derive consecutive n_{ fu }2cycles from Iα^{1}, say φ_{1}, φ_{2},...,${\phi}_{{n}_{fu}}$, that can join χ (α) cycles in α to a single one, where n_{ fu }= χ(α)  1. In other words, φ_{1}, φ_{2},...,${\phi}_{{n}_{fu}}$ function as χ(α)  1 fusions that transform genome α with χ(α) chromosomes into a genome, denoted by α', with a single chromosome. Clearly, we have α' = ${\phi}_{{n}_{fu}}$${\phi}_{{n}_{fu}1}$...φ_{1}α, Iα'^{1} = Iα^{1} φ_{1}φ_{2} ... ${\phi}_{{n}_{fu}}$, and f(Iα'^{1}) = f(Iα^{1}) + n_{ fu }. Hence, we can immediately claim the following.
Claim 1 α' = ${\phi}_{{n}_{fu}}$${\phi}_{{n}_{fu}1}$...φ_{1}α, Iα'^{1} = Iα^{1} φ_{1}φ_{2} ... ${\phi}_{{n}_{fu}}$, and f(Iα'^{1}) = f(Iα^{1}) + n_{ fu }, where n_{ fu }= χ(α)  1.
Without loss of generality, we now suppose that χ(I) > 1. Similarly as the discussion above, we can derive consecutive n_{ fi }2cycles from α'I^{1}, say ψ_{1}, ψ_{2}, ..., ${\psi}_{{n}_{fi}}$, such that they serve as the fusions to transform I with χ(I) chromosomes into a genome, denoted by I', with only one chromosome, where n_{ fi }= χ(I)  1 and α'I^{1} is the inverse of Iα'^{1} (i.e., α'I^{1} = (Iα'^{1})^{1}). Then we have I' = ${\psi}_{{n}_{fi}}$${\psi}_{{n}_{fi}1}$ ... ψ_{1}I (hence ψ_{1}ψ_{2} ... ${\psi}_{{n}_{fi}}$I' = I), α'I'^{1} = α'I^{1}ψ_{1}ψ_{2} ... ${\psi}_{{n}_{fi}}$ and f(α'I'^{1}) = f(α'I^{1}) + n_{ fi }. Conversely, we can use ${\psi}_{{n}_{fi}}$, ${\psi}_{{n}_{fi}1}$…, ψ_{1} as fissions to split I' with one chromosome into I with n_{ fi }chromosomes. Since α'I^{1} is the inverse of Iα'^{1}, it can be easily obtained from Iα'^{1} by just reversing the order of the numbers in each cycle of Iα'^{1} and hence f(α'I^{1}) = f(Iα'^{1}), which leads to f(α'I'^{1}) = f(Iα'^{1}) + n_{ fi }. As a result, we have f(I'α'^{1}) = f(α'I^{1}) = f(I α'^{1}) + n_{ fi }since I'α'^{1} = (α'I^{1})^{1}. Therefore, the following claim can be obtained.
Claim 2 ψ_{1}ψ_{2} ... ${\psi}_{{n}_{fi}}$I' = I, α'I'^{1} = α'I^{1}ψ_{1}ψ_{2} ... ${\psi}_{{n}_{fi}}$and f(I'α'^{1}) = f(Iα'^{1}) + n_{ fi }, where n_{ fi }= χ(I)  1.
Notice that both α' and I' now are the genomes with only one chromosome. Then based on the algorithm proposed by Lin et al. [14], we can find ${n}_{bi}=\frac{nf({I}^{\prime}{{\alpha}^{\prime}}^{1})}{2}$ blockinterchanges from I'α'^{1} to transform α' into I'. Certainly, these n_{ bi }blockinterchanges can be further expressed by a product of 2n_{ bi }2cycles, say ${\tau}_{1}^{1},{\tau}_{1}^{2},{\tau}_{2}^{1},{\tau}_{2}^{2},\dots ,{\tau}_{{n}_{bi}}^{1},{\tau}_{{n}_{bi}}^{2}$, such that every two consecutive 2cycles act as a blockinterchange in the process of transforming α' into I', where I'α'^{1} = ${\tau}_{{n}_{bi}}^{2}{\tau}_{{n}_{bi}}^{1}{\tau}_{{n}_{bi}1}^{2}{\tau}_{{n}_{bi}1}^{1}\dots {\tau}_{1}^{2}{\tau}_{1}^{1}$. Hence, we have the following claim immediately.
Claim 3 I' = ${\tau}_{{n}_{bi}}^{2}{\tau}_{{n}_{bi}}^{1}{\tau}_{{n}_{bi}1}^{2}{\tau}_{{n}_{bi}1}^{1}\dots {\tau}_{1}^{2}{\tau}_{1}^{1}$α'.
Now we let Φ = ψ_{1}ψ_{2}...${\psi}_{{n}_{fi}}{\tau}_{{n}_{bi}}^{2}{\tau}_{{n}_{bi}}^{1}\dots {\tau}_{1}^{2}{\tau}_{1}^{1}{\phi}_{{n}_{fu}}{\phi}_{{n}_{fu}1}\dots {\phi}_{1}$. Then the result of Φα = I (hence Φ = Iα^{1}) can be easily verified by Claims 1, 2 and 3 as follows.
In other words, Φ is a product of (n_{ fu }+ (n  f(I'α'^{1})) + n_{ fi }) 2cycles that can transform α into I. More clearly, Φ first uses φ_{1}, φ_{2},...,${\phi}_{{n}_{fu}}$ (acting as n_{ fu }fusions) to transform α into α', then uses ${\tau}_{1}^{1},{\tau}_{1}^{2},\dots ,{\tau}_{{n}_{bi}}^{1},{\tau}_{{n}_{bi}}^{2}$ (acting as n_{ bi }blockinterchanges) to transform α' into I', and finally uses ${\psi}_{{n}_{fi}}$, ${\psi}_{{n}_{fi}1}$,..., ψ_{1} (acting as n_{ fi }fissions) to transform I' into I. By Claims 1 and 2, we can show that (n_{ fu }+ (n  f(I'α'^{1})) + n_{ fi }= n  f(Iα^{1}) as follows.
n_{ fu }+ (n  f(I'α'^{1})) + n_{ fi }
= n_{ fu }+ (n  f(Iα'^{1}) + n_{ fi }) + n_{ fi } (by Claim 2)
= n_{ fu }+ (n  f(Iα^{1}) + n_{ fu } n_{ fi }) + n_{ fi } (by Claim 1)
= n  f(Iα^{1})
As mentioned before, χ(α)  1 and χ(I)  1 are the lower bounds of the numbers of fusions and fissions, respectively, required in any optimal series of rearrangement events for transforming α into I. Hence, the number of 2cycles in Φ that function as the fusions and fissions is minimum. Along with that Φ = Iα^{1} can be expressed as a product of n  f(Iα^{1}) 2cycles, we thus conclude that Φ is an optimal series of the events that transform α into I with first n_{ fu }fusions, then n_{ bi }blockinterchanges and finally n_{ fi }fissions, where n_{ fu }= χ(α) 1, ${n}_{bi}=\frac{nf({I}^{\prime}{{\alpha}^{\prime}}^{1})}{2}=\frac{nf(I{\alpha}^{1}){n}_{fu}{n}_{fi}}{2}$, and n_{ fi }= χ(I)  1.
Lemma 7 There is an optimal series of the events needed to transform α into I in a canonical order such that all n_{ fu }fusions come before all n_{ bi }blockinterchanges that come before all n_{ fi }fissions, where n_{ fu }= χ(α)  1, ${n}_{bi}=\frac{nf(I{\alpha}^{1}){n}_{fu}{n}_{fi}}{2}$and n_{ fi }= χ(I)  1.
Let us take α = (1, 2,10) (11, 8, 9, 3, 6) (7, 4, 5, 12) and I = (1, 2, 3) (4, 5) (6, 7, 8) (9, 10, 11, 12) for an example. It should be straightforward to see that $\mathcal{G}$(α,I) is a connected bipartite graph with χ(α) = 3 and χ(I) = 4, and Iα^{1} = (1, 11, 7, 9, 6) (3, 10) (4, 8, 12) and hence f(Iα^{2}) = 5, since two 1cycles (i.e., (2) and (5)) are not explicitly shown. First, we are to find two 2cycles φ_{1} and φ_{2} (since n_{ fu }= χ(α) 1 = 2) from Iα^{1} to transform genome α with three chromosomes into genome α' with exactly one chromosome. To this purpose, we let φ_{1} = (3,10) and φ_{2} = (4,8), since Iα^{1} = (1, 11, 7, 9, 6) (4, 12) (4, 8) (3, 10). Then by Claim 1, α' = φ_{2}φ_{1}α = (4, 5, 12, 7, 8, 9, 10, 1, 2, 3, 6, 11) and Iα'^{1} = Iα^{1}φ_{1}φ_{2} = (1, 11, 7, 9, 6) (4, 12). Next, we need to find three 2cycles ψ_{1}, ψ_{2} and ψ_{3} (since n_{ fi }= χ(I)  1 = 3) from α'I^{1}, which is equal to (Iα'^{1})^{1} = (6, 9, 7, 11, 1) (12, 4) = (1, 7, 11) (1, 9) (1, 6) (12, 4), to transform I into I' with only one chromosome. By letting ψ_{1} = (12, 4), ψ_{2} = (1, 6) and ψ_{3} = (1, 9), we have I' = ψ_{3}ψ_{2}ψ_{1} = (1, 2, 3, 6, 7, 8, 9, 10, 11, 4, 5, 12) and α'I'^{1} = α'I^{1}ψ_{1}ψ_{2}ψ_{3} = (1, 7, 11) according to Claim 2. Finally, we will find two 2cycles ${\tau}_{\text{1}}^{\text{1}}$ and ${\tau}_{\text{1}}^{\text{2}}$ (since n  f(Iα^{1})  n_{ fu } n_{ fi }= 12  5  2  3 = 2) from I'α'^{1} that act as a blockinterchange to transform α' into I', where I'α'^{1} = (α'I'^{1})^{1} = (11, 7, 1) = (11, 1) (11, 7). By letting ${\tau}_{\text{1}}^{\text{1}}$ = (11, 7) and ${\tau}_{\text{1}}^{\text{2}}$ = (11, 1), we have ${\tau}_{\text{1}}^{\text{2}}$${\tau}_{\text{1}}^{\text{1}}$α' = (11, 1) (11, 7) (4, 5, 12, 7, 8, 9, 10, 1, 2, 3, 6, 11) = (11, 4, 5, 12, 1, 2, 3, 6, 7, 8, 9, 10), which indeed equals I'. Consequently, we find an optimal series of events Φ= ψ_{1}ψ_{2}ψ_{3}${\tau}_{\text{1}}^{\text{2}}$${\tau}_{\text{1}}^{\text{1}}$φ_{2}φ_{1} that transform α into I (i.e., Φα = I).
Based on the idea above, we have designed Algorithm Sortingbyffbi (meaning sorting by fusions, fissions and blockinterchanges) to compute the genome rearrangement distance d(α, I) between two given circular multichromosomal genomes α and I, and also to generate an optimal scenario of the required rearrangement events in a canonical order. In Algorithm Sortingbyffbi, the purpose of Step 2.3.3 (respectively, Step 2.4.4) is to find two numbers x and y that are both in some cycle of γ = J_{ i }${\beta}_{i}^{1}$ (respectively, γ = ${{\beta}^{\prime}}_{i}{J}_{i}^{1}$), but in different cycles in β_{ i }. By Lemma 5, such x and y exist. In fact, they can be found using the following simple approach. For simplicity, let γ_{ k }= (${a}_{k}^{1},{a}_{k}^{2},\dots ,{a}_{k}^{{l}_{k}}$) be a cycle in γ that contains two numbers x and y such that they are in different cycles in β_{ i }. Then we only need to check whether ${a}_{k}^{1}$ and ${a}_{k}^{j}$, where 2 ≤ j ≤ l_{ k }, are in different cycles in β_{ i }or not. If so, we let x = ${a}_{k}^{1}$ and y = ${a}_{k}^{j}$. The reason is as follows. Suppose that both ${a}_{k}^{1}$ and ${a}_{k}^{j}$ for all 2 ≤ j ≤ l_{ k }are in the same cycle in β_{ i }. Then all of numbers ${a}_{k}^{1},{a}_{k}^{2},\dots ,{a}_{k}^{{l}_{k}}$ in γ_{ k }are in the same cycle in β_{ i }, which contradicts the above assumption that γ_{ k }contains x and y that are in different cycles in β_{ i }.
Algorithm sortingbyffbi
Input: Two circular multichromosomal genomes α and I.
Output: d(α,I) and a minimum series Φ of events required to transform α into I.
1: Find all connected components $\mathcal{C}$_{ 1 }, $\mathcal{C}$_{ 2 },...,$\mathcal{C}$_{ ω }in graph $\mathcal{G}$(α,I);
/* Denote by n_{ i } the number of genes in $\mathcal{C}$_{ i }. */
2: for each $\mathcal{C}$_{ i }, 1 ≤ i ≤ ω,do
/* Denote by β_{ i }(resp. J_{ i }) the collection of chromosomes in α (resp. I) whose corresponding vertices are in $\mathcal{C}$_{ i }. */
2.1: Compute J_{ i }${\beta}_{i}^{1}$ and let γ = J_{ i }${\beta}_{i}^{1}$;
2.2: n_{ fu }= χ(β_{ i })  1, n_{ fi }= χ(J_{ i }) 1, ${n}_{bi}=\frac{{n}_{i}f(\gamma ){n}_{fu}{n}_{fi}}{2}$ and δ_{ i }= n_{ fu }+ n_{ bi }+ n_{ fi };
2.3: if χ(β_{ i }) > 1 then /* To compute φ_{1}, φ_{2},...,${\phi}_{{n}_{fu}}$ */
2.3.1: for each cycle of β_{ i }do
Create a set to contain all the numbers in this cycle;
endfor
2.3.2: /* Let γ_{1}γ_{2} ... γ_{ p }be the cycle decomposition of the current γ and
let γ_{ q }= (${a}_{q}^{1},{a}_{q}^{2},\dots ,{a}_{q}^{{l}_{q}}$), where 1 ≤ q ≤ p and l_{ q }≥ 2 */
k = 1 and h = 2;
2.3.3: for j = 1 to n_{ fu }do
S = findset (${a}_{k}^{1}$);
while S = findset(${a}_{k}^{h}$)) do
if h <l_{ k }then h = h + 1; else k = k + 1, h = 2 and S = findset (${a}_{k}^{1}$);
endwhile
x = ${a}_{k}^{1}$ and y = ${a}_{k}^{h}$;
φ_{ j }= (x, y) and union(x, y);
endfor
2.3.4: ${{\beta}^{\prime}}_{i}$ = ${\phi}_{{n}_{fu}}$${\phi}_{{n}_{fu}1}$...φ_{1}β_{ i }and γ = γ φ_{1}φ_{2}...${\phi}_{{n}_{fu}}$; /* Currently, γ is J_{i}${{\beta}^{\prime}}_{i}^{1}$ */
endif
2.4: if χ(J_{ i }) > 1 then/* To compute ψ_{1}, ψ_{2},...,${\psi}_{{n}_{fi}}$ */
2.4.1: γ = γ^{1}/* New γ becomes ${{\beta}^{\prime}}_{i}{J}_{i}^{1}$ */
2.4.2: for each cycle of J_{ i } do
Create a set to contain all the numbers in this cycle;
endfor
2.4.3: /* Let γ_{1}γ_{2}...γ_{ p }be the cycle decomposition of the current γ and
let γ_{ q }= (${a}_{q}^{1},{a}_{q}^{2},\dots ,{a}_{q}^{{l}_{q}}$), where 1 ≤ q ≤ p and l_{ q }≥ 2 */
k = 1 and h = 2;
2.4.4: for j = 1 to n_{ fi }do
S = findset(${a}_{k}^{1}$);
while (S = findset(${a}_{k}^{h}$)) do
if h <I_{ k }then h = h + 1; else k = k + 1,h = 2 and S = findset(${a}_{k}^{1}$);
endwhile
x = ${a}_{k}^{1}$ and y = ${a}_{k}^{h}$;
ψ_{ j }= (x, y) and union(x, y);
endfor
2.4.5: = ${\psi}_{{n}_{fi}}$${\psi}_{{n}_{fi}1}$... ψ_{1}J_{ i }and γ = γ ψ_{1}ψ_{2}...${\psi}_{{n}_{fi}}$; /* Currently, γ is ${{\beta}^{\prime}}_{i}{{J}^{\prime}}_{i}^{1}$ */
endif
2.5: /* To compute ${\tau}_{1}^{1},{\tau}_{1}^{2},\dots ,{\tau}_{{n}_{bi}}^{1},{\tau}_{{n}_{bi}}^{2}$ */
2.5.1: γ = γ^{1};/* New γ becomes ${{J}^{\prime}}_{i}$${{\beta}^{\prime}}_{i}^{1}$ */
2.5.2: ${n}_{bi}=\frac{{n}_{i}f(\gamma )}{2}$;
2.5.3: for j = 1 to n_{ bi }do
Arbitrarily choose two adjacent elements x and y in γ;
/* Let ${{\beta}^{\prime}}_{i}$ = (a_{1},a_{2}...,${a}_{{n}_{i}}$) */
Circularly shift (a_{1}, a_{2},...,${a}_{{n}_{i}}$) such that a_{1} = x and assume y = a_{ k };
${\tau}_{j}^{\text{1}}$ = (x, y);
for h = 1 to n_{ i }do
index (a_{ h }) = h;
end for
Find two adjacent elements u and v in γ (x, y) such that
index(u) ≤ k  1 and index(v) ≥ k;
${\tau}_{j}^{\text{2}}$ = (u, v);
${{\beta}^{\prime}}_{i}$ = ${\tau}_{j}^{\text{2}}$${\tau}_{j}^{\text{1}}$${{\beta}^{\prime}}_{i}$ and γ = γ${\tau}_{j}^{\text{1}}$${\tau}_{j}^{\text{2}}$;
endfor
3: Let ${\Phi}_{i}={\psi}_{1}\dots {\psi}_{{n}_{fi}}{\tau}_{{n}_{bi}}^{2}{\tau}_{{n}_{bi}}^{1}\dots {\tau}_{1}^{2}{\tau}_{1}^{1}{\phi}_{{n}_{fu}}\dots {\phi}_{1}$ for each 1 ≤ i ≤ ω;.
4: Output d(α, I) = $\sum}_{i=1}^{\omega}{\delta}_{i$ and Φ = Φ_{1}Φ_{2}...Φ_{ ω };
Theorem 1 Given two circular multichromosomal genomes α and I over the same gene set E = {1,2,..., n}, the problem of computing the genome rearrangement distance between α and I using fusions, fissions and blockinterchanges can be solved and an optimal series of such events in a canonical order can be obtained in $\mathcal{O}$(n^{2}) time.
Proof. As discussed above, Algorithm Sortingbyffbi transforms α into I using the minimum number of fusions, fissions and blockinterchanges. Next, we follow to analyze its timecomplexity. Notice that given an undirected graph with p vertices and q edges, all the connected components in this graph can be found in $\mathcal{O}$(p + q) time using depthfirst search or breadthfirst search [31]. As a result, Step 1 can be done in $\mathcal{O}$(n^{2}) time for computing the connected components in the induced bipartite graph $\mathcal{G}$(α, I), since in the worst case, the number of edges in $\mathcal{G}$(α,I) is χ(α) × χ(I) and χ(α) = $\mathcal{O}$(n) and χ(I) = $\mathcal{O}$(n). In Step 2, there are ω outer iterations, each computing the minimum series of events needed to transform β_{ i }into J_{ i }, where 1 ≤ i ≤ ω. Clearly, Steps 2.1 costs $\mathcal{O}$(n_{ i }) time for computing J_{ i }${\beta}_{i}^{1}$ and Step 2.2 takes only a constant time. The time cost of Step 2.3 is mostly contributed from that of Step 2.3.3. There are n_{ fu }inner iterations in Step 2.3, each with the purpose of finding two numbers x and y that are both in the same cycle in γ and, however, in the different cycles in β_{ i }. In the worst case, Step 2.3 needs n_{ i } findset operations and n_{ fu }union operations to finish its overall process. Note that Step 2.3.1 can be implemented by initially creating a set for each number in gene(β_{ i }) and then performing n_{ i } χ(β_{ i }) union operations to generate χ(β_{ i }) sets with each corresponding to a cycle in β_{ i }, where χ(β_{ i }) = n_{ fu }+ 1. Hence, the total number of union operations is n_{ i } 1 in Step 2.3. In fact, these findset and union operations can be implemented in $\mathcal{O}$(n_{ i }) time using the socalled "static disjoint set union and find" algorithm proposed by Gabow and Tarjan [32]. In other words, Step 2.3 can cost only $\mathcal{O}$(n_{ i }) time. By the same principle, it can be verified that the time cost of Step 2.4 is $\mathcal{O}$(n_{ i }). As for Step 2.5, adopted from our previous work [14], it takes $\mathcal{O}$(n_{ bi }n_{ i }) time, where = $\mathcal{O}$n_{ i }). As a result, the time cost of Step 2 is $\mathcal{O}$(ζ n), where ζ is the maximum n_{ bi }among all iterations in Step 2 and ζ <n. Clearly, Steps 3 and 4 cost constant time. Therefore, the total timecomplexity of Algorithm Sortingbyffbi is $\mathcal{O}$(n^{2}).
Construction of orthologous genes
To analyze the rearrangement of three Vibrio genomes, we identified and constructed a table of orthologous genes that are putatively not involved in horizontal gene transfer (HGT) events by adopting the socalled symmetrical best hits (SymBets for short). In principle, two genes match and give SymBets if they are more similar to each other than they are to any other genes from the compared genomes [33, 34]. Detection of such SymBet genes is, yet arguably, the simplest and most suitable method for identification of probable orthologs for closely related genomes [33, 34]. Particularly, this prediction of orthologs holds to be applicable even when sequence similarity between the compared proteins is relatively low [33, 34]. For the purpose of excluding paralogous genes derived from lineagespecific gene duplications, we here considered only onetoone orthologous genes, which actually have been demonstrated as a major pattern in prokaryotic genome evolution [34]. Therefore, we used the following steps to identify an HGTfree table of onetoone orthologous genes from the three complete Vibrio genomes. First, the GenePlot [35] program offered by NCBI was utilized to find and construct a table of SymBet genes between each pair of Vibrio genomes. Next, after removing all onetomany or manytomany SymBets, the three resulting tables of SymBet genes were joined to give a new one of onetoone orthologous genes for all the three Vibrio genomes by using a rule as follows. If genes a (from genome A) and b (from genome B),b and c (from genome C), and c and a are all onetoone SymBet pairs, then a, b and c are considered as onetoone orthologous genes for the genomes A,B and C. In other words, the SymBet relationships among a, b and c result in a triangle. Finally, those genes that were involved in the putative HGT events detected and available in the Horizontal Gene Transfer Database [36] were then deleted from the table of onetoone orthologous genes.
Note
^{1}The Institute for Genome Research (TIGR) offers a website of Comprehensive Microbial Resource (CMR) that provides information on all of the publicly available and complete bacterial genomes.
References
 1.
Hannenhalli S, Pevzner PA: Transforming cabbage into turnip: polynomial algorithm for sorting signed permutations by reversals. J Assoc Comput Mach 1999, 46: 1–27.
 2.
Bader DA, Yan M, Moret BMW: A lineartime algorithm for computing inversion distance between signed permutations with an experimental study. J Comput Biol 2001, 8: 483–491.
 3.
Bafna V, Pevzner PA: Genome rearrangements and sorting by reversals. SIAM J Comput 1996, 25: 272–289.
 4.
Berman P, Hannenhalli S: Fast sorting by reversal. In Proceedings of the 7th Annual Symposium on Combinatorial Pattern Matching (CPM1996), Volume 1075 of Lecture Notes in Computer Science. Edited by: Hirschberg DS, Myers E. SpringerVerlag; 1996:168–185.
 5.
Berman P, Hannenhalli S, Karpinski M: 1.375approximation algorithm for sorting by reversals. In Proceedings of the 10th Annual European Symposium on Algorithms (ESA2002), Volume 2461 of Lecture Notes in Computer Science. Edited by: Mohring RH, Raman R. SpringerVerlag; 2002:200–210.
 6.
Caprara A: Sorting by reversal is difficult. In Proceedings of the 1th Annual International Conference on Research in Computational Molecular Biology (RECOMB1997). ACM Press; 1997:75–83.
 7.
Caprara A: Sorting permutations by reversals and Eulerian cycle decompositions. SIAM J Discrete Math 1999, 12: 91–110.
 8.
Christie DA: A 3/2approximation algorithm for sorting by reversals. In Proceedings of the 9th Annual ACMSIAM Symposium on Discrete Algorithms (SODA1998), ACM/SIAM 1998, 244–252.
 9.
Kaplan H, Shamir R, Tarjan RE: Faster and simpler algorithm for sorting signed permutations by reversals. SIAM J Comput 2000, 29: 880–892.
 10.
Kececioglu JD, Sankoff D: Exact and approximation algorithms for the inversion distance between two permutations. In Proceedings of the 4th Annual Symposium on Combinatorial Pattern Matching (CPM1993), Volume 684 of Lecture Notes on Computer Science. Edited by: Apostolico A, Crochemore M, Galil Z, Manber U. SpringerVerlag; 1993:87–105.
 11.
Bafna V, Pevzner PA: Sorting by transpositions. SIAM J Discrete Math 1998, 11: 221–240.
 12.
Walter MEMT, Dias Z, Meidanis J: Reversal and transposition distance of linear chromosomes. In Proceedings of String Processing and Information Retrieval (SPIRE1998). IEEE Computer Society; 1998:96–102.
 13.
Christie DA: Sorting by blockinterchanges. Inform Process Lett 1996, 60: 165–169.
 14.
Lin YC, Lu CL, Chang HY, Tang CY: An efficient algorithm for sorting by block interchanges and its application to the evolution of Vibrio species. J Comput Biol 2005, 12: 102–112.
 15.
Lu CL, Wang TC, Lin YC, Tang CY: ROBIN: a tool for genome rearrangement of blockinterchanges. Bioinformatics 2005, 21: 2780–2782.
 16.
Hannenhalli S: Polynomial algorithm for computing translocation distance between genomes. Discrete Appl Math 1996, 71: 137–151.
 17.
Kececioglu JD, Ravi R: Of mice and men: algorithms for evolutionary distances between genomes with translocation. In Proceedings of the 6th ACMSIAM Symposium on Discrete Algorithms (SODA1995). ACM/SIAM, San Francisco; 1995:604–613.
 18.
Hannenhalli S, Pevzner PA: Transforming men into mice (polynomial algorithm for genomic distance problem). In Proceedings of the 36th IEEE Symposium on Foundations of Computer Science (FOCS1995). IEEE Computer Society; 1995:581–592.
 19.
Meidanis J, Bias Z: Genome rearrangements distance by fusion, fission, and transposition is easy. In Proceedings of the 8th International Symposium on String Processing and Information Retrieval (SPIRE2001). Edited by: Navarro G. IEEE Computer Society; 2001:250–253.
 20.
Tillier ER, Collins RA: Genome rearrangement by replicationdirected translocation. Nat Genet 2000, 26: 195–197.
 21.
Yancopoulos S, Attie O, Friedberg R: Efficient sorting of genomic permutations by translocation, inversion and block interchange. Bioinformatics 2005, 21: 3340–3346.
 22.
FFBI: a tool of circular genome rearrangement by fusions, fissions and blockinterchanges[http://genome.life.nctu.edu.tw/FFBI/]
 23.
The COG database[http://www.ncbi.nlm.nih.gov/COG/]
 24.
Dorsch M, Lane D, Stackebrandt E: Towards a phylogeny of the genus Vibrio based on 16S rRNA sequences. Int J Syst Bacterial 1992, 42: 58–63.
 25.
KitaTsukamoto K, Oyaizu H, Nanba K, Simidu U: Phylogenetic relationships of marine bacteria, mainly members of the family Vibrionaceae ,determined on the basis of 16S rRNA sequences. Int J Syst Bacteriol 1993, 43: 8–19.
 26.
Okada K, Iida T, KitaTsukamoto K, Honda T: Vibrios commonly possess two chromosomes. J Bacteriol 2005, 187: 752–757.
 27.
Thompson JD, Higgins DG, Gibson TJ: CLUSTAL W: improving the sensitivity of progressive multiple sequence alignment through sequence weighting, positionspecific gap penalties and weight matrix choice. Nucleic Acids Res 1994, 22: 4673–4680.
 28.
Felsenstein J: PHYLIP: phylogeny inference package (version 3.2). Cladistics 1989, 5: 164–166.
 29.
Fraleigh JB: A First Course in Abstract Algebra. 7th edition. AddisonWesley; 2003.
 30.
Meidanis J, Dias Z: An alternative algebraic formalism for genome rearrangements. In Comparative Genomics: Empirical and Analytical Approaches to Gene Order Dynamics, Map Alignment and Evolution of Gene Families. Edited by: Sankoff D, Nadeau JH. Kluwer Academic Publisher; 2000:213–223.
 31.
Cormen TH, Leiserson CE, Rivest RL, Stein C: Introduction to Algorithms. 2nd edition. The MIT Press; 2001.
 32.
Gabow HN, Tarjan RE: A lineartime algorithm for a special case of disjoint set union. J Comput Syst Sci 1985, 30: 209–221.
 33.
Tatusov RL, Koonin EV, Lipman DJ: A genomic perspective on protein families. Science 1997, 278: 631–637.
 34.
Koonin EV: Orthologs, paralogs, and evolutionary genomics. Annu Rev Genet 2005, 39: 309–338.
 35.
 36.
GarciaVallve S, Guzman E, Montero MA, Romeu A: HGTDB: a database of putative horizontally transferred genes in prokaryotic complete genomes. Nucleic Acids Res 2003, 31: 187–189. [http://www.fut.es/~debb/HGT/]
Acknowledgements
The authors would like to thank anonymous referees for many constructive comments in the presentation of this paper. This work was supported in part by National Science Council of Republic of China under grants NSC 942213E009141 (C.L. Lu) and NSC 942113M009007 (H.T. Chiu).
Author information
Additional information
Authors' contributions
CLL and HTC contributed equally to this work. CLL conceived of the study, participated in the design and analysis of algorithm and drafted the manuscript. YLH participated in the algorithm design and software development, and carried out the bioinformatics experiments. TCW participated in the software development. HTC participated in the design and coordination of this study as well as in drafting the manuscript.
Authors’ original submitted files for images
Below are the links to the authors’ original submitted files for images.
Rights and permissions
About this article
Received
Accepted
Published
DOI
Keywords
 Orthologous Gene
 Genome Rearrangement
 Optimal Scenario
 Circular Chromosome
 Disjoint Cycle