 Research
 Open Access
Reconciling multiple genes trees via segmental duplications and losses
 Riccardo Dondi^{1},
 Manuel Lafond^{2}Email authorView ORCID ID profile and
 Celine Scornavacca^{3}
https://doi.org/10.1186/s1301501901396
© The Author(s) 2019
 Received: 5 November 2018
 Accepted: 23 February 2019
 Published: 20 March 2019
Abstract
Reconciling gene trees with a species tree is a fundamental problem to understand the evolution of gene families. Many existing approaches reconcile each gene tree independently. However, it is wellknown that the evolution of gene families is interconnected. In this paper, we extend a previous approach to reconcile a set of gene trees with a species tree based on segmental macroevolutionary events, where segmental duplication events and losses are associated with cost \(\delta \) and \(\lambda \), respectively. We show that the problem is polynomialtime solvable when \(\delta \le \lambda \) (via LCAmapping), while if \(\delta > \lambda \) the problem is NPhard, even when \(\lambda = 0\) and a single gene tree is given, solving a long standing open problem on the complexity of multigene reconciliation. On the positive side, we give a fixedparameter algorithm for the problem, where the parameters are \(\delta /\lambda \) and the number d of segmental duplications, of time complexity \(O\left(\lceil \frac{\delta }{\lambda } \rceil ^{d} \cdot n \cdot \frac{\delta }{\lambda }\right)\). Finally, we demonstrate the usefulness of this algorithm on two previously studied real datasets: we first show that our method can be used to confirm or raise doubt on hypothetical segmental duplications on a set of 16 eukaryotes, then show how we can detect whole genome duplications in yeast genomes.
Keywords
 Phylogenetics
 Gene trees
 Species trees
 Reconciliation
 Segmental duplications
 Fixedparameter tractability
 NPhardness
 Whole genome duplications
Introduction
It is nowadays well established that the evolution of a gene family can differ from that of the species containing these genes. This can be due to quite a number of different reasons, including gene duplication, gene loss, horizontal gene transfer or incomplete lineage sorting, to only name a few [22]. While this discongruity between the gene phylogenies (the gene trees) and the species phylogeny (the species tree) complicates the process of reconstructing the latter from the former , every cloud has a silver lining: “plunging” gene trees into the species tree and analyzing the differences between these topologies, one can infer the macroevolutionary events that shaped gene evolution. This is the rationale behind the species treegene tree reconciliation, a concept introduced in [13] and first formally defined in [24]. Understanding these macroevolutionary events allows us to better understand the mechanisms of evolution with applications ranging from orthology detection [9, 18, 19, 30] to ancestral genome reconstruction [11], and recently in dating phylogenies [5, 7].
It is wellknown that the evolution of gene families is interconnected. However, in current pipelines, each gene tree is reconciled independently with the species tree, even when posterior to the reconciliation phase the genes are considered as related, e.g. [11]. A more pertinent approach would be to reconcile the set of gene trees at once and consider segmental macroevolutionary events, i.e. events that concern a chromosome segment instead of a single gene.
Some work has been done in the past to model segmental gene duplications and three models have been considered: the ec (Episode Clustering) problem, the me (Minimum Episodes) problem [1, 14], and the mgd (Multiple Gene Duplication) problem [12]. The ec and mgd problems both aim at clustering duplications together by minimizing the number of locations in the species tree where at least one duplication occurred, with the additional requirement that a cluster cannot contain two gene duplications from a same gene tree in the mgd problem. The me problem is more biologicallyrelevant, because it aims at minimizing the actual number of segmental duplications (more details in "Reconciliation with segmental duplications" section). Most of the exact solutions proposed for the me problem [1, 20, 26] deal with a constrained version, since the possible mappings of a gene tree node are limited to given intervals, see for example [1, Def. 2.4]. In [26], a simple \(O^*(2^{k})\) time algorithm is presented for the unconstrained version (here \(O^*\) hides polynomial factors), where k is the number of speciation nodes that have descending duplications under the LCAmapping. This shows that the problem is fixedparameter tractable (FPT) in k. But since the LCAmapping maximizes the number of speciation nodes, there is no reason to believe that k is a small parameter, and so more practical FPT algorithms are needed. Recently, Delabre et al. [8] studied the problem of reconstructing the evolution of syntenic blocks. Their model allows segmental duplications, but is more constrained since it requires every gene family of every block to have evolved along the same tree.
In this paper, we extend the unconstrained ME model to gene losses and provide a variety of new algorithmic results. We allow weighing segmental duplication events and loss events by separate costs \(\delta \) and \(\lambda \), respectively. We show that if \(\delta \le \lambda \), then an optimal reconciliation can be obtained by reconciling each gene tree separately under the usual LCAmapping, even in the context of segmental duplications. On the other hand, we show that if \(\delta > \lambda \) and both costs are given, reconciling a set of gene trees while minimizing segmental gene duplications and gene losses is NPhard. The hardness also holds in the particular case that we ignore losses, i.e. when \(\lambda = 0\). This solves a long standing open question on the complexity of the reconciliation problem under this model (in [1], the authors already said “it would be interesting to extend the [...] model of Guigó et al. (1996) by relaxing the constraints on the possible locations of gene duplications on the species tree”. The question is stated as still open in [26]). The hardness holds also when only a single gene tree is given. On the positive side, we describe an algorithm that is practical when \(\delta \) and \(\lambda \) are not too far apart. More precisely, we show that multigene tree reconciliation is fixedparameter tractable in the ratio \(\delta /\lambda \) and the number d of segmental duplications, and can be solved in time \(O\left(\lceil \frac{\delta }{\lambda } \rceil ^{d} \cdot n \cdot \frac{\delta }{\lambda }\right)\). The algorithm has been implemented and tested and is freely available^{1} at https://github.com/AEVOlab/MultRec. We first evaluate the potential of multigene reconciliation on a set of 16 eukaryotes, and show that our method can find scenarios with less duplications than other approaches. While some previously identified segmental duplications are confirmed by our results, it casts some doubt on others as they do not occur in our optimal scenarios. We then show how the algorithm can be used to detect whole genome duplications in yeast genomes. Further work includes incorporating in the model segmental gene losses and segmental horizontal gene transfers, with a similar flavor than the heuristic method discussed in [10].
Preliminaries
Basic notions
For our purposes, a rooted phylogenetic tree \(T=(V(T),E(T))\) is an oriented tree, where V(T) is the set of nodes, E(T) is the set of arcs, all oriented away from r(T), the root. Unless stated otherwise, all trees in this paper are rooted phylogenetic trees. A forest \(F = (V(F),E(F))\) is a directed graph in which every connected component is a tree. Denote by t(F) the set of trees of F that are formed by its connected components. Note that a tree is itself a forest. In what follows, we shall extend the usual terminology on trees to forests.
For an arc (x, y) of F, we call x the parent of y, and y a child of x. If there exists a path that starts at x and ends at y, then x is an ancestor of y and y is a descendant of x. We say y is a proper descendant of x if \(y \ne x\), and then x is a proper ancestor of y. This defines a partial order denoted by \(y \le _{F} x\), and \(y <_F x\) if \(x\not =y\) (we may omit the F subscript if clear from the context). If none of \(x \le y\) and \(y \le x\) holds, then x and y are incomparable. The set of children of x is denoted ch(x) and its parent x is denoted par(x) [which is defined to be x if x itself is a root of a tree in t(F)]. For some integer \(k \ge 0\), we define \(par^k(x)\) as the kth parent of x. Formally, \(par^0(x) = par(x)\) and \(par^k(x) = par(par^{k  1}(x))\) for \(k > 0\). The number of children ch(x) of x is called the outdegree of x. Nodes with no children are leaves, all others are internal nodes. The set of leaves of a tree F is denoted by L(F). The leaves of F are bijectively labeled by a set \(\mathcal{L} (F)\) of labels. A forest is binary if \(ch(x)=2\) for all internal nodes x. Given a set of nodes X that belong to the same tree \(T \in t(F)\), the lowest common ancestor of X, denoted \(LCA_F(X)\), is the node z that satisfies \(x \le z\) for all \(x \in X\) and such that no child of z satisfies this property. We leave \(LCA_F(X)\) undefined if no such node exists (when elements of X belong to different trees of t(F)). We may write \(LCA_F(x,y)\) instead of \(LCA_F(\{x,y\})\). The height of a forest F, denoted h(F), is the number of nodes of a longest directed path from a root to a leaf in a tree of F (note that the height is sometimes defined as the number of arcs on such a path—here we use the number of nodes instead). Observe that since a tree is a forest, all the above notions also apply on trees.
Reconciliations
A reconciliation usually involves two rooted phylogenetic trees, a gene tree G and a species tree S, which we always assume to be both binary. In what follows, we will instead define reconciliation between a gene forest \({\mathcal {G}}\) and a species tree. Here \({\mathcal {G}}\) can be thought of as a set of gene trees. Each leaf of \({\mathcal {G}}\) represents a distinct extant gene, and \({\mathcal {G}}\) and S are related by a function \(s: {\mathcal{L}} ({\mathcal {G}}) \rightarrow {\mathcal{L}} (S)\), which means that each extant gene belongs to an extant species. Note that s does not have to be injective (in particular, several genes from a same gene tree G of \({\mathcal {G}}\) can belong to the same species) or surjective (some species may not contain any gene of \({\mathcal {G}}\)). Given \({\mathcal {G}}\) and S, we will implicitly assume the existence of the function s.
In a \({{\mathbb {D}}}{{\mathbb {L}}}\) reconciliation, each node of \({\mathcal {G}}\) is associated to a node of S and an event—a speciation (\({\mathbb {S}}\)), a duplication (\({\mathbb {D}}\)) or a contemporary event (\({\mathbb {C}}\))—under some constraints. A contemporary event \({\mathbb {C}}\) associates a leaf u of \({\mathcal {G}}\) with a leaf x of S such that \(s(u) = x\). A speciation in a node u of \({\mathcal {G}}\) is constrained to the existence of two separated paths from the mapping of u to the mappings of its two children, while the only constraint given by a duplication event is that evolution of \({\mathcal {G}}\) cannot go back in time. More formally:
Definition 1
 1
