Optimal Phylogenetic Networks with Constrained and Unconstrained Recombination

Optimal Phylogenetic Networks with Constrained and Unconstrained Recombination (The root-unknown case) Dan Gusfield UC Davis

Reconstructing the Evolution of Binary Bio-Sequences • Perfect Phylogeny (tree) model • Phylogenetic Networks (DAG) with recombination • Phylogenetic Networks with disjoint cycles: Galled-Trees • Phylogenetic Networks with unconstrained cycles: Blobbed-Trees • Combinatorial Structure and Efficient Algorithms

Geneological or Phylogenetic Networks • The major biological motivation comes from genetics and attempts to reconstruct the history of recombination in populations. • Also relates to phylogenetic-based haplotyping. • The algorithmic and mathematical results also have phylogenetic applications, for example in lateral gene transfer.

The Perfect Phylogeny Model for binary sequences sites 12345 Ancestral sequence 00000 1 4 Site mutations on edges 3 The tree derives the set M: 2 10100 5 10000 01011 01010 00010 10000 00010 01011 Extant sequences at the leaves

When can a set of sequences be derived on a perfect phylogeny with the all-0 root? Classic NASC: Arrange the sequences in a matrix. Then (with no duplicate columns), the sequences can be generated on a unique perfect phylogeny if and only if no two columns (sites) contain all four pairs: 0, 0 and 0, 1 and 1, 0 and 1, 1 This is the 4 -Gamete Test

A richer model 10100 10000 01011 01010 00010 10101 added pair 4, 5 fails the three gamete-test. The sites 4, 5 ``conflict”. 12345 00000 1 4 3 10100 10000 2 00010 5 0101101010 Real sequence histories often involve recombination.

Sequence Recombination 01011 10100 S P 5 Single crossover recombination 10101 A recombination of P and S at recombination point 5. The first 4 sites come from P (Prefix) and the sites from 5 onward come from S (Suffix).

Network with Recombination 10100 10000 01011 01010 00010 10101 new 12345 00000 1 4 3 2 10100 The previous tree with one recombination event now derives all the sequences. P 10000 5 10101 S 00010 5 0101101010

Multiple Crossover Recombination 4 -crossovers 2 -crossovers = ``gene conversion”

Elements of a Phylogenetic Network (single crossover recombination) • Directed acyclic graph. • Integers from 1 to m written on the edges. Each integer written only once. These represent mutations. • A choice of ancestral sequence at the root. • Every non-root node is labeled by a sequence obtained from its parent(s) and any edge label on the edge into it. • A node with two edges into it is a ``recombination node”, with a recombination point r. One parent is P and one is S. • The network derives the sequences that label the leaves.

A Phylogenetic Network 00010 a: 00010 10010 00000 4 3 1 00100 2 b: 10010 P c: 00100 S 3 01100 p d: 10100 5 00101 S 4 01101 f: 01101 e: 01100 g: 00101

Which Phylogenetic Networks are meaningful? Given M we want a phylogenetic network that derives M, but which one? A: A perfect phylogeny (tree) if possible. As little deviation from a tree, if a tree is not possible. Use as little recombination or gene-conversion as possible.

Minimizing recombinations • Any set M of sequences can be generated by a phylogenetic network with enough recombinations, and one mutation per site. This is not interesting or useful. • However, the number of (observable) recombinations is small in realistic sets of sequences. ``Observable” depends on n and m relative to the number of recombinations. • Two algorithmic problems: given a set of sequences M, find a phylogenetic network generating M, minimizing the number of recombinations (Hein’s problem). Find a network generating M that has some biologicallymotivated structural properties.

Minimization is NP-hard The problem of finding a phylogenetic network that creates a given set of sequences M, and minimizes the number of recombinations, is NPhard. (Wang et al 2000) They explored the problem of finding a phylogenetic network where the recombination cycles are required to be node disjoint, if possible. They gave a sufficient but not a necessary condition to recognize cases when this is possible. O(nm + n^4) time.

Recombination Cycles • In a Phylogenetic Network, with a recombination node x, if we trace two paths backwards from x, then the paths will eventually meet. • The cycle specified by those two paths is called a ``recombination cycle”.

Galled-Trees A recombination cycle in a phylogenetic network is called a “gall” if it shares no node with any other recombination cycle. A phylogenetic network is called a “galledtree” if every recombination cycle is a gall.

A galled-tree generating the sequences generated by the prior network. 4 3 1 p s a: 00010 3 b: 10010 c: 00100 d: 10100 2 p e: 01100 5 4 s f: 01101 g: 00101

Sales pitch for Galled-Trees Galled-trees represent a small deviation from true trees. There are sufficient applications where it is plausible that a galled tree exists that generates the sequences. Observable recombinations tend to be recent; block structure of human DNA; recombination is sparse, so the true history of observable recombinations may be a galled-tree. The number of recombinations is never more than m/2. Moreover, when M can be derived on a galled-tree, the number of recombinations used is the minimum possible over any phylogenetic network, even if multiple cross-overs at a recombination event are counted as a single recombination. A galled-tree for M is ``almost unique” - implications for reconstructing the correct history.

Old (Aug. 2003) Results • O(nm + n^3)-time algorithm to determine whether or not M can be derived on a galled-tree with all-0 ancestral sequence. • Proof that the galled-tree produced by the algorithm is a “nearly-unique” solution. • Proof that the galled-tree (if one exists) produced by the algorithm minimizes the number of recombinations used, over all phylogeneticnetworks with all-0 ancestral sequence.

The initial version is in the Proceedings of the 2003 IEEE CSB conference. The expanded journal version with proof of optimality is in “Optimal, Efficient Reconstruction of Phylogenetic Networks with Constrained Recombination” by Gusfield, Eddhu and Langley. J. Bioinformatics and Computational Biology, Vol. 2 No. 1 (2004) p. 173 -213

New work We derive the galled-tree results in a more general setting that addresses unconstrained recombination cycles and multiple crossover recombination. This also solves the problem of finding the ``most tree-like” network when a perfect phylogeny is not possible. In this algorithm, no ancestral sequence is known in advance.

Blobbed-trees: generalizing galled-trees • In a phylogenetic network a maximal set of intersecting cycles is called a blob. • Contracting each blob results in a directed, rooted tree, otherwise one of the “blobs” was not maximal. • So every phylogenetic network can be viewed as a directed tree of blobs - a blobbed-tree. The blobs are the non-tree-like parts of the network.

Every network is a tree of blobs. How do the tree parts and the blobs relate? How can we exploit this relationship? Ugly tangled network inside the blob.

The start of technical stuff

Site Conflicts A pair of sites (columns) of M that fail the 4 -gametes test are said to conflict. And each site in the pair is said to be conflicted. A site that is not in such a pair is unconflicted.

M 12345 a 00010 b 10010 c 00100 d 10100 e 01100 f 01101 g 00101 Conflict Graph 4 1 3 2 5 Two nodes are connected iff the pair of sites conflict, i. e, fail the 4 -gamete test. THE MAIN TOOL: We represent the pairwise conflicts in a conflict graph.

Simple Fact If sites two sites i and j conflict, then the sites must be together on some recombination cycle whose recombination point is between the two sites i and j. (This is a general fact for all phylogenetic networks. ) Ex: In the prior example, site 1 conflicts with 3 and 4; and site 2 conflicts with 5.

A Phylogenetic Network 00010 a: 00010 10010 00000 4 3 1 00100 2 b: 10010 P c: 00100 S 3 01100 p d: 10100 5 00101 S 4 01101 f: 01101 e: 01100 g: 00101

Simple Consequence of the simple fact All sites on the same (non-trivial) connected component of the conflict graph must be on the same blob in any blobbed-tree. Follows by transitivity. So we can’t subdivide a blob into a tree-like structure if it only contains sites from a single connected component of the conflict graph.

Key Result about Galls: For galls, the converse of the simple consequence is also true. Two sites that are in different (non-trivial) connected components cannot be placed on the same gall in any phylogenetic network for M. Hence, in a galled-tree T for M each gall contains all and only the sites of one (non-trivial) connected component of the conflict graph. All unconflicted sites can be put on edges outside of the galls. This was the key to the galled-tree solution.

A galled-tree generating the sequences generated by the prior network. 4 Conflict Graph 4 3 1 1 3 2 5 p s a: 00010 2 b: 10010 c: 00100 d: 10100 2 p e: 01100 5 4 s f: 01101 g: 00101

Reduced Galled-Trees • A galled-tree is called reduced if every gall only contains conflicted sites. • Theorem: If M can be derived on a galled-tree, it can be derived on a reduced galled-tree. • The number of recombination nodes in a reduced galled-tree equals the number of connected components of the conflict graph, which is the minimum number of recombinations possible in any galled-tree.