if u is a leaf of \({\mathcal {G}}\), then \(\alpha _e(u) = {\mathbb {C}}\) and \(\alpha _r(u) = s(u)\);
 2if u is an internal node of \({\mathcal {G}}\) with children \(u_1, u_2\), then exactly one of following cases holds:

\(\alpha _e(u) = {\mathbb {S}}\), \(\alpha _r(u) = LCA_S(\alpha _r(u_1),\alpha _r(u_2))\) and \(\alpha _r(u_1),\alpha _r(u_2)\) are incomparable;

\(\alpha _e(u) = {\mathbb {D}}\), \(\alpha _r(u_1) \le \alpha _r(u)\) and \(\alpha _r(u_2) \le \alpha _r(u)\)


if \(\alpha _e(u)={\mathbb {S}}\), then \(l_{\alpha }(u) = dist(\alpha _r(u), \alpha _r(u_1)) +\,dist(\alpha _r(u), \alpha _r(u_2))  2\);

if \(\alpha _e(u)={\mathbb {D}}\), then \(l_{\alpha }(u) = dist(\alpha _r(u), \alpha _r(u_1)) +\,dist(\alpha _r(u), \alpha _r(u_2))\).
Reconciliation with segmental duplications
Here we want to tackle the problem of reconciling several gene trees at the same time and counting segmental duplications only once. Given a set of duplications nodes \({\mathcal {D}} \in V({\mathcal {G}})\) occurring in a given node s of the species tree, it is easy to see that the minimum number of segmental duplications associated with s is the minimal number of parts in a partition of \({\mathcal {D}}\) in which each part does not contain comparable nodes. See Fig. 1(4) for an example. This number coincides [1] with \(h_{\alpha }(s) := h({\mathcal {G}}[\alpha , s])\), i.e. the height of the forest of the duplications in s. Now, denote \({\hat{d}}(\alpha ) = \sum _{s \in V(S)} h_{\alpha }(s)\). For instance in Fig. 1, under the mapping \(\mu \) in (2), we have \({\hat{d}}(\mu ) = 6\), because \(h_{\mu }(s) = 1\) for \(s \in \{A,B,C,E\}\) and \(h_{\mu }(F) = 2\). But under the mapping \(\alpha \) in (3), \({\hat{d}}(\alpha ) = 4\), since \(h_{\alpha }(A) = 1\) and \(h_{\alpha }(F) = 3\). Note that this does not consider losses though—the \(\alpha \) mapping has more losses than \(\mu \).
The cost of \(\alpha \) is \(cost^{SD}({\mathcal {G}}, S, \alpha ) = \delta \cdot {\hat{d}}(\alpha ) + \lambda \cdot l(\alpha )\). If \({\mathcal {G}}\) and S are unambiguous, we may write \(cost^{SD}(\alpha )\). We have the following problem:
Problem 1
Most parsimonious reconciliation of a set of trees with segmental duplications (MPRSTSD)
Instance: a species tree S, a gene forest \({\mathcal {G}}\), costs \(\delta \) for duplications and \(\lambda \) for losses.
Output: a reconciliation \(\alpha \) of \({\mathcal {G}}\) in S such that \(cost^{SD}({\mathcal {G}}, S, \alpha )\) is minimum.
Note that, when \(\lambda =0\), \(cost^{SD}\) coincides with the unconstrained ME score defined in [26] (where it is called ME under the FHS model).
Properties of multigene reconciliations
We finish this section with some additional terminology and general properties of multigene reconciliations that will be useful throughout the paper. The next basic result states that in a reconciliation \(\alpha \), we should set the events of internal nodes to \({\mathbb {S}}\) whenever it is allowed.
Lemma 1
Let \(\alpha \) be a reconciliation for \({\mathcal {G}}\) in S, and let \(u \in V({\mathcal {G}})\) such that \(\alpha _e(u) = {\mathbb {D}}\). Let \(\alpha '\) be identical to \(\alpha \), with the exception that \(\alpha '_e(u) = {\mathbb {S}}\), and suppose that \(\alpha '\) satisfies the requirements of Definition 1. Then \(cost^{SD}(\alpha ') \le cost^{SD}(\alpha )\).
Proof
Observe that changing \(\alpha _e(u)\) from \({\mathbb {D}}\) to \({\mathbb {S}}\) cannot increase \({\hat{d}}(\alpha )\). Moreover, as \(dist(\alpha '_r(u), \alpha '_r(u_1))\) and \(dist(\alpha '_r(u), \alpha '_r(u_2))\) are the same as in \(\alpha \) for the two children \(u_1\) and \(u_2\) of u, by definition of duplications and losses this decreases the number of losses by 2. Thus \(cost^{SD}(\alpha ') \le cost^{SD}(\alpha )\), and this inequality is strict when \(\lambda > 0\). \(\square \)
Since we are looking for a most parsimonious reconciliation, by Lemma 1 we may assume that for an internal node \(u \in V({\mathcal {G}})\), \(\alpha _e(u) = {\mathbb {S}}\) whenever allowed, and \(\alpha _e(u) = {\mathbb {D}}\) otherwise. Therefore, \(\alpha _e(u)\) is implicitly defined by the \(\alpha _r\) mapping. To alleviate notation, we will treat \(\alpha \) as simply as a mapping from \(V({\mathcal {G}})\) to V(S) and thus write \(\alpha (u)\) instead of \(\alpha _r(u)\). We will assume that the events \(\alpha _e(u)\) can be deduced from this mapping \(\alpha \) by Lemma 1.
Therefore, treating \(\alpha \) as a mapping, we will say that \(\alpha \) is valid if for every \(v \in V({\mathcal {G}})\), \(\alpha (v) \ge \alpha (v')\) for all descendants \(v'\) of v. We denote by \(\alpha [v \rightarrow s]\) the mapping obtained from \(\alpha \) by remapping \(v \in V({\mathcal {G}})\) to \(s \in V(S)\), i.e. \(\alpha [v \rightarrow s](w) = \alpha (w)\) for every \(w \in V({\mathcal {G}}) \setminus \{v\}\), and \(\alpha [v \rightarrow s](v) = s\). Since we are assuming that \({\mathbb {S}}\) and \({\mathbb {D}}\) events can be deduced from \(\alpha \), the LCAmapping is now unique: we denote by \(\mu : V({\mathcal {G}}) \rightarrow V(S)\) the LCAmapping, defined as \(\mu (v) = s(v)\) if \(v \in L({\mathcal {G}})\), and otherwise \(\mu (v) = LCA_{S}( \mu (v_1), \mu (v_2) )\), where \(v_1\) and \(v_2\) are the children of v. Note that for any valid reconciliation \(\alpha \), we have \(\alpha (v) \ge \mu (v)\) for all \(v \in V({\mathcal {G}})\). We also have the following, which will be useful to establish our results.
Lemma 2
Let \(\alpha \) be a mapping from \({\mathcal {G}}\) to S. If \(\alpha (v) > \mu (v)\), then v is a \({\mathbb {D}}\) node under \(\alpha \).
Proof
Let \(v_1\) and \(v_2\) be the two children of v. If \(\alpha (v) \ne LCA_{S}(\alpha (v_1), \alpha (v_2))\), then v must be a duplication, by the definition of \({\mathbb {S}}\) events. The same holds if \(\alpha (v_1)\) and \(\alpha (v_2)\) are not incomparable. Thus assume \(\alpha (v) = LCA_{S}(\alpha (v_1), \alpha (v_2)) > \mu (v)\) and that \(\alpha (v_1)\) and \(\alpha (v_2)\) are incomparable. This implies that one of \(\alpha (v_1)\) or \(\alpha (v_2)\) is incomparable with \(\mu (v)\), say \(\alpha (v_1)\) w.l.o.g. But \(\mu (v_1) \le \alpha (v_1)\), implying that \(\mu (v_1)\) is also incomparable with \(\mu (v)\), a contradiction to the definition of \(\mu = LCA_S(\mu (v_1), \mu (v_2))\). \(\square \)
Lemma 3
Let \(\alpha \) be a mapping from \({\mathcal {G}}\) to S, and let \(v \in V({\mathcal {G}})\). Suppose that there is some proper descendant \(v'\) of v such that \(\alpha (v') \ge \mu (v)\). Then v is a duplication under \(\alpha \).
Proof
If \(\alpha (v) = \mu (v)\), we get \(\mu (v) \le \alpha (v') \le \alpha (v) = \mu (v)\), and so \(\alpha (v') = \mu (v)\). We must then have \(\alpha (v'') = \mu (v)\) for every node \(v''\) on the path between \(v'\) and v. In particular, v has a child \(v_1\) with \(\alpha (v) = \alpha (v_1)\) and thus v is a duplication. If instead \(\alpha (v) > \mu (v)\), then v is a duplication by Lemma 2. \(\square \)
Lemma 4
(Shiftdown lemma) Let \(\alpha \) be a mapping from \({\mathcal {G}}\) to S, let \(v \in V({\mathcal {G}})\), let \(s \in V(S)\) and \(k > 0\) be such that \(par^k(s) = \alpha (v)\). Suppose that \(\alpha [v \rightarrow s]\) is a valid mapping. Then \(l(\alpha [v \rightarrow s]) \le l(\alpha )  k\).
Proof
Let \(v_1\) and \(v_2\) be the children of v, and denote \(t := \alpha (v), t_1 := \alpha (v_1)\) and \(t_2 := \alpha (v_2)\). Moreover denote \(\alpha ' := \alpha [v \rightarrow s]\). Let P be the set of nodes that appear on the path between s and t, excluding s but including t (note that s is a proper descendant of t but an ancestor of both \(t_1\) and \(t_2\), by the validity of \(\alpha '\)). For instance in Fig. 2, \(P = \{x,t\}\). Observe that \(P = k\). Under \(\alpha \), there is a loss for each node of P on both the \((v, v_1)\) and \((v, v_2)\) branches. (noting that v is a duplication by Lemma 2). These 2k losses are not present under \(\alpha '\). On the other hand, there are at most k losses that are present under \(\alpha '\) but not under \(\alpha \), which consist of one loss for each node of P on the (par(v), v) branch (in the case that v is not the root of its tree—otherwise, no such loss occurs). This proves that \(l(\alpha ') \le l(\alpha )  k\). \(\square \)
A “proofbypicture” of the above Lemma appears in Fig. 2. When shifting down the mapping of v, some losses that appear left and right of v get “combined” on the branch above it.
The computational complexity of the mprstsd problem
We separate the study of the complexity of the mprstsd problem into two subcases: when \(\lambda \ge \delta \) and when \(\lambda < \delta \).
The case of \(\lambda \ge \delta \)
The following theorem states that, when \(\lambda \ge \delta \), the MPR (ie the LCAmapping) is a solution to the mprstsd problem.
Theorem 1
Let \({\mathcal {G}}\) and S be an instance of mprstsd, and suppose that \(\lambda \ge \delta \), Then the LCAmapping \(\mu \) is a reconciliation of minimum cost for \({\mathcal {G}}\) and S. Moreover if \(\lambda > \delta \), \(\mu \) is the unique reconciliation of minimum cost.
Proof
Let \(\alpha \) be a mapping of \({\mathcal {G}}\) into S of minimum cost. Let \(v \in V({\mathcal {G}})\) be a minimal node of \({\mathcal {G}}\) with the property that \(\alpha (v) \ne \mu (v)\) (i.e. all proper descendants \(v'\) of v satisfy \(\alpha (v) = \mu (v)\)). Note that v must exists since, for every leaf \(l \in \mathcal{L} ({\mathcal {G}})\), we have \(\alpha (l) = \mu (l)\). Because \(\alpha (v) \ge \mu (v)\), it follows that \(\alpha (v) > \mu (v)\). Denote \(s = \mu (v)\) and \(t = \alpha (v)\). Then there is some \(k \ge 1\) such that \(t = par^k(s)\). Consider the mapping \(\alpha ' = \alpha [v \rightarrow s]\). This possibly increases the sum of duplications by 1, so that \({\hat{d}}(\alpha ') \le {\hat{d}}(\alpha ) + 1\). But by the Shiftdown lemma, \(l(\alpha ') \le l(\alpha )  1\). Thus we have at most one duplication but save at least one loss.
If \(\lambda > \delta \), this contradicts the optimality of \(\alpha \), implying that v cannot exist and thus that \(\alpha = \mu \). This proves the uniqueness of \(\mu \) in this case.
If \(\delta = \lambda \), then \(\delta {\hat{d}}(\alpha ') + \lambda l(\alpha ') \le \delta {\hat{d}}(\alpha ) + \lambda l(\alpha )\). By applying the above transformation successively on the minimal nodes v that are not mapped to \(\mu (v)\), we eventually reach the LCAmapping \(\mu \) with an equal or better cost than \(\alpha \). \(\square \)
The case of \(\delta > \lambda \)
We show that, in contrast with the \(\lambda \ge \delta \) case, the mprstsd problem is NPhard when \(\delta > \lambda \) and the costs are given as part of the input. More specifically, we show that the problem is NPhard when one only wants to minimize the sum of duplication heights, i.e. \(\lambda = 0\). Note that if \(\lambda > 0\) but is small enough, the effect will be the same and the hardness result also holds—for instance, putting \(\delta = 1\) and say \(\lambda < \frac{1}{2V({\mathcal {G}})V(S)}\) ensures that even if a maximum number of losses appears on every branch of \({\mathcal {G}}\), it does not even amount to the cost of one duplication.
We briefly outline the main ideas of the reduction. We start from the Vertex Cover problem, where we are given a graph G and must find a subset of vertices \(V' \subseteq V(G)\) of minimum size such that each edge has at least one endpoint in \(V'\). The species tree S and the forest \({\mathcal {G}}\) are constructed so that, for each vertex \(v_i \in V(G)\), there is a gene tree \(A_i\) in \({\mathcal {G}}\) with a long path of duplications, all of which could either be mapped to a species called \(y_i\) or another species \(z_i\). We make it so that mapping to \(y_i\) introduces one more duplication than mapping to \(z_i\), hence we have to “pay” for each \(y_i\). On the other hand, we construct a tree \(C_h\) in \({\mathcal {G}}\) for each edge in \(e_h = v_iv_j \in E(G)\) such that if, in \(A_i\) or \(A_j\), we chose to map to either \(y_i\) or \(y_j\), then \(C_h\) can “reuse” these \(y_i\) or \(y_j\) duplications with no extra cost. However if we did not choose either \(y_i\) or \(y_j\), \(C_h\) will introduce a large number of new duplications. We are therefore forced to pick a minimum number of \(y_i\)’s to “cover” all the \(C_h\) trees.
The construction
We will first need particular trees as described by the following Lemma. These trees guarantee that a prescribed set of leaves L are at distance exactly k from the root, and any two of the leaves in L have their LCA at distance at least k / 2. Recall that for a tree T and \(u,v \in V(T)\), \(dist_T(u,v)\) denotes the number of edges on the path between u and v in T (we write dist(u, v) for short). A caterpillar is a binary rooted tree in which each internal node has one child that is a leaf, with the exception of one internal node which has two such children.
Lemma 5
Let \(k \ge 8\) be an integer, and let L be a given set of at most k labels. Then there exists a rooted tree T with leaf set \(L'\) with \(L \subseteq L'\) such that for any \(l \in L\), \(dist(l, r(T)) = k\) and for any two distinct \(l_1, l_2 \in L\), \(dist(l_1, LCA_{T}(l_1, l_2)) \ge k/2\). Moreover, T can be constructed in polynomial time with respect to k.
Proof
Let p be the smallest integer such that \(k \le 2^p\). First consider a fully balanced binary tree T on \(2^{p}\) leaves, so that each leaf is at distance p from the root. Then replace each leaf by the root of a caterpillar of height \(k  p + 1\) (hence in the caterpillars the longest roottoleaf path has \(k  p\) edges). The resulting tree is T. Choose k of these caterpillars, and in each of them, assign a distinct label of L to a deepest leaf (any of the two). Thus \(dist(l, r(T)) = k\) for each \(l \in L\), and clearly T can be built in polynomial time. We also have \(dist(l_1, LCA_{S}(l_1, l_2)) \ge k  p\) for each distinct \(l_1, l_2 \in L\). As \(2^p \le 2k\), we have \(p \le \log (2k)\). This implies \(k  p \ge k  \log (2k) \ge k/2\) for \(k \ge 8\). \(\square \)
In the following, we will assume that \(\lambda = 0\) and \(\delta = 1\). We reduce the Vertex Cover problem to that of finding a mapping of minimum cost for given \({\mathcal {G}}\) and S. Recall that in the decision version of Vertex Cover, we are given a graph \(G = (V, E)\) and an integer \(\beta < n\) and are asked if there exists a subset \(V' \subseteq V\) with \(V' \le \beta \) such that every edge of E has at least one endpoint in \(V'\). For such a given instance, denote \(V = \{v_1, \ldots , v_n\}\) and \(E = \{e_1, \ldots , e_m\}\) (so that \(n = V\) and \(m = E\)). The ordering of the \(v_i\)’s and \(e_j\)’s can be arbitrary, but must remain fixed for the remainder of the construction.
We proceed with the construction of the set of gene trees \({\mathcal {G}}\). Most of the trees of \({\mathcal {G}}\) consist of a subset of the nodes of S, to which we graft additional leaves to introduce duplications—some terminology is needed before proceeding. For \(w \in V(S)\), deleting w consists in removing w and all its descendants from S, then contracting the possible resulting degree two vertex (which was the parent of w if \(p(w) \ne r(S)\)). If this leaves a root with only one child, we contract the root with its child. For \(X \subseteq \mathcal{L} (S)\), keeping X consists of successively deleting every node that has no descendant in X until none remains (the tree obtained by keeping X is sometimes called the restriction of S to X).

The A trees Let \(A = \{A_1, \ldots , A_n\}\), one tree for each vertex of G. For each \(i \in [n]\), obtain \(A_i\) by first taking a copy of S, then deleting all the special \(y_iz_i\) leaves. Then on the resulting \(z_iy_i\) branch, graft 10K leaves labeled \(q_i\). Then delete the child of \(z_i\) that is also an ancestor of \(r_i\) (removing \(z_i\) in the process). Figure 3 bottomleft might be helpful.
As a result, under the LCAmapping \(\mu \), \(A_i\) has a path of 10K duplications mapped to \(y_i\). One can choose whether to keep this mapping in \(A_i\), or to remap these duplications to \(z_i\).

The B trees Let \(B = \{B_1, \ldots , B_n\}\). For \(i \in [n]\), \(B_i\) is obtained from S by deleting all except \(10K  2\) of the special \(r_iz_i\) leaves, and grafting a leaf labeled \(r_i\) on the edge between \(r_i\) and its parent, thereby creating a single duplication mapped to \(r_i\) under \(\mu \).

The C trees Let \(C = \{C_1, C_2, \ldots , C_m\}\), where for \(h \in [m]\), \(C_h\) corresponds to edge \(e_h\). Let \(v_i, v_j\) be the two endpoints of edge \(e_h = \{v_i,v_j\}\), where \(i < j\). To describe \(C_h\), we list the set of leaves that we keep from S. Keep all the leaves of the \(P_i\) subtree of S, and keep a subset of the special leaves defined as follows:

Keep 9K of the special \(x_iy_i\) leaves;

Keep \((j  i  1)10K + K\) of the special \(y_iz_i\) leaves;

Keep 9K of the special \(x_jy_j\) leaves;

Keep all the special \(y_jz_j\) leaves.
No other leaves are kept. Next, in the tree obtained by keeping the aforementioned list of leaves, for each \(k \in [K]\) we graft, on the edge between \(c_{i,k}\) and its parent, another leaf labeled \(c_{i,k}\). Thus \(C_h\) has K duplications, all located at the bottom of the \(P_i\) subtree. This concludes our construction.

Now consider a \(C_h\) tree, \(e_h = \{v_i,v_j\}\). Under the LCAmapping, the \(c_{i,k}\) duplications at the bottom enforce an additional K duplications. This can be avoided by, say, mapping all these duplications to the same species. For instance, we could remap all these duplications to some \(y_k\) species of S. But in this case, because of Lemma 3, every node v of \(C_h\) above a \(c_{i,k}\) duplication for which \(\mu (v) \le y_k\) will become a duplication. This will create a duplication subtree D in \(C_h\) with a large height, and our goal will be to “reuse” the duplications we chose in the \(A_{k'}\) and \(B_{k'}\) trees. As it turns out, this “reuse” of duplications will be feasible only if some \(y_i\) or \(y_j\) has duplication height 10K. If this does not occur, any attempt at mapping the \(c_{i,k}\) nodes to a common species will induce a chain reaction of too many duplications created above.
The complete proof is somewhat technical. The interested reader can find it in Appendix.
Theorem 2
The mprstsd problem is NPhard for \(\lambda = 0\) and for given \(\delta > \lambda \).
The above hardness supposes that \(\delta \) and \(\lambda \) can be arbitrarily far apart. This leaves open the question of whether MPRSTSD is NPhard when \(\delta \) and \(\lambda \) are fixed constants—in particular when \(\delta = 1 + \epsilon \) and \(\lambda = 1\), where \(\epsilon < 1\) is some very small constant. We end this section by showing that the above hardness result persists event if only one gene tree is given. The idea is to reduce from the MPRSTSD shown hard just above. Given a species tree S and a gene forest \({\mathcal {G}}\), we make \({\mathcal {G}}\) a single tree by incorporating a large number of speciations (under \(\mu \)) above the root of each tree of \({\mathcal {G}}\) (modifying S accordingly), then successively joining the roots of two trees of \({\mathcal {G}}\) under a common parent until \({\mathcal {G}}\) has only one tree.
Theorem 3
The mprstsd problem is NPhard for \(\lambda = 0\) and for given \(\delta > \lambda \), even if only one gene tree is given as input.
Proof
Notice the following: in any mapping \(\alpha \) of T, the \(k  1\) internal nodes of the \(T'\) caterpillar must be duplications mapped to \(r(S')\), so that \(h_{\alpha }(r(S')) \ge k  1\). Also note that under the LCA mapping \(\mu \) for T and \(S'\), the only duplications other than those \(k  1\) mentioned above occur in the \(G_i\) subtrees. The \((\Rightarrow )\) is then easy to see: given \(\alpha \) such that \(cost({\mathcal {G}}, S, \alpha ) \le t\), we set \(\alpha '(v) = \alpha (v)\) for every node v of T that is also in \({\mathcal {G}}\) (namely the nodes of \(G_1, \ldots , G_k\)), and set \(\alpha '(v) = \mu (v)\) for every other node. This achieves a cost of \(t + k  1\).
As for the \((\Leftarrow )\) direction, suppose that \(cost^{SD}(T', S', \alpha ') \le t + k  1\) for some mapping \(\alpha '\). Observe that under the LCAmapping in T, each root of each \(G_i\) subtree has a path of \(2(t + k)\) speciations in its ancestors. If any node in a \(G_i\) subtree of T is mapped to \(r(S')\), then all these speciations become duplications (by Lemma 3), which would contradict \(cost^{SD}(T', S', \alpha ') \le t + k  1\). We may thus assume that no node belonging to a \(G_i\) subtree is mapped to \(r(S')\). Since \(h_{\alpha '}(r(S')) \ge k  1\), this implies that the restriction of \(\alpha '\) to the \(G_i\) subtrees has cost at most t.
More formally, consider the mapping \(\alpha ''\) from \({\mathcal {G}}\) to \(S'\) in which we put \(\alpha ''(v) = \alpha '(v)\) for all \(v \in V({\mathcal {G}})\). Then \(cost^{SD}({\mathcal {G}}, S', \alpha '') \le cost^{SD}(T, S', \alpha ')  (k  1) \le t\), because \(\alpha ''\) does not contain the top \(k  1\) duplications of \(\alpha '\), and cannot introduce longer duplication paths than in \(\alpha '\).
We are not done, however, since \(\alpha ''\) is a mapping from \({\mathcal {G}}\) to \(S'\), and not from \({\mathcal {G}}\) to S. Consider the set \(Q \subseteq V({\mathcal {G}})\) of nodes v of \({\mathcal {G}}\) such that \(\alpha ''(v) \in \overline{V(S)} := V(S') \setminus V(S)\). We will remap every such node to r(S) and show that this cannot increase the cost. Observe that if \(v \in Q\), then every ancestor of v in \({\mathcal {G}}\) is also in Q. Also, every node in Q is a duplication (by invoking Lemma 2).
Consider the mapping \(\alpha ^*\) from \({\mathcal {G}}\) to \(S'\) in which we put \(\alpha ^*(v) = \alpha ''(v)\) for all \(v \notin Q\), and \(\alpha ^*(v) = r(S)\) for all \(v \in Q\). It is not difficult to see that \(\alpha ^*\) is valid.
An FPT algorithm
In this section, we show that for costs \(\delta > \lambda \) and a parameter \(d > 0\), if there is an optimal reconciliation \(\alpha \) of cost \(cost^{SD}({\mathcal {G}}, S)\) satisfying \({\hat{d}}(\alpha ) \le d\), then \(\alpha \) can be found in time \(O(\lceil \frac{\delta }{\lambda }\rceil ^d \cdot n \cdot \frac{\delta }{\lambda })\).
In what follows, we allow mappings to be partially defined, and we use the \(\bot \) symbol to indicate undetermined mappings. The idea is to start from a mapping in which every internal node is undetermined, and gradually determine those in a bottomup fashion. We need an additional set of definitions. We will assume that \(\delta> \lambda > 0\) (although the algorithm described in this section can solve the \(\lambda = 0\) case by setting \(\lambda \) to a very small value).
We say that the mapping \(\alpha : V({\mathcal {G}}) \rightarrow V(S) \cup \{ \bot \}\) is a partial mapping if \(\alpha (l) = s(l)\) for every leaf \(l \in \mathcal{L} ({\mathcal {G}})\), and it holds that whenever \(\alpha (v) \ne \bot \), we have \(\alpha (v') \ne \bot \) for every descendant \(v'\) of v. That is, if a node is determined, then all its descendants also are. This also implies that every ancestor of a \(\bot \)node is also a \(\bot \)node. A node \(v \in V({\mathcal {G}})\) is a minimal \(\bot \)node (under \(\alpha \)) if \(\alpha (v) = \bot \) and \(\alpha (v') \ne \bot \) for each child \(v'\) of v. If \(\alpha (v) \ne \bot \) for every \(v \in V({\mathcal {G}})\), then \(\alpha \) is called complete. Note that if \(\alpha \) is partial and \(\alpha (v) \ne \bot \), one can already determine whether v is an \({\mathbb {S}}\) or a \({\mathbb {D}}\) node, and hence we may say that v is a speciation or a duplication under \(\alpha \). Also note that the definitions of \({\hat{d}}(\alpha ), l(\alpha )\) and \(h_{\alpha }(s)\) extend naturally to a partial mapping \(\alpha \) by considering the forest induced by the nodes not mapped to \(\bot \).
If \(\alpha \) is a partial mapping, we call \(\alpha '\) a completion of \(\alpha \) if \(\alpha '\) is complete, and \(\alpha (v) = \alpha '(v)\) whenever \(\alpha (v) \ne \bot \). Note that such a completion always exists, as in particular one can map every \(\bot \)node to the root of S (such a mapping must be valid, since all ancestors of a \(\bot \)node are also \(\bot \)nodes, which ensures that \(r(S) = \alpha '(v) \ge \alpha '(v')\) for every descendant \(v'\) of a newly mapped \(\bot \)node v). We say that \(\alpha '\) is an optimal completion of \(\alpha \) if \(cost^{SD}(\alpha ')\) is minimum among every possible completion of \(\alpha \). For a minimal \(\bot \)node v with children \(v_1\) and \(v_2\), we denote \(\mu _{\alpha }(v) = LCA_S(\alpha (v_1), \alpha (v_2))\), i.e. the lowest species of S to which v can possibly be mapped to in any completion of \(\alpha \). Observe that \(\mu _{\alpha }(v) \ge \mu (v)\). Moreover, if v is a minimal \(\bot \)node, then in any completion \(\alpha '\) of \(\alpha \), \(\alpha '[v \rightarrow \mu _{\alpha }(v)]\) is a valid mapping. A minimal \(\bot \)node v is called a lowest minimal \(\bot \) node if, for every minimal \(\bot \)node w distinct from v, either \(\mu _{\alpha }(v) \le \mu _{\alpha }(w)\) or \(\mu _{\alpha }(v)\) and \(\mu _{\alpha }(w)\) are incomparable.
The following Lemma forms the basis of our FPT algorithm, as it allows us to bound the possible mappings of a minimal \(\bot \)node.
Lemma 6
Let \(\alpha \) be a partial mapping and let v be a minimal \(\bot \) node. Then for any optimal completion \(\alpha ^*\) of \(\alpha \), \(\alpha ^*(v) \le par^{\lceil \delta /\lambda \rceil }(\mu _{\alpha }(v))\).
Proof
Let \(\alpha ^*\) be an optimal completion of \(\alpha \) and let \(\alpha ' := \alpha ^*[v \rightarrow \mu _{\alpha }(v)]\). Note that \({\hat{d}}(\alpha ') \le {\hat{d}}(\alpha ^*) + 1\). Now suppose that \(\alpha ^*(v) > par^{\lceil \delta /\lambda \rceil }(\mu _{\alpha }(v))\). Then by the Shiftdown lemma, \(l(\alpha ^*)  l(\alpha ') > \lceil \delta /\lambda \rceil \ge \delta /\lambda \). Thus \(cost^{SD}(\alpha ^*)  cost^{SD}(\alpha ') > \delta + \lambda (\delta /\lambda ) = 0\). This contradicts the optimality of \(\alpha ^*\). \(\square \)
A node \(v \in V({\mathcal {G}})\) is a required duplication (under \(\alpha \)) if, in any completion \(\alpha '\) of \(\alpha \), v is a duplication under \(\alpha '\). We first show that required duplications are easy to find.
Lemma 7
Let v be a minimal \(\bot \) node under \(\alpha \), and let \(v_1\) and \(v_2\) be its two children. Then v is a required duplication under \(\alpha \) if and only if \(\alpha (v_1) \ge \mu (v)\) or \(\alpha (v_2) \ge \mu (v)\).
Proof
Suppose that \(\alpha (v_1) \ge \mu (v)\), and let \(\alpha '\) be a completion of \(\alpha \). If \(\alpha '(v) = \alpha '(v_1)\), then v is a duplication by definition. Otherwise, \(\alpha '(v) > \alpha '(v_1) = \alpha (v_1) \ge \mu (v)\), and v is a duplication by Lemma 2. The case when \(\alpha (v_2) \ge \mu (v)\) is identical.
Conversely, suppose that \(\alpha (v_1) < \mu (v)\) and \(\alpha (v_2) < \mu (v)\). Then \(\alpha (v_1)\) and \(\alpha (v_2)\) must be incomparable descendants of \(\mu (v)\) (because otherwise if e.g. \(\alpha (v_1) \le \alpha (v_2)\), then we would have \(\mu (v) = LCA_S(\mu (v_1), \mu (v_2)) \le LCA_S(\alpha (v_1), \alpha (v_2)) = \alpha (v_2),\) whereas we are assuming that \(\alpha (v_2) < \mu (v)\)). Take any completion \(\alpha '\) of \(\alpha \) such that \(\alpha '(v) = \mu (v)\). To see that v is a speciation under \(\alpha '\), it remains to argue that \(\alpha '(v) = \mu (v) = LCA_{S}(\alpha (v_1), \alpha (v_2))\). Since \(\mu (v)\) is an ancestor of both \(\alpha (v_1)\) and \(\alpha (v_2)\), we have \(LCA_{S}(\alpha (v_1), \alpha (v_2)) \le \mu (v)\). We also have \(\mu (v) = LCA_{S}(\mu (v_1), \mu (v_2)) \le LCA_{S}(\alpha (v_1), \alpha (v_2))\), and equality follows. \(\square \)
Lemmas 8 and 9 allow us to find minimal \(\bot \)nodes of \({\mathcal {G}}\) that are the easiest to deal with, as their mapping in an optimal completion can be determined with certainty.
Lemma 8
Let v be a minimal \(\bot \) node under \(\alpha \). If v is not a required duplication under \(\alpha \), then \(\alpha ^*(v) = \mu _{\alpha }(v)\) for any optimal completion \(\alpha ^*\) of \(\alpha \).
Proof
Let \(v_1, v_2\) be the children of v, and let \(\alpha ^*\) be an optimal completion of \(\alpha \). Since v is not a required duplication, by Lemma 7 we have \(\alpha (v_1) < \mu (v)\) and \(\alpha (v_2) < \mu (v)\) and, as argued in the proof of Lemma 7, \(\alpha (v_1)\) and \(\alpha (v_2)\) are incomparable. We thus have that \(\mu _{\alpha }(v) = \mu (v)\). Then \(\alpha ^*[v \rightarrow \mu (v)]\) is a valid mapping, and v is a speciation under this mapping. Hence \({\hat{d}}(\alpha ^*[v \rightarrow \mu (v)]) \le {\hat{d}}(\alpha ^*)\). Then by the Shiftdown lemma, this new mapping has fewer losses, and thus attains a lower cost than \(\alpha ^*\). \(\square \)
Lemma 9
Let v be a minimal \(\bot \) node under \(\alpha \), and let \(\alpha _{v} := \alpha [v \rightarrow \mu _{\alpha }(v)]\). If \({\hat{d}}(\alpha ) = {\hat{d}}(\alpha _{v})\), then \(\alpha ^*(v) = \mu _{\alpha }(v)\) for any optimal completion \(\alpha ^*\) of \(\alpha \).
Proof
Let \(\alpha ^*\) be an optimal completion of \(\alpha \). Denote \(s := \mu _{\alpha }(v)\), and assume that \(\alpha ^*(v) > s\) (as otherwise, we are done). Let \(\alpha ' = \alpha ^*[v \rightarrow s]\). We have that \(l(\alpha ') < l(\alpha ^*)\) by the Shiftdown lemma. To prove the Lemma, we then show that \({\hat{d}}(\alpha ') \le {\hat{d}}(\alpha ^*)\). Suppose otherwise that \({\hat{d}}(\alpha ') > {\hat{d}}(\alpha ^*)\). As only v changed mapping to s to go from \(\alpha ^*\) to \(\alpha '\), this implies that \(h_{\alpha '}(s) > h_{\alpha ^*}(s)\) because of v. Since under \(\alpha ^*\), no ancestor of v is mapped to s, it must be that under \(\alpha '\), v is the root of a subtree T of height \(h_{\alpha '}(s)\) of duplications in s. Since T contains only descendants of v, it must also be that \(h_{\alpha _{v}}(s) = h_{\alpha '}(s)\) (here \(\alpha _v\) is the mapping defined in the Lemma statement). As we are assuming that \(h_{\alpha '}(s) > h_{\alpha ^*}(s)\), we get \(h_{\alpha _{v}}(s) > h_{\alpha ^*}(s)\). This is a contradiction, since \(h_{\alpha ^*}(s) \ge h_{\alpha }(s) = h_{\alpha _{v}}(s)\) (the left inequality because \(\alpha ^*\) is a completion of \(\alpha \), and the right equality by the choice of \(\alpha _v\)). Then \(l(\alpha ') < l(\alpha ^*)\) and \({\hat{d}}(\alpha ') \le {\hat{d}}(\alpha ^*)\) contradicts the fact that \(\alpha ^*\) is optimal. \(\square \)
 [C1]:

v is not easy;
 [C2]:

for all duplication nodes w (under \(\alpha \) with \(\alpha (w) \ne \bot \)), either \(\alpha (w) \le \mu _{\alpha }(v)\) or \(\alpha (w)\) is incomparable with \(\mu _{\alpha }(v)\).
Lemma 10
Suppose that \(\alpha \) is a clean partial mapping, and let \(\alpha ^*\) be an optimal completion of \(\alpha \). Let v be a lowest minimal \(\bot \) node under \(\alpha \), and let \(s := \alpha ^*(v)\). Then for every minimal \(\bot \) node w such that \(\mu _{\alpha }(w) \le s\), we have \(\alpha ^*(w) = s\).
Proof
Denote \(\alpha ' := \alpha [v \rightarrow s]\). Suppose first that \(s = \mu _{\alpha }(v)\). Note that since \(\alpha \) is clean, v is not easy, which implies that \(h_{\alpha '}(s) = h_{\alpha }(s) + 1\). Since v is a lowest minimal \(\bot \)node, if w is a minimal \(\bot \)node such that \(\mu _{\alpha }(w) \le s\), we must have \(\mu _{\alpha }(w) = s\), as otherwise v would not have the ‘lowest’ property. Moreover, because v and w are both minimal \(\bot \)nodes under the partial mapping \(\alpha \), one cannot be the ancestor of the other and so v and w are incomparable. This implies that mapping w to s under \(\alpha '\) cannot further increase \(h_{\alpha '}(s)\) (because we already increased it by 1 when mapping v to s). Thus \({\hat{d}}(\alpha ') = {\hat{d}}(\alpha '[w \rightarrow s])\), and w is easy under \(\alpha '\) and must be mapped to s by Lemma 9. This proves the \(\alpha ^*(v) = \mu _{\alpha }(s)\) case.
Now assume that \(s > \mu _{\alpha }(v)\), and let w be a minimal \(\bot \)node with \(\mu _{\alpha }(w) \le s\). Let us denote \(s' := \alpha ^*(w)\). If \(s' = s\), then we are done. Suppose that \(s' < s\), noting that \(h_{\alpha ^*}(s') > 0\) (because w must be a duplication node, due to \(\alpha \) being clean). If \(s' = \mu _{\alpha }(v)\), then w is also a lowest minimal \(\bot \)node. In this case, using the arguments from the previous paragraph and swapping the roles of v and w, one can see that v is easy in \(\alpha [w \rightarrow s']\) and must be mapped to \(s' < s\), a contradiction. Thus assume \(s' > \mu _{\alpha }(v)\). Under \(\alpha ^*\), for each child \(v'\) of v, we have \(\alpha ^*(v') \le \mu _{\alpha }(v) < s'\), and for each ancestor \(v''\) of v, we have \(\alpha ^*(v'') \ge \alpha ^*(v) = s > s'\). Therefore, by remapping v to \(s'\), v is the only duplication mapped to \(s'\) among its ancestors and descendants. In other words, because \(h_{\alpha ^*}(s') > 0\), we have \({\hat{d}}(\alpha ^*[v \rightarrow s']) \le {\hat{d}}(\alpha ^*)\). Moreover by the Shiftdown lemma, \(l(\alpha ^*[v \rightarrow s']) < l(\alpha ^*)\), which contradicts the optimality of \(\alpha ^*\).
The remaining case is \(s' > s\). Note that \(h_{\alpha ^*}(s) > 0\) (because v must be a duplication node, due to \(\alpha \) being clean). Since it holds that v is a minimal \(\bot \)node, that \(\alpha \) is clean and that \(s > \mu _{\alpha }(v)\), it must be the case that \(\alpha \) has no duplication mapped to s (by the second property of cleanness). In particular, w has no descendant that is a duplication mapped to s under \(\alpha \) (and hence under \(\alpha ^*\)). Moreover, as \(s' = \alpha ^*(w) > s\), w has no ancestor that is a duplication mapped to s. Thus \({\hat{d}}(\alpha ^*[w \rightarrow s]) \le {\hat{d}}(\alpha ^*)\), and the Shiftdown lemma contradicts the optimality of \(\alpha ^*\). This concludes the proof. \(\square \)
We are finally ready to describe our algorithm. We start from a partial mapping \(\alpha \) with \(\alpha (v) = \bot \) for every internal node v of \({\mathcal {G}}\). We gradually “fillup” the \(\bot \)nodes of \(\alpha \) in a bottomup fashion, maintaining a clean mapping at each step and ensuring that each decision leads to an optimal completion \(\alpha ^*\). To do this, we pick a lowest minimal \(\bot \)node v, and “guess” \(\alpha ^*(v)\) among the \(\lceil \delta /\lambda \rceil \) possibilities. This increases some \(h_{\alpha }(s)\) by 1. For each such guess s, we use Lemma 10 to map the appropriate minimal \(\bot \)nodes to s, then take care of the easy nodes to obtain another clean mapping. We repeat until we have either found a complete mapping or we have a duplication height higher than d. An illustration of a pass through the algorithm is shown in Fig. 5.
The complexity follows from the fact that the algorithm creates a search tree of degree \(\lceil \delta /\lambda \rceil \) of depth at most d. The main technicality is to show that the algorithm maintains a clean mapping before each recursive call.^{2}
Theorem 4
Algorithm 1 is correct and finds a minimum cost mapping \(\alpha ^*\) satisfying \({\hat{d}}(\alpha ^*) \le d\), if any, in time \(O(\lceil \frac{\delta }{\lambda } \rceil ^d \cdot n \cdot \frac{\delta }{\lambda })\).
Proof
We show by induction over the depth of the search tree that, in any recursive call made to Algorithm 1 with partial mapping \(\alpha \), the algorithm returns the cost of an optimal completion \(\alpha ^*\) of \(\alpha \) having \({\hat{d}}(\alpha ^*) \le d\), or \(\infty \) if no such completion exists—assuming that the algorithm receives a clean mapping \(\alpha \) as input. Thus in order to use induction, we must also show that at each recursive call done on line 15, \(\alpha '\) is a clean mapping. We additionally claim that the search tree created by the algorithm has depth at most d. To show this, we will also prove that every \(\alpha '\) sent to a recursive call satisfies \({\hat{d}}(\alpha ') = {\hat{d}}(\alpha ) + 1\).
The base cases of lines 3–5 are trivial. For the induction step, let v be the lowest minimal \(\bot \)node chosen on line 7. By Lemma 6, if \(\alpha ^*\) is an optimal completion of \(\alpha \) and \(s = \alpha ^*(v)\), then \(\mu _{\alpha }(v) \le s \le par^{\lceil \delta /\lambda \rceil }(\mu _{\alpha }(v))\). We try all the \(\lceil \delta /\lambda \rceil \) possibilities in the forloop on line 9. The forloop on line 11 is justified by Lemma 10, and the forloop on line 13 is justified by Lemmas 8 and 9. Assuming that \(\alpha '\) is clean on line 15, by induction the recursive call will return the cost of an optimal completion \(\alpha ^*\) of \(\alpha '\) having \({\hat{d}}(\alpha ^*) \le d\), if any such completion exists. It remains to argue that for every \(\alpha '\) sent to a recursive call on line 15, \(\alpha '\) is clean and \({\hat{d}}(\alpha ') = {\hat{d}}(\alpha ) + 1\) .
Let us first show that such a \(\alpha '\) is clean, for each choice of s on line 9. There clearly cannot be an easy node under \(\alpha '\) after line 13, so we must show C2, i.e. that for any minimal \(\bot \)node w under \(\alpha '\), there is no duplication z under \(\alpha '\) satisfying \(\alpha '(z) > \mu _{\alpha '}(w)\). Suppose instead that \(\alpha '(z) > \mu _{\alpha '}(w)\) for some duplication node z. Let \(w_0\) be a descendant of w that is a minimal \(\bot \)node in \(\alpha \) (note that \(w_0 = w\) is possible). We must have \(\mu _{\alpha '}(w) \ge \mu _{\alpha }(w_0)\). By our assumption, we then have \(\alpha '(z) > \mu _{\alpha '}(w) \ge \mu _{\alpha }(w_0)\). Then z cannot be a duplication under \(\alpha \), as otherwise \(\alpha \) itself could not be clean (by C2 applied on z and \(w_0\)). Thus z is a newly introduced duplication in \(\alpha '\), and so z was a \(\bot \)node under \(\alpha \). Note that Algorithm 1 maps \(\bot \)nodes of \({\mathcal {G}}\) one after another, in some order \((z_1, z_2, \ldots , z_k)\). Suppose without loss of generality that z is the first duplication node in this ordering that gets mapped to \(\alpha '(z)\). There are two cases: either \(\alpha '(z) \ne s\), or \(\alpha '(z) = s\).
Suppose first that \(\alpha '(z) \ne s\). Lines 10 and 11 can only map \(\bot \)nodes to s, and line 13 either maps speciation nodes, or easy nodes that become duplications. Thus when \(\alpha '(z) \ne s\), we may assume that z falls into the latter case, i.e. z is easy before being mapped, so that mapping z to \(\alpha '(z)\) does not increase \(h_{\alpha '}(\alpha '(z))\). Because z is the first \(\bot \)node that gets mapped to \(\alpha '(z)\), this is only possible if there was already a duplication \(z_0\) mapped to \(\alpha '(z)\) in \(\alpha \). This implies that \(\alpha (z_0) = \alpha '(z) > \mu _{\alpha }(w_0)\), and that \(\alpha \) was not clean (by C2 applied on \(z_0\) and \(w_0\)). This is a contradiction.
We may thus assume that \(\alpha '(z) = s\). This implies \(\mu _{\alpha '}(w) < \alpha '(z) = s\). If w was a minimal \(\bot \)node in \(\alpha \), it would have been mapped to s on line 11, and so in this case w cannot also be a minimal \(\bot \)node in \(\alpha '\), as we supposed. If instead w was not a minimal \(\bot \)node in \(\alpha \), then w has a descendant \(w_0\) that was a minimal \(\bot \)node under \(\alpha \). We have \(\mu _{\alpha }(w_0) \le \mu _{\alpha '}(w) < s\), which implies that \(w_0\) gets mapped to s on line 11. This makes \(\mu _{\alpha '}(w) < s\) impossible, and we have reached a contradiction. We deduce that z cannot exist, and that \(\alpha '\) is clean.
It remains to show that \({\hat{d}}(\alpha ') = {\hat{d}}(\alpha ) + 1\). Again, let s be the chosen species on line 9. Suppose first that \(s = \mu _{\alpha }(v)\). Then \(h_{\alpha [v \rightarrow s]}(s) = h_{\alpha }(s) + 1\), as otherwise v would be easy under \(\alpha \), contradicting its cleanness. In this situation, as argued in the proof of Lemma 10, each node w that gets mapped to s on line 11 or on line 13 is easy, and thus cannot further increase the height of the duplications in s. If \(s > \mu _{\alpha }(v)\), then \(h_{\alpha [v \rightarrow s]}(s) = 1 = h_{\alpha }(s) + 1\), since by cleanness no duplication under \(\alpha \) maps to s. Here, each node w that gets mapped on line 11 has no descendant nor ancestor mapped to s, and thus the height does not increase. Noting that remapping easy nodes on line 13 cannot alter the duplication heights, we get in both cases that \({\hat{d}}(\alpha [v \rightarrow s]) = {\hat{d}}(\alpha ) + 1\). This proves the correctness of the algorithm.
As for the complexity, the algorithm creates a search tree of degree \(\lceil \delta /\lambda \rceil \) and of depth at most d. Each pass can easily seen to be feasible in time \(O(\delta /\lambda \cdot n)\) (with appropriate preparsing to compute \(\mu _{\alpha }(v)\) in constant time, and to decide if a node is easy or not in constant time as well), and so the total complexity is \(O(\lceil \delta /\lambda \rceil ^d n \cdot \frac{\delta }{\lambda })\). \(\square \)
Experiments
We also reanalyzed the data set of yeast species described in [2]. First, we selected from the data set the 2379 gene trees containing all 16 species and refined unsupported branches using the method described in [16] and implemented in ecceTERA [15] with a bootstrap threshold of 0.9 and \(\delta =\lambda =1\). Using our method with \(\delta =1.5\), \(\lambda =1\) we were able to detect the ancient genome duplication in Saccharomyces cerevisiae already established using synteny information [17], with 216 gene families supporting the event. Other nodes with a signature of segmental duplication are nodes 7, 6, 13 and 2 (refer to Fig. 7) with respectively 190, 157, 148 and 136 gene families supporting the event. It would be interesting to see if the synteny information supports these hypotheses.
Note on the costs: As for any other parsimony method, the costs associated to the events in Algorithm 1 (i.e. \(\delta \) and \(\lambda \)) have to be fixed by the user and cannot be estimated. Several possibilities exist to estimate these costs. A method for cost estimation for \({{\mathbb {D}}}{{\mathbb {L}}}\) reconciliation, based on reducing large fluctuations in ancient genome sizes, has been proposed in [6]. Another possibility is to use costs estimated by a maximum likelihood method, as done in [23], where the costs estimated in [27, 28] were used. Uncertainty in the input costs can also be tackled by considering several cost vectors and paretooptimal reconciliations, or even via sampling; examples of these approach can be found respectively in [29] in the context of \(\mathbb {DTL}\) reconciliations (i.e. \({{\mathbb {D}}}{{\mathbb {L}}}\) reconciliations where horizontal gene transfers are also allowed) and in [4], where Boltzmann sampling is used in the context of evolution of gene adjacencies, a problem related to reconciliations.
Conclusion
We have presented an approach for the reconciliation of a set of gene trees and a species tree, based on segmental macroevolutionary events, where segmental duplication events and losses are associated with cost \(\delta \) and \(\lambda \), respectively. We have shown that the problem is polynomialtime solvable when \(\delta \le \lambda \), since LCAmapping is already an optimal solution. When \(\delta > \lambda \) the problem is NPhard, even when \(\lambda = 0\) and a single gene tree is given. This result solves a long standing open problem on the complexity of the reconciliation of a set of gene trees with a species tree. Moreover, we have given a fixedparameter algorithm of time complexity \(O\left(\lceil \frac{\delta }{\lambda } \rceil ^{d} \cdot n \cdot \frac{\delta }{\lambda }\right)\), where d is the number of segmental duplications, that has been tested on real data, showing its effectiveness.
This work poses a variety of questions that deserve further investigation. The complexity of the problem when \(\delta /\lambda \) is a constant remains an open problem. Moreover, our FPT algorithm can handle data sets with a sum of duplication height of about \(d = 30\), but in the future, one might consider whether there exist fast approximation algorithms for MPRSTSD in order to attain better scalability. Other future directions include a multivariate complexity analysis of the problem, in order to understand whether it is possible to identify other parameters that are small in practice. Finally, we plan to extend the experimental analysis to other data sets, for instance for the detection of whole genome duplications in plants.
To our knowledge, this is the first publicly available reconciliation software for segmental duplications.
There is a subtlety to consider here. What we have shown is that if there exists a mapping \(\alpha \) of minimum cost \(cost^{SD}({\mathcal {G}}, S)\) with \({\hat{d}}(\alpha ) \le d\), then the algorithm finds it. It might be that a reconciliation \(\alpha \) satisfying \({\hat{d}}(\alpha ) \le d\) exists, but that the algorithm returns no solution. This can happen in the case that \(\alpha \) is not of cost \(cost^{SD}({\mathcal {G}}, S)\).
Note that for this data set we used a high value for \(\frac{\delta }{\lambda }\) since, because of the sampling strategy, we expect that all relevant genes have been sampled (recall that in ExactMGD, \(\lambda \) is implicitly set to 0).
Declarations
Authors' contributions
RD, ML and CS all participated in writing the manuscript and establishing the theoretical results. ML implemented the algorithm, CS ran the experiments. All authors read and approved the final manuscript.
Acknowledgements
The authors would like to thank Mukul Bansal for providing the eukaryotes data set for the experiments section.
Competing interests
The authors declare that they have no competing interests.
Availability of data and materials
The implementation is available at https://github.com/AEVOlab/MultRec. The data used in the experiments is available on demand.
Consent for publication
Not applicable.
Ethics approval and consent to participate
Not applicable.
Funding
The work of ML was supported by a postdoctoral fellowship from the Natural Sciences and Engineering Research Council (NSERC), Canada.
Publisher’s Note
Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.
Open AccessThis article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.
Authors’ Affiliations
References
 Bansal MS, Eulenstein O. The multiple gene duplication problem revisited. Bioinformatics. 2008;24(13):132–8.View ArticleGoogle Scholar
 Butler G, Rasmussen MD, Lin MF, Santos Manuel AS, Sakthikumar S, Munro CA, Rheinbay E, Grabherr M, Forche A, Reedy JL, et al. Evolution of pathogenicity and sexual reproduction in eight candida genomes. Nature. 2009;459(7247):657.View ArticleGoogle Scholar
 Cedric C, Nadia EM. New perspectives on gene family evolution: losses in reconciliation and a link with supertrees. In: Annual international conference on research in computational molecular biology. Springer; 2009. p. 46–58.Google Scholar
 Chauve C, Ponty Y, Zanetti JPP. Evolution of genes neighborhood within reconciled phylogenies: an ensemble approach. In: Brazilian symposium on bioinformatics. Springer; 2014. p. 49–56.Google Scholar
 Chauve C, Rafiey A, Davin AA, Scornavacca C, Veber P, Boussau B, Szollosi G, Daubin V, Tannier E. Maxtic: fast ranking of a phylogenetic tree by maximum time consistency with lateral gene transfers. bioRxiv, 2017. https://www.biorxiv.org/content/early/2017/11/07/127548.
 David LA, Alm EJ. Rapid evolutionary innovation during an archaean genetic expansion. Nature. 2011;469(7328):93.View ArticleGoogle Scholar
 Davín AA, Tannier E, Williams TA, Boussau B, Daubin V, Szöllősi GJ. Gene transfers can date the tree of life. Nat Ecol Evol. 2018;2(5):904.View ArticleGoogle Scholar
 Delabre M, ElMabrouk N, Huber KT, Lafond M, Moulton V, Noutahi E, Castellanos MS. Reconstructing the history of syntenies through superreconciliation. In: RECOMB international conference on comparative genomics. Springer; 2018. 179–195.Google Scholar
 Dondi R, Lafond M, ElMabrouk N. Approximating the correction of weighted and unweighted orthology and paralogy relations. Algorith Mol Biol. 2017;12(1):1–4. https://doi.org/10.1186/s130150170096x.View ArticleGoogle Scholar
 Duchemin W. Phylogeny of dependencies and dependencies of phylogenies in genes and genomes. Ph.D. thesis, Universit de Lyon, 2017.Google Scholar
 Duchemin W, Anselmetti Y, Patterson M, Ponty Y, Bérard S, Chauve C, Scornavacca C, Daubin V, Tannier E. Decostar: reconstructing the ancestral organization of genes or genomes using reconciled phylogenies. Genome Biol Evol. 2017;9(5):1312–9.View ArticleGoogle Scholar
 Fellows M, Hallett M, Stege U. On the multiple gene duplication problem. In: International symposium on algorithms and computation. Springer; 1998. p. 348–357.Google Scholar
 Goodman M, Czelusniak J, Moore GW, RomeroHerrera AE, Matsuda G. Fitting the gene lineage into its species lineage, a parsimony strategy illustrated by cladograms constructed from globin sequences. Syst Biol. 1979;28(2):132–63.View ArticleGoogle Scholar
 Guigo R, Muchnik I, Smith TF. Reconstruction of ancient molecular phylogeny. Mol Phylogenet Evol. 1996;6(2):189–213.View ArticleGoogle Scholar
 Jacox E, Chauve C, Szllsi GJ, Ponty Y, Scornavacca C. Eccetera: comprehensive gene treespecies tree reconciliation using parsimony. Bioinformatics. 2016;32(13):2056–8. https://doi.org/10.1093/bioinformatics/btw105.View ArticlePubMedGoogle Scholar
 Jacox E, Weller M, Tannier E, Scornavacca C. Resolution and reconciliation of nonbinary gene trees with transfers, duplications and losses. Bioinformatics. 2017;33(7):980–7.PubMedGoogle Scholar
 Kellis M, Birren BW, Lander ES. Proof and evolutionary analysis of ancient genome duplication in the yeast saccharomyces cerevisiae. Nature. 2004;428(6983):617.View ArticleGoogle Scholar
 Lafond M, Dondi R, ElMabrouk N. The link between orthology relations and gene trees: a correction perspective. Algorith Mol Biol. 2016;11:4. https://doi.org/10.1186/s1301501600677.View ArticleGoogle Scholar
 Lafond M, Miardan MM, Sankoff D. Accurate prediction of orthologs in the presence of divergence after duplication. Bioinformatics. 2018;34(13):i366–75.View ArticleGoogle Scholar
 Luo CW, Chen MC, Chen YC, Yang RW, Liu HF, Chao KM. Lineartime algorithms for the multiple gene duplication problems. IEEE/ACM Trans Comput Biol Bioinf. 2011;8(1):260–5.View ArticleGoogle Scholar
 Ma B, Li M, Zhang L. From gene trees to species trees. SIAM J Comput. 2000;30(3):729–52.View ArticleGoogle Scholar
 Maddison WP. Gene trees in species trees. Syst Biol. 1997;46(3):523–36.View ArticleGoogle Scholar
 Nguyen TH, Ranwez V, Pointet S, Chifolleau AMA, Doyon JP, Berry V. Reconciliation and local gene tree rearrangement can be of mutual profit. Algorith Mol Biol. 2013;8(1):12.View ArticleGoogle Scholar
 Page RDM. Maps between trees and cladistic analysis of historical associations among genes, organisms, and areas. Syst Biol. 1994;43(1):58–77.Google Scholar
 Page RDM, Cotton JA. Vertebrate phylogenomics: reconciled trees and gene duplications. Pac Symp Biocomput. 2002;7:536–47.Google Scholar
 Paszek J, Gorecki P. Efficient algorithms for genomic duplication models. IEEE/ACM Trans Comput Biol Bioinf. 2017;15(5):1515–24.Google Scholar
 Szöllősi GJ, Boussau B, Abby SS, Tannier E, Daubin V. Phylogenetic modeling of lateral gene transfer reconstructs the pattern and relative timing of speciations. In: Proceedings of the national academy of sciences, 2012. p. 201202997.Google Scholar
 Szöllősi GJ, Daubin V. Modeling gene family evolution and reconciling phylogenetic discord. In: Evolutionary genomics. Springer; 2012. p. 29–51.Google Scholar
 To TH, Jacox E, Ranwez V, Scornavacca C. A fast method for calculating reliable event supports in tree reconciliations via pareto optimality. BMC Bioinf. 2015;16(1):384.View ArticleGoogle Scholar
 Ullah I, Sjöstrand J, Andersson P, Sennblad B, Lagergren J. Integrating sequence evolution into probabilistic orthology analysis. Syst Biol. 2015;64(6):969–82.View ArticleGoogle Scholar