• A reduced galled-tree for M (if one exists) minimizes the number of recombinations used over any phylogenetic network for M. The initial proof is based on a lower bound method. • Another proof (together with Dean Hickerson): The minimum number of recombinations needed in any phylogenetic network for M is at least the number of non-trivial connected components of the conflict graph for M.

The main new fact For any set of sequences M, there is a blobbed-tree T(M) that derives M, where each blob contains all and only the sites in one non-trivial connected component of the conflict graph. The unconflicted sites can always be put on edges outside of any blob. Moreover, the tree part of T(M) is unique. This is bit weaker than the result for galled-trees: it replaces “must” with “can”.

My proof is direct and self-contained, but the result may also be derivable from well-established results about splits, for example, from facts known about Bunemann graphs. Thanks to Mike Steel for pointing this out to me.

What does T(M) look like? How do we construct it? How do we exploit it? Let T’(M) be the tree where each blob in T(M) has been contracted to a single node. Finding T’(M) is easy from M. Then we expand each node in T’(M) to get T(M).

How to find T’(M) For a connected component C in the conflict graph, let M[C] be M restricted to the sites in C. Consider each distinct sequence S in M[C] as defining a binary character (split) partitioning those rows of M[C] that contain S, and those that do not. Let M’ be the binary matrix obtained from these characters, over all connected components in the conflict graph.

Fact: M’ satisfies the 4 -gamete test and so has a unique undirected perfect phylogeny, which is T’(M). Then we need to expand each node in T’(M) that corresponds to a blob for connected component C to some network B so that the external nodes in B have the labels in M[C], when restricted to C.

T(M) v sites on B from C w B’ Leaf set Tv Key point: For any leaf in Tv, and only at those leaves, the states of the sites in C are the same as at v. B: blob C: component of the conflict graph whose sites are on B

Example: 4 3 v 1 134 010 p s a: 00010 3 b: 10010 c: 00100 d: 10100 2 p e: 01100 5 4 s f: 01101 g: 00101

More Formally: Let M[C] be the matrix M restricted to the sites in C. Let S[C] be a sequence S restricted to the sites in C. Key point: Each distinct (non-zero) sequence in M[C] must be the sequence S[C] for some sequence S labeling a branching node v on B. So, although we don’t know much about the interior of B, or the arrangement of the exterior nodes (braching off of B), we know precisely their number and the sequences (restricted to sites of C) that label the exterior nodes on B. And we know that the states of the non-C sites are identical at each node in B.

M 134 a 001 b 101 c 010 d 110 e 010 f 010 g 01 0 12345 a 00010 b 10010 c 00100 d 10100 e 01100 f 01101 g 00101 4 C 1 3 000 4 B 3 001 010 1 Matrix M[C] is Matrix M restricted to the columns in C. a 101 b p s 2 110 c, e, f, g d

Punch Line • Each distinct non-zero sequence S[C] in M[C] corresponds to an edge e branching off of B, and exactly defines the set of sequences in M (and hence leaves in T(M)) that must be below e in T(M). • If S[C] is the all-zero sequence, then the leaf labeled S must connect to B through the coalesenct node of B. These two facts allow the construction of the tree part of T(M) in O(nm) time, starting from M. It also defines the complete label of each exterior node.

M 12345 a 00010 b 10010 c 00100 d 10100 e 01100 f 01101 g 00101 134 a 001 b 101 c 010 d 110 e 010 f 010 g 01 0 25 00 00 10 11 01 a b c d e f g M[C 1] M[C 2] 4 C 2 C 1 1 3 B 1 2 5 B 2 So the paths to every leaf pass through the blob B 1, but only the paths to e, f, g pass through gall B 2. The path from B 1 to B 2 exits B 1 at the node whose C 1 -restricted label is 010.

M 1 2 3 4 3 3 3 134 a 001 b 101 c 010 d 110 e 010 f 010 g 01 0 M[C 1] 12345 a 00010 b 10010 c 00100 d 10100 e 01100 f 01101 g 00101 5 5 6 7 8 a b c d e f g 25 00 00 10 11 01 M[C 2] 4 C 2 C 1 1 3 a b c d e f g 2 5 12345678 1000 01001000 00101000 00011000 00100100 00100001 M’

Algorithmically • Finding the tree part of the blobbed-tree is easy. • Determining the sequences labeling the exterior nodes on any blob is easy. • Determining a “good” structure inside a blob B is the problem of generating the sequences of the exterior nodes of B. • It is easy to test whether the exterior sequences on B can be generated with only a single (possibly multiple-crossover) recombination. The original galled-tree problem is now just the problem of testing whether one single-crossover recombination is sufficient for each blob.

The main open question • Is there always a blobbed-tree where each blob has all and only the sites of a single connected component of the conflict graph, which minimizes the number of recombinations over all possible phylogenetic networks for M? • Proof attempt by splitting blob B into B’ and B”. The catch is the possibility of a recombination node in B that is needed in both B’ and B”.

If Yes, Then any method that computes a lower bound on the number of needed recombinations should be applied separately on sequences defined by the sites on each connected component, and the results added together. This may be true even if the desired result is not. It is worth trying to prove this.

How to arrange the sites on a blob Given a single connected component of the conflict graph with k sites, how do we arrange those k sites to generate the required sequences, using only one (multiple-crossover) recombination, or using a multiple-crossover recombination with the fewest cross-overs?

Arranging the sites We will describe an O(n^3) time method to arrange all of the blobs. O(n^2) time is possible with a more complex method when only single-crossover recombinations are allowed.

A needed fact in words Let Q be a gall for the sites on connected-component C of the conflict graph. Let M[C] be the matrix M restricted to the sites on C. Let LQ[C] be the sequences labeling the nodes of Q, restricted to the sites on C. Claim: The two sets of sequences are identical, i. e. , M[C] = LQ[C].

M 134 a 001 b 101 c 010 d 110 e 010 f 010 g 01 0 12345 a 00010 b 10010 c 00100 d 10100 e 01100 f 01101 g 00101 LQ[C] are the node labels on Q restricted to the sites in C 4 C 1 3 LQ[C] 4 Q 3 001 010 1 Matrix M[C] is Matrix M restricted to the columns in C. Fact: M[C] = LQ[C] a 101 b p s 2 110 c, e, f, g d

The idea for arranging the sites of C on B: via a short movie

4 Q 3 001 010 1 a 101 b p s 2 110 c, e, f, g d

4 Q 3 001 010 1 a 101 b 110 d c, e, f, g

B 4 001 1 a 3 010 101 b: 101 110 c, e, f, g: 010 d Blob B minus the recombination node is a perfect phylogeny for M[C] minus the recombinant sequence; all sites are on one or two paths from the root; and the two end sequences of those paths can recombine at point r to recreate the recombinant sequence.

The point If we remove the recombinant node from B, we have a phylogenetic tree (no cycles) for the remaining sequences and hence a perfect phylogenetic tree for the sequences in M[C] minus the recombinant sequence of B. The sites in this tree are on one or two paths. Moreover, the two end sequences on that perfect phylogeny can recombine to create the removed recombinant sequence.

The algorithm for arranging a blob B for C 1. Form the matrix M[C]. 2. For each row of M[C], remove the row, see if there is a perfect phylogeny for the remaining rows. If yes, see if the sites are in one or two paths, and the end sequences can generate the removed row by a recombination. Fact: Every row that works gives a permitted

Optimality Theorem: The minimum number of recombination nodes in any phylogenetic network for M is at least the number of non-trivial connected components of the conflict graph. Hence, when the sequences on each blob on T(M) can be generated with a single recombination node, the blobbed-tree minimizes the number of recombination nodes over all phylogenetic networks.

The number of arrangements on a gall (all-0 ancestral sequence) Following the algorithmic approach above, one can prove that the number of arrangements of any gall is at most three, and this happens only if the gall has two sites. If the gall has more than two sites, then the number of arrangements is at most two. If the gall has four or more sites, with at least two sites on each side of the recombination point (not the side of the gall) then the arrangement is forced and unique.

What is the most tree-like Network? Assume no node with only one child. Contracting each blob to a single node results in a tree. The number of edges in that tree is a measure of how tree-like the network is. The larger the number of edges, the more tree-like is the network. Then T(M) is the most tree-like network.

OPEN QUESTIONS • Main open question - stated earlier • PPH problem when the haplotypes were derived on a galled-tree or blobbed tree rather than a tree. • Better lower bounds on the number of recombinations based on a finer examination of the structure of the conflict graph.