Christian Baudet

Béatrice Donati

Blerina Sinaimeri

Pierluigi Crescenzi

Christian Gautier

Catherine Matias

Marie-France Sagot
email: marie-france.sagot@inria.fr

M.-F Sagot

Marie-France Sagot

Cophylogeny Reconstruction via an Approximate Bayesian Computation

ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

event-based methods then assign a cost to each type of event in order to find a reconstruction of minimum cost. The main problem with this approach is that the cost of the events strongly influences the reconciliation obtained. Some earlier approaches attempt to avoid this problem by finding a Pareto set of solutions and hence by considering event costs under some minimisation constraints.

To deal with this problem, we developed an algorithm, called Coala, for estimating the frequency of the events based on an approximate Bayesian computation approach. The benefits of this method are twofold: (1) it provides more confidence in the set of costs to be used in a reconciliation, and (2) it allows estimation of the frequency of the events in cases where the dataset consists of trees with a large number of taxa.

We evaluate our method on simulated and on biological datasets. We show that in both cases, for the same pair of host and parasite trees, different sets of frequencies for the events lead to equally probable solutions. Moreover, often these solutions differ greatly in terms of the number of inferred events. It appears crucial to take this into account before attempting any further biological interpretation of such reconciliations. More generally, we also show that the set of frequencies can vary widely depending on the input host and parasite trees. Indiscriminately applying a standard vector of costs may thus not be a good strategy.

(Keywords: cophylogeny, host/parasite systems, likelihood-free inference, approximate Bayesian computation.)

Cophylogeny is the reconstruction of ancient relationships among ecologically linked groups of organisms from their phylogenetic information. The study of host-parasite systems has a long history and has been already well addressed in the literature (Page 1994b;[START_REF] Huelsenbeck | Statistical tests of host-parasite cospeciation[END_REF][START_REF] Charleston | Jungles: a new solution to the host/parasite phylogeny reconciliation problem[END_REF][START_REF] Paterson | Analytical approaches to measuring cospeciation of host and parasites: Through a glass, darkly[END_REF][START_REF] Merkle | Reconstruction of the cophylogenetic history of related phylogenetic trees with divergence timing information[END_REF]Conow et al. 2010, for example). It also has broad applications throughout biology. For instance, the same mathematical model can be applied to gene-species associations [START_REF] Hallett | Efficient algorithms for lateral gene transfer problems[END_REF]Doyon et al. 2011a,b;[START_REF] Tofigh | Simultaneous identification of duplications and lateral gene transfers[END_REF][START_REF] Bansal | Efficient algorithms for the reconciliation problem with gene duplication, horizontal transfer and loss[END_REF]. Hence, any single method for host/parasite associations that is developed could be applicable to both situations. Lately indeed, there have been attempts to introduce a general framework that incorporates all existing models [START_REF] Wieseke | Unifying parsimonious tree reconciliation[END_REF]).

Our work is particularly focused on reconstructing the coevolutionary history of host-parasite systems. Specifically, we are given a host tree H, a parasite tree P , and a function ϕ mapping the leaves of P to the leaves of H. In general, four main macro-evolutionary events are assumed to be recovered: (a) cospeciation, when the parasite diverges in correspondence to the divergence of a host species; (b) duplication, when the parasite diverges "without the stimulus of host speciation" [START_REF] Paterson | Analytical approaches to measuring cospeciation of host and parasites: Through a glass, darkly[END_REF]; (c) host-switching, when a parasite switches, or jumps from one host species to another independent of any host divergence; and (d) loss, which can describe three different and undistinguishable situations: (i) speciation of the host species independently of the parasite, which then follows just one of the new host species due to factors such as, for instance, geographical isolation; (ii) cospeciation of host and parasite, followed by extinction of one of the new parasite species and; (iii) same as (ii) with failure to detect the parasite in one of the two new host species. These events are depicted in Figure 1.

A parsimonious solution for reconciling the phylogenetic trees for hosts on one side, and parasites on the other, simply assigns a cost to each of the four types of events and then seeks to minimise the total cost of the mapping. If host switches are forbidden, exact solutions can be found in time linear in the size of the trees [START_REF] Goodman | Fitting the gene lineage into its species lineage, a parsimony strategy illustrated by cladograms constructed from globin sequences[END_REF]Page 1994a;[START_REF] Mirkin | A biologically consistent model for comparing molecular phylogenies[END_REF]Guigó et al. 1996, for example). If timing information is available, e.g. if we happen to know the order in which speciation events occurred in the host phylogeny, then any proposed reconciliation must also respect the temporal constraints imposed by the available timing information. Host switches are thus restricted to occur only between co-existing species. When co-existence relationships are known for all host species, the reconciliation problem can again easily be solved using dynamic programming, this time polynomially in the size of the trees [START_REF] Libeskind-Hadas | On the computational complexity of the reticulate cophylogeny reconstruction problem[END_REF][START_REF] Conow | Jane: a new tool for the cophylogeny reconstruction problem[END_REF][START_REF] Drinkwater | An improved node mapping algorithm for the cophylogeny reconstruction problem[END_REF]. However, when timing information is not available, the difficulty of separating between compatible and incompatible switches makes the reconciliation problem NP-hard [START_REF] Ovadia | The cophylogeny reconstruction problem is NP-complete[END_REF][START_REF] Tofigh | Simultaneous identification of duplications and lateral gene transfers[END_REF]. A number of algorithms have been developed that allow for solutions that are biologically unfeasible, that is, solutions where some of the switches induce a contradictory timing ordering for the internal vertices of the host tree (Doyon et al. 2011c).

In this case, the algorithms are able to generate optimal solutions in polynomial time. For the fastest existing ones, see for example [START_REF] Bansal | Efficient algorithms for the reconciliation problem with gene duplication, horizontal transfer and loss[END_REF].

Clearly in all situations, the choice for the cost values is crucial in the solution(s) found. Indeed, arbitrarily choosing a cost vector may lead to solutions where the events in the optimal solutions do not necessary reflect the reality [START_REF] Charleston | Recent results in cophylogeny mapping[END_REF], for example, describes a study on the distribution of the events in optimal reconciliations). From a biological point of view, reasonable cost values for an event-based reconciliation are not easily chosen. It is also natural to think that the frequency of the events is not constant across datasets. Thus, different pairs of host/parasite phylogenies might be associated with different cost events. Moreover, our results show that for the same pair of host and parasite trees, different reconciliations -in the sense of presenting a different set of frequencies for the events -may constitute equally probable solutions. It is thus crucial to take this into account before attempting any further biological interpretation of such reconciliations. Some approaches [START_REF] Charleston | Treemap 3b[END_REF][START_REF] Libeskind-Hadas | Pareto-optimal phylogenetic tree reconciliation[END_REF] attempt to choose the costs of the events by adopting some minimisation constraints and by focusing on Pareto optimal solutions. As indicated in [START_REF] Ronquist | Tangled trees: phylogeny, cospeciation, and coevolution chap. Parsimony analysis of coevolving species associations[END_REF], if each event is associated with a cost that is inversely related to its likelihood (the more likely is the event, the smaller its cost), then the most parsimonious reconstruction will also, in some sense, be the most likely explanation of the observed data. Likelihood-based approaches should in general be preferred to parsimony-based methods as they remove the subjective step of cost parameter choice and rely instead on a simultaneous inference of parameter values and events. Some work has been done along these lines, for instance in testing for coevolution [START_REF] Huelsenbeck | Statistical tests of host-parasite cospeciation[END_REF][START_REF] Huelsenbeck | A Bayesian framework for the analysis of cospeciation[END_REF]. This however excluded duplications and tended to over-estimate the number of host switches. Instead, in [START_REF] Szöllősi | Efficient exploration of the space of reconciled gene trees[END_REF] all four types of events are considered, but the method was developed with the objective of reconstructing a species tree starting from multiple gene trees. The aim is similar in [START_REF] Arvestad | Bayesian gene/species tree reconciliation and orthology analysis using MCMC[END_REF] but the type of approach is different and the model again incomplete as in [START_REF] Huelsenbeck | Statistical tests of host-parasite cospeciation[END_REF][START_REF] Huelsenbeck | A Bayesian framework for the analysis of cospeciation[END_REF], this time not allowing for host switches. The likelihood approach adopted in [START_REF] Huelsenbeck | Statistical tests of host-parasite cospeciation[END_REF][START_REF] Huelsenbeck | A Bayesian framework for the analysis of cospeciation[END_REF][START_REF] Szöllősi | Efficient exploration of the space of reconciled gene trees[END_REF]) moreover presents the inconvenience of being computationally intensive.

The huge space of possible solutions is also an issue, for instance, in population genetics for reconstructing the evolutionary history of a set of individuals. Since the early work of [START_REF] Pritchard | Population growth of human Y chromosomes: a study of Y chromosome microsatellites[END_REF], the literature from this domain has seen classical Monte Carlo methods and their variants being replaced by Approximate Bayesian Computation (ABC), a set of more efficient statistical techniques [START_REF] Beaumont | Approximate Bayesian computation in population genetics[END_REF]. In complex models, likelihood calculation is often unfeasible or computationally prohibitive. ABC methods, also called likelihood-free inference methods, bypass this issue while remaining statistically well-founded. For more details, we refer to the review of [START_REF] Marin | Approximate Bayesian computational methods[END_REF] as well as the convergence results in [START_REF] Fearnhead | Constructing summary statistics for approximate Bayesian computation: semi-automatic approximate Bayesian computation[END_REF].

Following these ideas, we developed an algorithm, called Coala (Coala stands for "COevolution Assessment by a Likelihood-free Approach", and is also the Portuguese spelling for Koala, the arboreal herbivorous marsupial native to Australia), for estimating the frequency of the events based on a likelihood-free approach. Given a pair of "known" host and parasite trees and a prior probability distribution associated with the events, Coala simulates the temporal evolution of a set of species (the parasites) following the evolution of another set (the hosts) as represented by the latter's known phylogenetic tree.

In this way, it generates under different parameter values a number of simulated multi-labelled parasite trees which are then compared to the known parasite tree. The ABC principle is to keep the parameter values (event probabilities) giving rise to parasite trees that are "close" to the known one. The output of the algorithm is then a distribution on such parameter values that is a surrogate of the posterior probability for the events which would best explain the observed data.

To the best of our knowledge, the only other method that might be compared to ours is the parameter adaptive approach CoRe-Pa [START_REF] Merkle | A parameter-adaptive dynamic programming approach for inferring cophylogenies[END_REF]. In this case, the space of cost vectors is explored either by sampling such vectors at random assuming a uniform distribution model or by using a more sophisticated approach, the so-called Nelder-Mead simplex method [START_REF] Nelder | A simplex method for function minimization[END_REF]. The first appears to be the option by default in CoRe-Pa. In both cases, the function to minimise is the difference between the probabilities directly computed from the cost vector chosen and the actual relative frequencies observed during the reconstruction using such vector. This choice may appear somewhat circular as one would expect that, since reconstruction is driven by the cost vector, the frequency of the events thus reconstructed not only would, but indeed should agree with it.

Method

General framework

The method we propose relies on an approximate Bayesian computation (ABC).

This belongs to a family of likelihood-free Bayesian inference algorithms that attempt to estimate posterior densities for problems where the likelihood is unknown a priori. Given a set of observed data D 0 and starting with a prior distribution π on the parameter space Θ of the model, the objective is to estimate the parameter values θ ∈ Θ that could lead to the observed data using a Bayesian framework. More precisely, the Bayesian paradigm consists in finding the posterior given D 0 defined as:

p(θ|D 0 ) = p(D 0 |θ)π(θ) p(D 0 ) .
If the likelihood function p(D 0 |θ) cannot be derived, then a likelihood-free approximation can be used to estimate this posterior distribution and thus the parameter values. In general, a likelihood-free computation involves a chain of parameter proposals and only accepts a set of parameter values on condition that the model with these values generates data that satisfy a performance criterion with respect to the observed data [START_REF] Sisson | Sequential Monte Carlo without likelihoods[END_REF][START_REF] Sisson | Sequential Monte Carlo without likelihoods[END_REF]. Strict acceptance (or inversely rejection) is based on whether the generated data D S perfectly matches the observed data D 0 . In cases where the probability of perfectly matching the data is very small, a tolerance d(D s , D 0 ) ≤ is adopted to relax the rejection policy, where d is a distance measure. In either case, this is called the fitting criterion. Note that this fitting criterion often relies only on a summary statistic instead of the full datasets D S and D 0 . Moreover, for complex models where the prior and posterior densities are believed to be sufficiently different, the acceptance rate is very low and then the use of a likelihood-free Sequential Monte Carlo (SMC) search that involves many iterations leads to a more appropriate strategy. SMC is also preferred among other possible methods as it is flexible, easy to implement, parallelisable and applicable to general settings [START_REF] Del Moral | An adaptive sequential Monte Carlo method for approximate Bayesian computation[END_REF].

The ABC-SMC algorithms approximate the posterior distribution by using a large set of randomly chosen parameter values. Over sufficiently many iterations and under suitable conditions, the stationary distribution of the Markov chain will approach the distribution of p(θ|d(D S , D 0 ) ≤ ), which will converge to the posterior density p(θ|D 0 ) if the statistics used to compare the generated data with the real one are sufficient and is small enough. In our case, the observed data are a pair of host and parasite trees, denoted by H and P respectively, and a list of associations between parasite and host leaves. The parameter vector of the model is composed of the probabilities of each one of four events corresponding to respectively: speciation of the parasite together with a speciation of its host (called cospeciation); speciation of the parasite without concomittant speciation of the host (called duplication); switch (also known as jump) of the parasite to another host (called host switch, which is further assumed to be without loss on the original host); and speciation of the host without concomitant speciation of the parasite, and thus loss of the parasite for one of the new host species (called loss). We thus have that θ stands for a vector of four probabilities p c , p d , p s , p l . Note that each node in the host tree either matches a node in the parasite tree or represents a loss, giving rise to the four possible events. For this reason, the parameter θ is constrained such that p c + p d + p s + p l = 1 (see Section "Parasite tree generation algorithm" for more details).

Starting from the host and respecting the probabilities of the events specified in a given parameter vector θ i , we generate M parasite trees, where M ≥ 1.

Once a parasite tree P is thus simulated, it can be compared to the real parasite tree P by computing a distance between the two. For a given parameter vector θ i , we can then produce a distance summary of the generated trees, and use this as a criterion in the ABC rejection method. The latter selects the parameter vector(s) that approximate the observed data within a given tolerance threshold.

The ABC-SMC procedure allows us to refine the list of accepted probability vectors by sampling a vector θ i , introducing a small perturbation to it to produce a vector θ i , and then collecting a new distance summary for θ i .

The list of vectors output in the final step of the algorithm defines the posterior distribution of the coevolutionary event probabilities for the given pair H and P . Table 1 shows a summary of the notation used throughout this work.

Parasite tree generation algorithm

The Duplication-Transfer-Loss (DTL) model.-To simulate the coevolutionary history of the two input phylogenies, we rely on the event-based model presented in [START_REF] Tofigh | Simultaneous identification of duplications and lateral gene transfers[END_REF], and later further analysed in [START_REF] Bansal | Efficient algorithms for the reconciliation problem with gene duplication, horizontal transfer and loss[END_REF].

A rooted phylogenetic tree is a leaf-labelled tree that models the evolution of a set of taxa from their most recent common ancestor (placed at the root). The internal vertices of the tree correspond to the speciation events. The tree is rooted so a direction is intrinsically assumed that corresponds to the direction of increasing evolutionary time.

Henceforth, by a phylogenetic tree T , we mean a rooted tree with labelled leaves where every vertex has in-degree 1 and out-degree 2 except for the leaves, which have out-degree 0. For such a tree T , the set of vertices is denoted by V (T ), the set of arcs by A(T ), and the set of leaves by L(T ). The root of T is denoted by r(T ). Given an arc a = (v, w) ∈ A(T ), going from v to w, we call its head, denoted by h(a), the vertex w and its tail, denoted by t(a), the vertex v. For a vertex v ∈ V (T ), we define the set of descendants of v, denoted by Des(v), as the set of vertices in the subtree of T rooted at v (including v). Similarly, the set of ancestors of v, denoted by Anc(v), is the set of vertices in the unique path from the root of T to v (including the end points). For a vertex v ∈ V (T ) different from the root, we call its parent, denoted by par(v), the vertex x for which there is the arc (x, v) ∈ A(T ). We denote by mrca(v, w) the most recent common ancestor of v, w in T . Finally, we denote by ≥ the partial order induced by the ancestry relation in the tree. Formally, for x, y ∈ V (T ), we say that x ≥ y if x ∈ Anc(y). If neither

x ∈ Anc(y) nor y ∈ Anc(x), the vertices are said to be incomparable.

Let H, P be the phylogenetic trees for the host and parasite species respectively. We define ϕ as a function from the leaves of P to the leaves of H that represents the association between currently living host species and parasites. These associations are part of the input of our algorithm, together with the trees themselves. In our model, we allow each parasite to be related to one and only one host, while a host can be related to zero, one, or more than one parasite. More formally, ϕ is thus a function which needs not be surjective nor injective.

A reconciliation γ is a function γ : V (P ) → V (H) that is an extension of ϕ. In particular γ partitions the set V (P ) into three sets Σ, ∆, and Γ which correspond to the vertices of P associated with, respectively, cospeciations, duplications, and host switches.

The reconciliation γ also defines a subset Ξ of A(P ) which corresponds to the arcs associated with host switches.

Given a reconciliation γ, the following holds [START_REF] Tofigh | Simultaneous identification of duplications and lateral gene transfers[END_REF][START_REF] Charleston | Biological Evolution and Statistical Physics[END_REF]:

1. For any p ∈ L(P ), γ(p) = ϕ(p) (γ extends ϕ).
2. For any internal vertex p ∈ V (P ) -L(P ) with children p 1 and p 2 :

(a) mrca(γ(p), γ(p i )) ≥ γ(p i ), for i = 1, 2 (a child cannot be mapped to an ancestor of the parent).

(b) mrca(γ(p), γ(p 1 )) = γ(p) or mrca(γ(p), γ(p 2 ))) = γ(p) (one of the two children is mapped to the subtree rooted at the parent).

3. For any (p 1 , p 2 ) ∈ Ξ ⇔ mrca(γ(p 1 ), γ(p 2 )) ∈ {γ(p 1 ), γ(p 2 )} (the arc (p 1 , p 2 ) is an arc denoting a host switch).

4. For any p ∈ V (P ) -L(P ) with children p 1 and p 2 :

(a) p ∈ Γ ⇔ (p, p 1 ) ∈ Ξ or (p, p 2 ) ∈ Ξ (p is associated with a host switch). (b) p ∈ ∆ ⇔ mrca(γ(p 1 ), γ(p 2 )) ∈ {γ(p 1 ), γ(p 2
)} (the children are mapped to comparable vertices and p is associated with a duplication event).

(c) p ∈ Σ ⇔ mrca(γ(p 1 ), γ(p 2 )) = γ(p) and γ(p 1 ) and γ(p 2 ) are incomparable (p is associated with a cospeciation event).

The losses are identified by a multi-set (generalisation of a set where the elements are allowed to appear more than once) Λ whose elements are in V (H) containing all the vertices h ∈ V (H) that are in the path between the image of a vertex p ∈ V (P ) and the image of one of its children. The images themselves are not included in the count, except for the duplication event, where one of the images is included.

The triple S = H, P, γ is said to be a reconciliation. Given a vector c c , c d , c s , c l of non-negative real values that correspond to the cost of each type of event, the cost of a reconciliation is equal to

c c |Σ| + c d |∆| + c s |Γ| + c l |Λ|.
Finally, a reconciliation is said to be acyclic or time feasible if there exists a total order on V (H) ∪ V (P ) that is consistent with the two partial orders induced by H and P and respects all temporal constraints imposed by both tree topologies and by the set of host switch events. For a detailed definition of a time-feasible scenario, we refer to [START_REF] Stolzer | Inferring duplications, losses, transfers and incomplete lineage sorting with nonbinary species trees[END_REF].

Evolution of parasites.-The evolution of the parasites is simulated by following the evolution of the hosts traversing the phylogenetic tree H from the root to the leaves, and progressively constructing the phylogenetic tree for the parasites. During this process, a single parasite vertex can be in two different states: mapped or unmapped. At the moment of its creation, a new vertex v is unmapped and is assigned a temporary position on an arc a of the host tree H. We denote this position by v, a . From this position, we can decide to map v to a vertex w of H (all coevolutionary events except for loss), or, in the case of a loss, to move v to another position. In the first case, v is always mapped to the vertex h(a)

that is the head of the arc a. We denote this mapping by [v : w] with w = h(a).

Since in all three non-loss cases (cospeciation, duplication, and host switch), the parasite is supposed to speciate and two children are created for v, denoted by v 1 and v 2 .

Their positioning along arcs of the host then depends on which of the three events took place. In the case of a loss, no child for v is created (at this step) since there is no parasite speciation, and v is just moved to one of the two arcs outgoing from h(a) chosen randomly.

Notice however that, in order to avoid confusing a loss with another event (for instance, a cospeciation), some precautions must be taken, as explained more specifically in the next paragraph concerning the simulation of a loss event.

These choices, together with the general framework for our parasite tree generation method, are provided next.

Starting the generation.-The generation of the simulated parasite tree P starts with the creation of its root vertex Proot . This vertex is positioned before the root of H on the arc a = (ρ, H root ). This allows the simulation of events that happened in the parasite tree before the most recent common ancestor of all host species in H. The evolutionary events.-For any vertex v of P that is not yet mapped and whose position is v, a (Fig. 2 b)), we choose to apply one among the four allowed operations, depending on the probability of each event. In what follows, we denote by a 1 , a 2 the arcs outgoing from the head h(a) of the arc a.

• Cospeciation (Fig. 2 c)): We apply the mapping [v : h(a)] and we create the vertices v 1 and v 2 as children of v. We position them as follows: v 1 , a 1 and v 2 , a 2 . This operation is executed with probability p c .

• Duplication (Fig. 2 d)): We apply the mapping [v : h(a)] and we create the vertices v 1 and v 2 as children of v. Both v 1 and v 2 are positioned on a. This operation is executed with probability p d .

• Host switch (Fig. 2 e)): We apply the mapping [v : h(a)] and we create the vertices v 1 and v 2 as children of v. We then randomly choose one of the two children and position it on a. Finally, we randomly choose an arc a that does not violate the time feasibility of the reconstruction so far [START_REF] Stolzer | Inferring duplications, losses, transfers and incomplete lineage sorting with nonbinary species trees[END_REF]). If such an arc does not exist, it is not possible for a host switch to take place. In this case, we choose between the three remaining events with probability p i /(p c + p d + p l ) with i ∈ {c, d, l}.

Otherwise, we position v 2 on a . This operation is executed with probability p s .

• Loss (Fig. 2 f)): This operation is executed with probability p l and consists of randomly choosing an arc outgoing from the head h(a) of a and positioning v on it.

Observe that we are considering only losses resulting from lineage sorting. It would be interesting to incorporate extinction or failure to detect infection but this would require the addition of new parameters, thus making the model more complex to analyse. However, if v was created by a duplication event and is being processed for the first time, we have to verify if its sibling vertex v was already processed and also suffered a loss. In this case, v must be positioned on the same arc a where v was positioned. This procedure is adopted to avoid later mappings where a duplication followed by two losses would be confused with a cospeciation.

We also assume that no evolutionary event takes place whenever a leaf of H is reached. This means that, if v is positioned on an arc incoming to a leaf, then v is mapped to the leaf and no further operation is executed. Hence, the generation of P terminates when all the created vertices are mapped (i.e. have reached a leaf of the host tree). Finally, the leaves of the parasite tree P are labelled according to their mapping to the leaves of the host tree. Observe that as more than one parasite can be mapped to the same host, P is a multi-labelled tree (that is, trees whose leaf labels need not be unique). Finally, some combinations of host switches can introduce an incompatibility due to the temporal constraints imposed by the host and parasite trees, as well as by the reconciliation itself.

During the generation of the parasite tree, we always allow only for host switches that do not violate the time-feasibility constraints. For the criteria enabling to assess time-feasibility, we refer to [START_REF] Stolzer | Inferring duplications, losses, transfers and incomplete lineage sorting with nonbinary species trees[END_REF].

Note that in this model, we do not use information about edge lengths. This is a positive aspect of the method in the sense that branch lengths are not always easy to determine with accuracy. In contrast, we cannot simulate the "null events" (parasite doing nothing in the host tree). Moreover, for now, we do not simulate "failure to diverge" which describes a situation where a host speciates while the parasite does not but continues to inhabit both of the two new species of hosts. Despite the importance of this event, mathematically speaking it is not clear how to include it in the cophylogenetic reconciliation model since we have to allow the association of a parasite to multiple hosts.

The ideas presented by [START_REF] Drinkwater | An improved node mapping algorithm for the cophylogeny reconstruction problem[END_REF] for the improvement of node mapping algorithms may help on the simulation of the "failure to diverge" event in future work.

Since the simulation model is restricted to the events of cospeciation, duplication, host switch, and loss, the probabilities of these four events sum up to one.

Cophylogeny parameter estimation algorithm

Prior distribution π.-The parameter θ = p c , p d , p s , p l lives in the simplex S 3 (the p's are positive and sum to one). It is then standard to sample θ from a Dirichlet distribution which is a family of continuous multivariate probability distributions parameterised by a vector α of positive real numbers that determine the shape of the distribution [START_REF] Gelman | Bayesian data analysis[END_REF].

In our simulations, we adopt a uniform Dirichlet distribution (namely α = (1, 1, 1, 1)) that corresponds to sampling uniformly from the simplex S 3 . This is often used when there is no previous knowledge favouring one component (e.g. coevolutionary event) of θ over another. However, the method we implemented allows the user to specify other prior distributions when such knowledge is available.

Choice of summary statistic and fitting criterion.-The ABC inference method is based on the choice of a summary statistic that describes the data while performing a dimension reduction task. The latter is used to evaluate the quality of agreement (similarity) between the simulated datasets (the generated parasite trees) and the observed (the real parasite tree). In our case, the summary statistic will be based on the measured distances between the generated parasite trees and the real one.

The distance of each simulated tree to the real parasite tree is therefore informative about the quality of the vector that generated it. Hence, the distance that will be used must take into account: (i) how well does the simulated tree represent the set of trees generated by a given vector, and (ii) how topologically similar is the simulated tree to the real parasite tree.

Concerning the first point, the intuition is as follows. In our model, when generating a parasite tree, the expected frequency of an event should be close to the corresponding probability value of the parameter vector used to generate the tree. To this purpose, for a

given vector θ = p c , p d , p s , p l and for each simulated tree P θ that was generated according to this vector, we kept track of the number of events obs = o c , o d , o s , o l associated with this simulation. We compared the observed number of events to the expected exp = e c , e d , e s , e l . Observe that the expected number of events can be easily calculated using the size of the parasite tree P and the vector θ. A tree is a good representative if the observed number of events is near to the expected. More formally, for a real parasite tree P , a vector θ = p c , p d , p s , p l , and a simulated parasite tree P θ for which the observed number of events are obs = o c , o d , o s , o l , we define a measure D 1 (P, P θ ) as follows:

D 1 (P, P θ ) = 1 4 × i∈{c,d,s,l} |e i -o i | max{e i , o i } .
As concerns point (ii), we use a metric for comparing phylogenetic trees. There is a wide literature on distances for phylogenetic trees [START_REF] Felsenstein | Inferring phylogenies[END_REF]. Our choice was driven by the need to have one that can be computed efficiently and accurately.

Unfortunately, many of the distances used in biology are also NP-hard to compute [START_REF] Waterman | On the similarity of dendrograms[END_REF][START_REF] Hein | Reconstructing evolution of sequences subject to recombination using parsimony[END_REF][START_REF] Baroni | Bounding the number of hybridisation events for a consistent evolutionary history[END_REF], whereas some of the fastest, like for instance, the Robinson-Foulds distance [START_REF] Robinson | Comparison of phylogenetic trees[END_REF] which can be calculated in linear time [START_REF] Day | Optimal algorithms for comparing trees with labeled leaves[END_REF], are poorly distributed and thus not good enough discriminators [START_REF] Steel | Distributions of tree comparison metrics -some new results[END_REF][START_REF] Bryant | Computing the distribution of a tree metric[END_REF]. Moreover, many efficient-to-compute distances are not robust to small changes (such as in the position of a single leaf) in one of the two trees.

Recall that in our method the leaves of the parasite tree P are labelled according to their mapping to the leaves of the host tree and that more than one parasite can be mapped to the same host. Hence, we are interested in distances between multi-labelled trees.

In our context, the distance that best meets the requirement of efficiency and accuracy appears for now to be the maximum agreement area cladogram (MAAC) [START_REF] Ganapathy | Pattern identification in biogeography[END_REF]. This is a generalisation for multi-labelled trees of the well-known maximum agreement subtree [START_REF] Finden | Obtaining common pruned trees[END_REF][START_REF] Farach-Colton | On the agreement of many trees[END_REF] and it corresponds to the number of leaves in the largest isomorphic subtree that is common to two (multi-labelled) trees. Clearly this isomorphism takes into account the labels of the trees. The MAAC distance can be calculated in O(n 2 ) time where n is the size of the largest input tree [START_REF] Ganapathy | Pattern identification in biogeography[END_REF].

We use a normalised version of MAAC that takes into account also the number of leaves in common between the two trees. More formally, for two trees P and P with leaf sets L(P ) and L(P ) respectively, we define the measure D 2 (P, P ) as follows:

D 2 (P, P ) =        1 - M AAC(P, P ) |L(P ) ∩ L(P )| if L(P ) ∩ L(P ) = ∅ 1 otherwise.
Observe that the intersection operation involves multi-sets. We recall that a multi-set is a generalisation of a set where the elements are allowed to appear more than once, hence the operations take into account their multiplicity in the following way: if the multiplicity of an element e in a multi-set A is given by [e](A), then [e](A ∩ B) is given by min{

[e](A), [e](B)}.
Finally, we propose a distance that is based on these two components D 1 and D 2 .

For a real parasite P , a vector θ = p c , p d , p s , p l , and a simulated parasite tree P θ , we define the distance d(P, P θ ) as follows:

d(P, P θ ) = β 1 D 1 (P, P θ ) + β 2 D 2 (P, P θ ), with β 1 + β 2 = 1.
According to our experiments (see Supplementary Material, http://dx.doi.org/10.5061/dryad.9q5fp), the most appropriate values are β 1 = 0.7 and β 2 = 0.3 but this can be set by the user. The main drawback of this distance is that it is not a metric; however it achieves good results with respect to discriminating the trees as observed in our experiments.

In Coala, we implemented two other distances, both of which are variations of the MAAC. A user can choose the most appropriate one depending on the case. In this paper, we show only the results for the two-component distance, as this had the most discriminating power (data not shown).

Given a parameter vector θ = p c , p d , p s , p l , we generate M trees and for each of them we consider the distance of P from the real parasite tree P . From this set of distances, D, we produce a summary, denoted by S(θ), that characterises the set of trees generated with the parameter vector. In our experiments, we choose S(θ) as the average of all the produced distances.

The summary S(θ) is the value that is used in the rejection/acceptance step of the ABC method.

Finally, it is worth noting that while the choice of a summary statistic (or equivalently here a summary tree distance) is independent from the generation process (coevolution model), such a choice may have a deep impact on the performance and the results of the method. This is one of the main issues with ABC-related methods. Some recent works have attempted to improve this step [START_REF] Fearnhead | Constructing summary statistics for approximate Bayesian computation: semi-automatic approximate Bayesian computation[END_REF]. From the experiments done however, we can already see that the two-component distance seems to be a good enough discriminator.

Approximate Bayesian Computation -Sequential Monte Carlo procedure.-The ABC-SMC procedure is composed of a sequence of R > 1 rounds. For each of these rounds, we define a tolerance value τ r (1 ≤ r ≤ R) which determines the percentage of parameter vectors to be accepted. Associated with a tolerance value τ r , we have a threshold r which is the largest value of the summary statistic associated with the accepted parameter vectors.

• Initial round (r = 1):

Draw an initial set of N parameter vectors {θ i 1 } (1≤i≤N ) from the prior π. Then,

for each θ i 1 generate M trees { Pj (θ i 1 )} (1≤j≤M ) . Select Q = τ 1 × N parameter
vectors θ 1 that have the smallest S(θ 1 ), thus defining the threshold 1 and the set A 1 of accepted parameter vectors.

• Following rounds (2 ≤ r ≤ R):

1. Sample a parameter vector θ from the set A (r-1) . Create a parameter vector θ by perturbing θ . The perturbation is performed by adding to each coordinate of θ a randomly chosen value in [-0.01, +0.01] and normalising it. 3. Based on the set S r , select τ r × Q parameter vectors θ r that have the smallest S(θ r ), thus defining the threshold r and the set A r of accepted parameters.

The final set A R of accepted parameter vectors is the result of the ABC-SMC procedure and characterises the list of parameter vectors that may explain the evolution of the pair of host and parasite trees given as input.

Let us observe that, since in all our experiments we are assuming a uniform prior distribution and also are performing the perturbations in a uniform way, the weights induced by the proposals appear to be uniform [START_REF] Beaumont | Adaptive approximate Bayesian computation[END_REF]). However, in the case of a different prior, weights should be used in the process in order to correct the posterior distribution according to the perturbation made.

Clustering the results

Coala implements a hierarchical clustering procedure to group the final list of accepted parameter vectors. The basic process of a hierarchical clustering is as follows. At the beginning, each parameter vector forms a single cluster. Then at each step, the pair of clusters that have the smallest distance to each other are merged to form a new cluster.

The distance that we use between the vectors θ = p c , p d , p s , p l and θ = q c , q d , q s , q l is the χ 2 distance, which is a weighted Euclidean distance defined as follows:

d(θ, θ ) = i∈{c,d,s,l} 2 × (p i -q i ) 2 (p i + q i ) .
At the end of this process, we have a single cluster containing all the items represented as a tree (hierarchical cluster tree or dendrogram) showing the relationship among all the original items. As we make no assumptions concerning the space of the vectors we are dealing with, we chose to apply a more general but still efficient method, introduced in [START_REF] Langfelder | Defining clusters from a hierarchical cluster tree: the Dynamic Tree Cut package for R[END_REF], to select the branches to be cut in the dendrogram.

The method proceeds in two steps. Starting with the complete dendrogram, it first identifies preliminary clusters that satisfy some criteria: for example they contain a certain minimum number of objects (to avoid spurious divisions), any two clusters are at least some distance apart, etc. [START_REF] Langfelder | Defining clusters from a hierarchical cluster tree: the Dynamic Tree Cut package for R[END_REF], for more details). In a second step, all the items that have not been assigned to any cluster are tested for sufficient proximity to preliminary clusters; if the nearest cluster is close enough, the item is assigned to that cluster, otherwise the item remains clustered according to the complete dendrogram.

Finally, once the vectors are split into clusters, we associate to each one a representative parameter vector. To define each coordinate of the "consensus" parameter vector, we take the mean value of the respective coordinate in all the parameter vectors which are inside the cluster. We then normalise the "consensus" coordinates to sum to one.

Experimental Results and Discussion

We evaluated our method in two different ways. First we designed a self-test to

show that the principle underlying it is sound and to test it on simulated datasets.

We then extended the evaluation to four real examples that correspond to biological datasets from the literature. This choice was dictated by: (1) the availability of the trees and of their leaf mapping; and (2) the desire to, again, cover for situations as widely different as possible in terms of the events supposed to have taken place during the host-parasite coevolution. As a matter of fact, the first point drove the choice more than the second: there are not so many examples available from which it is easy to extract the tree and/or leaf mapping and that are big enough to represent meaningful datasets on which to test Coala. All four examples were also analysed in the original paper from which they were extracted by one or more of the existing algorithms that search for a most parsimonious (possibly cyclic) reconciliation (i.e. for a reconciliation of minimum cost).

Except in one case, which is a heuristic strategy and therefore does not guarantee optimality of the solution, all existing algorithms need to receive as input the cost of the events, which is thus established a priori and drives the conclusions on the results obtained.

Finally, we applied Coala to a biological dataset of our own, representing the coevolution of bacteria from the Wolbachia genus and the various arthropods that host them. This dataset was selected because of its size: the trees have each 387 leaves.

Experimental parameters

All datasets were processed by Coala configured with the same parameters. For each dataset, we generated N = 2000 parameter vectors in the first round. For each of the vectors, we generated M = 1000 parasite trees using our method. We required these trees to have a size at most twice the one of the real parasite tree, otherwise the tree was discarded as being too different from the original. If a given vector did not generate M such trees in 5000 trials, then the vector was immediately associated with a distance equal to 1 which indicated that it represented the real data badly.

We used the average of all the 1000 distances produced as a fitting criterion in the rejection/acceptance step of the ABC method. The tolerance value used in the first round was τ 1 = 0.1. For the remaining rounds 2 ≤ i ≤ R, we defined τ i = 0.25. Notice that τ 1 × N = 200 defines the size Q of the quantile set which must be produced in each new round. Thus, after the last round, we have τ R × Q = 50 accepted vectors. These vectors are grouped into clusters and a representative vector is associated with each cluster as explained in the Section "Clustering the results".

We ran the experiments using R = 3 and R = 5 rounds. The number of rounds is an important parameter, which defines the characteristics of the list of accepted parameter vectors.

However, observe that a high number of rounds will tend to overfit the data and thus hide a possible variability in the list of accepted vectors that could provide significant alternatives for explaining the studied pair of trees.

Since we are interested in exploring different alternatives for each dataset, we present only the results which were obtained after running Coala for 3 rounds. The results involving 5 rounds may be found in the Supplementary Material.

Simulated datasets

We first evaluated our model on simulated data. Clearly, in order to do this, we have to generate the phylogenies for the hosts and parasites whose coevolution is being studied in such a way that the probability of each event is known. The basic idea is that if we are able to select a "typical" (or representative) parasite tree P θ that is generated starting from a host tree H and a given probability vector θ, Coala should be able to list values close to θ among the vectors accepted in the last round.

It is important to observe that many different probability vectors can explain the same pair of trees. We will therefore consider it acceptable if Coala produces clusters that are relatively close to θ.

Generating simulated datasets.-Due to the high variability of the parasite trees which can be simulated given a host tree H and a vector θ, the task of choosing the most "typical" tree can be hard. To simplify this task and select a typical tree, we impose two conditions which must be observed by the simulated tree. The first one requires that the candidate tree should have a size close to the median for all the trees which are simulated using H and θ. The second condition requires that the observed number of events of a candidate tree should be very close to the expected number given θ.

In practical terms, we execute the following procedure: in order to get realistic datasets we choose a real host tree H (see Supplementary Material for more information).

Then, given a probability vector θ and H, we generate 2000 parasite trees using our model, without imposing any limit on the size of the generated trees. We then compute the median size of all generated trees and we filter out those whose size is far from this value (difference greater than 1 or 2 leaves from the median value). Finally, we select as typical tree P θ the one that shows the smallest χ 2 distance between the vector θ and the vector of observed frequencies of events.

We generated in this way 9 datasets (H,P ) associated with the following 9 probability vectors: θ 1 = 0.70, 0.10, 0.10, 0.10 , θ 2 = 0.80, 0.15, 0.01, 0.04 , θ 3 = 0.75, 0.01, 0.16, 0.08 , θ 4 = 0.70, 0.05, 0.02, 0.23 , θ 5 = 0.60, 0.20, 0.00, 0.20 , θ 6 = 0.55, 0.00, 0.20, 0.25 , θ 7 = 0.45, 0.10, 0.15, 0.30 , θ 8 = 0.40, 0.20, 0.10, 0.30 and θ 9 = 0.30, 0.20, 0.40, 0.10 (see the Supplementary material for more details). The choice of vectors was done with the aim to cover different patterns of probability. All datasets were generated with the same host tree H of 36 leaves.

Self-Test.-As concerns the self-test, we designed the following procedure. Let P θ denote the simulated parasite tree chosen in correspondence of the probability vector θ, as explained in the previous section. We recall that the host tree H remains the same during all the self-test experiments. For a pair of host and parasite trees (H, P θ ), we ran Coala 50 times. In each run j, we computed the quality q j that corresponded to how well the method was able to recover the target vector θ used for generating the dataset P θ . To do this, for each run j, we considered the representative vectors of the clusters produced as output. We computed the χ 2 distance for each of the representative vectors to the target vector θ and set q j to the smallest value among them.

Figure 3 shows the distribution of the quality values which were obtained at the end of each round (from 2 to 5) for the simulated datasets θ 1 , θ 3 , θ 4 , and θ 7 (the results for the remaining datasets can be found in the Supplementary Material). Figure 4 shows the histograms of the event probabilities observed for the 50 parameter vectors with smallest χ 2 distance at the end of each round for dataset θ 3 (again, the results for the remaining datasets are available in the Supplementary Material).

Up to a certain level of cospeciation probability (≥ 0.50), our results (Figure 3)

show that in the rounds 2 and 3, Coala is able to select parameter vectors that are close to the target probability vector. Looking to the histograms of these two rounds, we can observe that in most of the runs, the closest parameter vector has low χ 2 distance to the target. After the third round, this tendency changes and the closest parameter vectors show high χ 2 distances indicating that Coala is mainly selecting vectors which are far from the target one.

Since Coala is based on an ABC-SMC approach, the accepted vectors in one round have summary statistics (i.e. average distance) smaller than the defined in the previous round. This means that at each new round, Coala is selecting parameter vectors that have more probability of explaining the pair of trees given as input because their simulated parasite trees are, on average, closer to the real one.

Although we try to choose the best representative parasite tree P for each pair (H, θ), we cannot guarantee that θ is the best explanation for the association between H and P . Even so, Coala was able to select parameter vectors that are close to the target probability vector in the first rounds. Figure 4 shows the histograms of the event probabilities observed among the 50 parameter vectors with smallest χ 2 distance at the end of each round for dataset θ 3 , and confirms these observations. We can see that at round 2, the median and mean event probabilities (solid and dotted vertical lines respectively) are very close to the target value (dashed vertical line). When we increase the number of rounds, the distance between the median/mean probabilities and the target values increases.

When we decrease the cospeciation probability to values smaller than 0.50, Coala selects very few vectors which are close to the target vector. When the cospeciation probability decreases while the duplication and host switch probabilities increase, the variability of the tree topologies observed increases exponentially. Due to this, selecting a typical tree becomes an almost impossible task and this may explain the obtained results.

Increasing the number of simulated trees to compute the summary statistic might enable us to improve the quality of the results. However, this would require a much longer execution time.

Biological datasets extracted from the literature

To evaluate Coala on biological datasets, we extracted four pairs of host and parasite trees from the literature. However, due to space issues, in this work we present only two of them. A description and the results obtained on the additional biological datasets can be found in the Supplementary Material. Before presenting and discussing the datasets, we provide details on how we performed the analyses.

Each dataset was processed by Coala as described in the Section "Experimental parameters". Table 2 shows the representative parameter vectors obtained for each one of the datasets and Figure 5 the histograms of the event probabilities of the list of accepted vectors obtained at the end of the third round.

In order to compare our results to the existing literature, we transformed each one of the representative parameter vectors p c , p d , p s , p l into a vector of costs that was then used to compute optimal reconciliations between the host and parasite trees given as input.

The transformation was done by defining c i = -ln p i , with i ∈ {c, d, s, l}, which is based on a commonly accepted idea that the cost of an event is inversely related to its probability [START_REF] Charleston | Jungles: a new solution to the host/parasite phylogeny reconciliation problem[END_REF][START_REF] Ronquist | Tangled trees: phylogeny, cospeciation, and coevolution chap. Parsimony analysis of coevolving species associations[END_REF][START_REF] Huelsenbeck | Statistical tests of host-parasite cospeciation[END_REF], for example). Indeed, if p i is equal to 1, then we expect all the events to be of type i, thus the cost of the corresponding event must be 0. Similarly, if p i is equal to 0, we expect that event i never happens, and thus the cost must be assigned to +∞.

To the best of our knowledge, the only methods that enumerate all optimal reconciliations are CoRe-Pa [START_REF] Merkle | A parameter-adaptive dynamic programming approach for inferring cophylogenies[END_REF], Notung [START_REF] Stolzer | Inferring duplications, losses, transfers and incomplete lineage sorting with nonbinary species trees[END_REF] and

Eucalypt [START_REF] Donati | Eucalypt: Efficient tree reconciliation enumerator[END_REF]. However CoRe-Pa in some cases misses solutions, probably because it considers some additional constraints. Notung does not allow cospeciation costs different from zero and the remaining event costs must be described by integer values. We thus present the results of Eucalypt which allows the configuration of all event costs and accepts real numbers.

Table 3 shows, for each dataset, the vector of costs (c c , c d , c s , c l ) produced by transforming the representative parameter vectors obtained after the third round (Table 2).

Column Opt indicates the cost of the optimal solution and columns #c, #d, #s, #l the numbers of each event type which are observed among the enumerated scenarios. Finally, columns #A and #C indicate, respectively, the total number of acyclic and cyclic scenarios. Running Coala on this dataset, we obtain 3 non-singleton clusters which are quite different from each other (Table 2). Cluster 0 is formed by a single accepted vector which did not cluster with any other because it is too far apart. Cluster 1 shows probabilities of 0.46, 0.26 and 0.28, respectively, for cospeciation, duplication and loss. After transforming these into costs (Table 3) Figure 2) and considers the coevolution of hantaviruses with their insectivore and rodent hosts. The host tree consists of a total of 34 hosts (28 rodents and 6 insectivores) and the parasite tree includes 42 hantaviruses. It was strongly believed that hantaviruses cospeciated with rodents since their phylogenetic trees have topological similarities with three consistently well-defined clades [START_REF] Hughes | Evolutionary diversification of protein-coding genes of hantaviruses[END_REF][START_REF] Plyusnin | Virus evolution and genetic diversity of hantaviruses and their rodent hosts[END_REF][START_REF] Nemirov | Hantaviruses: co-evolution with natural hosts[END_REF][START_REF] Jackson | A cophylogenetic perspective of RNA-virus evolution[END_REF]. The authors show that to support this hypothesis, the evolutionary rate of the RNA sequences of the hantaviruses should be several orders of magnitude smaller than the rates which are normally observed in RNA viruses that replicate with RNA-dependent RNA polymerase [START_REF] Hanada | A large variation in the rates of synonymous substitution for RNA viruses and its relationship to a diversity of viral infection and transmission modes[END_REF].

By analysing the cophylogenetic reconciliations, the authors show that scenarios with more than 20 cospeciations are statistically non-significant. To explain the topological congruences, the authors point to the fact that host-switching followed by pathogen speciation can generate congruence between trees, particularly when pathogens preferentially switch among closely related hosts. Based on this fact and on the observed patterns of amino acid replacement observed in these viruses (compatible with host-specific adaptation), the authors conclude that the coevolutionary history of these hosts and parasites is the result of a recent history of preferential host-switching and local adaptation.

Looking at Table 2, we can observe that Clusters 1, 2, and 3 have representative vectors with zero probability for host switch events: Cluster 1 has a very high cospeciation probability (0.85), while Clusters 2 and 3 have probability values which are almost equally distributed among cospeciation, duplication and loss events. After transforming these vectors into costs (Table 3), we obtain scenarios with a high number of cospeciations which is considered non-significant by [START_REF] Ramsden | Hantavirus evolution in relation to its rodent and insectivore hosts[END_REF].

Differently from the others, Cluster 4 shows a vector with host switch probability higher than the probabilities of duplication and loss. When converted into costs (Table 3), this generates time-consistent scenarios with 17 cospeciations, 5 duplications, 19 host switches and 4 losses, a result much closer to the explanation given by [START_REF] Ramsden | Hantavirus evolution in relation to its rodent and insectivore hosts[END_REF]. These results reinforce the idea that, although Coala is able to identify vectors which can explain a pair of trees, having a prior knowledge of the dynamics of the interactions of the two groups of species is important to identify the clusters that better explain their coevolution.

Wolbachia and their arthropod hosts dataset

Wolbachia is a large, phylogenetically diverse monophyletic genus of intracellular bacteria that are currently considered the most abundant endosymbionts in arthropods. In insects alone, it is estimated that over 65% of the species are infected by Wolbachia. The dataset used in this paper corresponds to Wolbachia species that were detected in an extensive set of arthropods collected from 4 young, isolated islands (less than 5 Myr old) [START_REF] Simões | Wolbachia detection: an assessment of standard pcr protocols[END_REF][START_REF] Simões | Diversity and dynamics of Wolbachia-host associations in arthropods from the Society archipelago[END_REF]. The trees are a subset of those discussed in [START_REF] Simões | Wolbachia detection: an assessment of standard pcr protocols[END_REF][START_REF] Simões | Diversity and dynamics of Wolbachia-host associations in arthropods from the Society archipelago[END_REF], where we retained only those parasites which were associated with a unique host, the hosts diverge by at least 2% at the level of the CO1 genes that were used for reconstructing their phylogenetic tree and the Wolbachia sequences (corresponding to the fbpA gene) differ by at least one SNP. Each resulting tree is composed of 387 leaves. The initial results presented in [START_REF] Simões | Diversity and dynamics of Wolbachia-host associations in arthropods from the Society archipelago[END_REF] seemed to indicate that host switches might be quite frequent even among hosts that are physiologically and molecularly very distinct and thus phylogenetically distant.

The Wolbachia-arthropods dataset was also processed by Coala as described in the Section "Experimental parameters". Table 4 shows the three clusters which were obtained at the end of the third round. All these clusters have significantly high cospeciation probabilities (> 0.77). The first cluster has a very low duplication probability and a host switch probability around 0.5. The two other clusters point to a relatively high duplication probability and low level of host switches. The difference between them is related to the probability of losses, which is around 0.14 for Cluster 2 and zero for Cluster 3.

Cluster 1 goes in the direction of what was presented in [START_REF] Simões | Diversity and dynamics of Wolbachia-host associations in arthropods from the Society archipelago[END_REF] where the author suggested that in the last 3 Ma, there were many transfers of Wolbachia, including between different arthropod orders, i.e. over large phylogenetic distances. Clusters 2 and 3 point to an opposite scenario.

Similarly to the analysis performed for the small biological datasets, we transformed each one of the representative parameter vectors into a vector of costs that was then used to compute optimal reconciliations between the host and parasite trees given as input.

What is most striking with the results obtained for this dataset is the absolutely huge number of optimal reconciliations that can be derived for all clusters. Since the total number of solutions makes impossible the enumeration of all the results, for this dataset, we therefore only computed the costs of the optimal solutions and the total number of solutions. Additionally, for each cluster, we sampled 10000 solutions and we checked for the presence of acyclic solutions. Table 5 summarises the results obtained.

For the small sampling that we performed, we were able to find feasible (acyclic) solutions only with the cost vector produced with the event probabilities of Cluster 3.

However, the results obtained with all the other four datasets used here lead us to suggest that the number of feasible solutions might quite possibly remain large.

Conclusions

We have developed an automated method that, starting from two phylogenies representing sets of host and parasite species, allows extraction of information about the costs of the events in a most probable reconciliation. It is clear that within a parsimony-based approach, an optimal solution strictly depends on the specific values attributed to these costs. However, there seldom is enough information for assigning those values a priori. Indeed, we observe in the results we obtained on a diverse selection of datasets that the costs inferred by our simulations may be very different across datasets, thus motivating the use of estimated instead of fixed costs. Such costs may even differ widely for a same pair of host-parasite trees, as is observed for the Wolbachia-arthropods dataset.

These costs are inversely related to their likelihood, and so to their expected frequency. For this reason, providing information on the frequencies of the events is an important issue, in particular in the cases where the reconciliation methods fail to find a solution. The latter can happen, for instance, if all the optimal solutions that are identified by the existing reconciliation algorithms are biologically unfeasible due to the presence of cycles, since finding an acyclic reconciliation is an NP-hard problem. In addition, if the host and parasite trees are large (for instance, on the order of hundreds of taxa), these cases cannot be handled by the existing reconciliation algorithms in the sense that there are too many solutions to test for acyclicity.

As a future work, we first plan to refine the model used for the reconciliation problem, including more biological information and making it more realistic. In particular, we could include information about the distance of the allowed host switches (for instance if we expect a host switch to rarely happen between species that are too far from each other), or allow the mapping of the leaves to be an association instead of a function (thus addressing the cases where a parasite can be found in more than one host species).

Moreover, we should also consider the case where the input phylogenies are not fully resolved, meaning that the trees are not binary.

A more efficient exploration of the parameter space is another important future issue that would significantly increase the efficiency of our procedure and also allow to handle larger trees.

It is important to observe that most studies on cophylogeny assume that the phylogenies of the organisms are correct. Clearly, this may affect the results observed. It would therefore be interesting to be able to infer the cophylogenetic reconciliation directly from sequence data.

Finally, the accuracy of the results obtained by our method depends on the choice of the metric used for comparing trees. Designing new metrics that can be computed efficiently while still capturing the similarity for multi-labelled, not fully resolved trees is therefore another important future issue which we believe is also interesting per se.

Tables Table 1: Notation

N

  Figure 2 a) depicts this initial configuration.

  2. Generate M trees { Pj (θ )} (1≤j≤M ) and compute S(θ ). If S(θ ) ≤ (r-1) , add θ into the quantile set S r . If |S r | < Q, return to Step 1.

Dataset 1 -

 1 Flavobacterial endosymbionts and their insect hosts.-This dataset was extracted from the work of[START_REF] Rosenblueth | Evolutionary relationships of flavobacterial and enterobacterial endosymbionts with their scale insect hosts (hemiptera: Coccoidea)[END_REF] and is composed of a pair of host and parasite trees which have each 17 species (see Supplementary Material). The parameter adaptive approach of CoRe-Pa[START_REF] Merkle | A parameter-adaptive dynamic programming approach for inferring cophylogenies[END_REF]) was used to infer the more appropriate cost vectors for analysing this dataset. Nine such vectors were produced.However, only one, c c = 0.088, c d = 0.325, c s = 0.339, c l = 0.248 , was associated with a feasible reconciliation in the sense that host switches happened between contemporary species only (the branch length was used to infer this information). Since CoRe-Pa can produce unfeasible (i.e. cyclic) solutions during the parameter adaptive approach, Rosenblueth et al. decided to complement their study with Jane 3[START_REF] Conow | Jane: a new tool for the cophylogeny reconstruction problem[END_REF], which uses a genetic algorithm approach to produce only acyclic reconciliations. They thus started with the only cost vector obtained by CoRe-Pa associated with a feasible reconciliation, however transforming it into integer numbers (a requirement of the software), and then gradually changed the costs until a feasible reconciliation was produced (again using branch-length information). This procedure resulted in the cost vectorc c = 1, c d = 1, c s = 1, c l = 2and a reconciliation with 9 cospeciations, 0 duplication, 7 host switches and 1 loss, the same as obtained by CoRe-Pa.

  , the obtained reconciliation scenarios have 11 cospeciations, 2 duplications, 3 host switches and 11 losses. Clusters 2 and 3 show very low duplication probability. While Cluster 2 exhibits intermediate values for the remaining probabilities, Cluster 3 has a very high cospeciation probability value (0.91) and low host switch (0.06)and loss (0.02). Due to the low duplication value, these clusters show the same reconciliation scenario: 9 cospeciations, 0 duplications, 7 host switches and 1 loss, which is identical to the one proposed by[START_REF] Rosenblueth | Evolutionary relationships of flavobacterial and enterobacterial endosymbionts with their scale insect hosts (hemiptera: Coccoidea)[END_REF].Dataset 2 -Rodents and Hantaviruses.-This dataset is taken from[START_REF] Ramsden | Hantavirus evolution in relation to its rodent and insectivore hosts[END_REF]

FiguresFigure 1 :Figure 2 :

 12 Figures

Table 2 :

 2 Function from the leaves of P to the leaves of H. It represents the associations between currently living host species and parasites. γ Function from the vertices of P to the vertices of H. It represents the reconciliation between H and P and extends ϕ. Σ, ∆, Γ Sets of parasite vertices associated with, respectively, cospeciation, duplication, and host switch events. Ξ Set containing arcs of the parasite tree that are associated to host switch events. Λ Multi-set containing all vertices h ∈ V (H) that are associated to loss events. Probability of the event i, where i ∈ {c, d, s, l}. c i Cost of the event i, where i ∈ {c, d, s, l}. o i Number of observed events of the type i, where i ∈ {c, d, s, l}. Note: c = cospeciation, d = duplication, s = host switch, and l = loss. Representative probability vectors produced by Coala at Round 3.

	otation Description
	H	Host tree
	P	Parasite tree
	ϕ
	D 0	Observed data.
	D S	Generated data.
	Θ	Parameter space.
	θ	Parameter value.
	P	Simulated parasite tree.
	p i

Table 3 :

 3 Event vectors obtained by transforming the probability vectors (Table2) into cost vectors. c, #d, #s, and #l denote the number of each event type which are observed among the enumerated scenarios. #A and #C indicate, respectively, the total number of acyclic and cyclic scenarios.

	Dataset Cluster	c c	c d	c s	c l	Opt	#c #d #s #l #A #C
		0	3.517 13.816 0.584 0.885 14.044 1	0 15 2 2944 0
	1	1 2	0.775 1.355 7.824 1.270 48.664 11 2 0.591 8.517 1.310 1.736 16.217 9 0	3 11 7 1	2 1	0 0
		3	0.094 4.160 2.844 4.154 24.892 9	0	7	1	1	0
		1	0.161 2.496 9.210 2.717 153.544 22 11 8 18	0	12
	2	2 3	0.748 1.592 9.210 1.130 105.393 22 19 0 52 1.436 1.053 8.112 0.884 97.548 22 19 0 52	1 1	0 0
		4	0.545 6.266 1.265 1.996 72.588 17 5 19 4	4	0
	Note: #

Table 4 :

 4 Representative probability vectors produced by Coala, at the end of the third round, while processing the Wolbachia-arthropods datasets.

	Cluster	p c	p d	p s	p l	#vectors
	1	0.866 0.006 0.055 0.073	26
	2	0.771 0.078 0.010 0.141	22
	3	0.964 0.022 0.014 0.000	2

Table 5 :

 5 Total number of solutions obtained by transforming the probability vectors (Table4) into cost vectors for Wolbachia-arthropods datasets.

	Cluster	c c	c d	c s	c l	Opt	Solutions Acyclic solutions
	1	0.144 5.116 2.899	2.623	917.475	5.4 × 10 43	No
	2	0.260 2.551 4.595	1.961	1407.877 9.8 × 10 40	No
	3	0.037 3.817 4.269 13.816 1375.725 1.6 × 10 51	Yes

Availability

The Coala program is available at http://coala.gforge.inria.fr/ and runs on any machine with Java 1.6 or higher. The Eucalypt program is available at http://eucalypt.gforge.inria.fr/.

Data available from the Dryad Digital Repository at http://dx.doi.org/10.5061/dryad.9q5fp.

Running time

The experiments were executed at the IN2P3 Computing Center (http://cc.in2p3.fr/).

For the simulated datasets, each pair of trees was processed with 3 threads for speeding up the simulation process. The time necessary to complete 5 rounds for all 50 runs varied from 1 to 2 days depending on the size of the trees. For the biological datasets 1 to 4, we also used 3 threads. The observed execution times for 5 rounds were between a couple of hours for the smallest dataset (Dataset 1) and one day for Dataset 4. Due to its size, the dataset Wolbachia-arthropods was processed with 150 threads and it required approximately 8 days to complete 5 rounds.

Supplementary Material

Funding

The research leading to these results was funded by the Figure 3: For each simulated dataset, we ran Coala 50 times and, at the end of each round (from 2 to 5), we took note of the cluster whose representative parameter vector had the smallest χ 2 distance to the probability vector used to generate the simulated dataset. The histograms show the distribution of the smallest χ 2 distance observed on each one of the 50 runs at the end of each round (for the simulated datasets


Sylvain Malacria

Jonathan Aceituno

Philip Quinn

Géry Casiez

Andy Cockburn

Nicolas Roussel

Géry Casiez

Push-Edge and Slide-Edge: Scrolling by Pushing Against the Viewport Edge

Keywords: Edge-scroll, rate control, position control, inertia, autoscroll. ACM Classification H.5.2 [User Interfaces] Interaction styles

The cat was sitting in the kitchen when the lie-no decided, rather than being

ABSTRACT

Edge-scrolling allows users to scroll a viewport while simultaneously dragging near or beyond a window's edge. Common implementations rely on rate control, mapping the distance between the pointer and the edge of the viewport to the scrolling velocity. While ubiquitous in operating systems, edge-scrolling has received little attention, even though previous works suggest that (1) rate control may be suboptimal for isotonic pointing devices like mice and trackpads and (2) space beyond the window's edge might be scarce, limiting scrolling control. To address these problems, we developed Push-edge scrolling (and Slide-edge scrolling, its inertial variant), two novel position-based techniques that allow scrolling by 'pushing' against the viewport edge. A controlled experiment shows that our techniques reduce overshoots and offer performance improvements by up to 13% over traditional edge-scrolling.

INTRODUCTION

Scrolling while dragging an object is a common interaction problem when working within the confines of a viewport. Dragging is achieved by moving an input device in a special Permission to make digital or hard copies of all or part of this work for personal or classroom use is granted without fee provided that copies are not made or distributed for profit or commercial advantage and that copies bear this notice and the full citation on the first page. Copyrights for components of this work owned by others than ACM must be honored. Abstracting with credit is permitted. To copy otherwise, or republish, to post on servers or to redistribute to lists, requires prior specific permission and/or a fee. Request permissions from permissions@acm.org. CHI 2015, April 18 -23 2015, Seoul, Republic of Korea Copyright is held by the owner/author(s). Publication rights licensed to ACM. ACM 978-1-4503-3145-6/15/04...$15.00 http://dx.doi.org/10.1145/2702123.2702132 state (e.g., holding down a mouse button) -which allows actions such as moving file icons and selecting regions of text. Dragging tasks are straightforward when the start and end locations are both contained within the viewport, but many activities require scrolling to bring the end location into view -for example, selecting several paragraphs of text, or manipulating data in a large spreadsheet. The dragging state precludes using scrollbars or other on-screen navigation tools, and device-specific mechanisms, like the mouse wheel, are not always responsive during dragging.

To allow scrolling while dragging, interfaces support edgescroll (also called autoscroll) that triggers a scrolling motion when a dragged pointer enters a control area near the edge of the viewport. For example, in Figure 1a the user first presses the trackpad button with the thumb, then uses another finger to drag the pointer downwards and select a text region. When the pointer enters the control area (typically the region between the window border and the screen edge), the document starts scrolling until the pointer leaves this area or the user exits the dragging mode. The most common approach to control edge scroll is rate-based: inside the control area, scrolling velocity is based on the pointer's distance from the viewport edge [START_REF]Macintosh human interface guidelines[END_REF].

While edge-scroll is commonplace and designs have been disclosed [START_REF] Kwatinetz | Scrolling contents of a window[END_REF][START_REF] Li | Scrolling a target window during a drag and drop operation[END_REF], we are unaware of any empirical work, which is surprising because current techniques suffer from two major limitations. First, prior research on pointing and scrolling suggests that rate control may be inappropriate for edge-scrolling [START_REF] Zhai | Human performance in six degree of freedom input control[END_REF][START_REF] Zhai | Improving browsing performance: A study of four input devices for scrolling and pointing tasks[END_REF]. Second, the size of the control area is inconsistent and possibly limited, for instance when the window is maximised or with small displays, offering suboptimal control.

Based on these limitations, we introduce push-edge scrolling (Figure 1b) and its inertial variant, slide-edge scrolling. These novel position-based techniques block the pointer as it reaches the edge of the viewport, and transform subsequent pointing device movements that would otherwise move the pointer into scrolling displacements. Experimental results demonstrate that push-edge and slide-edge scrolling outperform rate-based techniques in a text selection task.

LIMITATIONS IN CURRENT EDGE-SCROLL TECHNIQUES

Edge-scrolling should allow users to quickly scroll to close or distant target locations, but there are two main reasons to question the adequacy of current designs with respect to these goals: rate control may not be superior to position control with common input devices, and window placement can limit scrolling control because it constrains the control area.

Mouse and Trackpad Rate Control

Prior research suggests that rate-based methods for scrolling control are most desirable when the user's primary task is reading because they facilitate a smooth and continuous scrolling speed [START_REF] Zhai | Improving browsing performance: A study of four input devices for scrolling and pointing tasks[END_REF], or when input is provided through an isometric input device such as a trackpoint [START_REF] Zhai | Human performance in six degree of freedom input control[END_REF]. In contrast, isotonic devices, such as mouse or trackpads, are less suited to rate control because they lack a self-centring effect, requiring fine coordination to control movement and to explicitly get back to zero point to stop it [START_REF] Kim | A comparison of position and rate control for telemanipulations with consideration of manipulator system dynamics[END_REF]. In summary, there is support for exploring position-based alternatives to rate-based edge-scroll. Such a mechanism for normal scrolling was patented [START_REF] Hinckley | Positional scrolling[END_REF], scrolling a view depending on the distance between the cursor and an activation point. However, the distance-based nature of its transfer function makes it hard to precisely control scrolling velocity. Moreover, scrolling has to be explicitly activated by pressing a button, unlike edge-scrolling which is as a direct extension of the dragging action.

Control Area Size and Location

Current edge-scroll methods rely on the availability of a control area outside the viewport. Scrolling begins when a dragged pointer enters this area, and the scroll velocity increases as the pointer moves away from the viewport's edge. Therefore, the actual size of the control area limits the maximum scrolling velocity in rate control. If the control area lies inside the viewport, it has to be small enough to avoid unintended scrolling, at the expense of the granularity of control. A common alternative is to expand the control area from the viewport edge so that it takes up the whole display space beyond. However, available screen space is likely scarce when maximised windows or small displays, making the control area smaller. Either way, a small control area is likely to make rate-based edge scrolling hard to control for the user.

When using position control, it is possible to overcome this problem by blocking the pointer when it crosses the edge [START_REF] Accot | More than dotting the i's -foundations for crossing-based interfaces[END_REF] so that every subsequent movement will scroll the view instead of moving the pointer, effectively having a zero-sized control area that does not impede scrolling control.

PUSH-EDGE AND SLIDE-EDGE SCROLLING

Push-edge scrolling is designed to overcome the above limitations and permit rapid and accurate scrolling while a dragging state is maintained. It operates by 'capturing' the cursor at the viewport edge, and mapping subsequent device signals that would otherwise move the cursor as input to a scrolling transfer function (Figure 1b) -further movement away from the viewport results in scrolling displacements of the view, while the cursor remains captured in a fixed position at the edge. The captured cursor is released either by moving it back towards the viewport, by releasing the dragging state, or when scrolling reached the end of the document.

However, the range of a pointing movement -from a pixel to a screen -is much narrower than the range of a scrolling movement -from a few lines to hundreds of pages [START_REF] Cockburn | Improving scrolling devices with document length dependent gain[END_REF][START_REF] Hinckley | Quantitative analysis of scrolling techniques[END_REF]. Thus, using raw pointing device displacements only might not scale to long distance scrolling. Consequently, push-edge scrolling converts motor-space device displacements d (in metres) into display-space scrolling displacements D (in pixels) by multiplying d with a unitless gain produced by a sigmoid transfer function [START_REF] Roussel | Giving a hand to the eyes: Leveraging input accuracy for subpixel interaction[END_REF] and converting the result in pixels using the pixel density of the display. The transfer function transforms the device velocity v (in metres per second) using a linear interpolation between constants g min and g max in the interval (v 1 , v 2 ), and is clamped elsewhere.

Since the physical range of motion is different between devices, we used the results of a pilot study to calibrate the transfer function differently for a trackpad (g min = 1, g max = 15, v 1 = 0.2 m/s, v 2 = 0.5 m/s) and mouse (g min = 1, g max = 10,

v 1 = 0.1 m/s, v 2 = 0.5 m/s).
Because the transfer function is not rate-based, it requires repeatedly lifting and repositioning the mouse or finger (clutching) to scroll very long distances. To reduce clutching-induced effort, we also designed slide-edge scrolling, a variant of pushedge that applies inertia to continue scrolling after the user stops operating the pointing device based on a simulation of residual momentum [START_REF] Quinn | Touch scrolling transfer functions[END_REF]. When the raw device velocity at liftoff v l exceeds 0.15 m/s, friction progressively decreases velocity V over time according to the function V (t) = V l ⇥ e 4t -where t is the time elapsed since liftoff, in seconds.

EXPERIMENT

We conducted an experiment to compare performance, amount of control, and perceived workload between push-edge (Push), slide-edge (Slide), and a traditional rate-based edge-scroll (Rate) using a mouse and a trackpad in a one-dimensional text selection task. The primary hypothesis (H1) was that selection time would be lower for Push and Slide than for Rate; secondary hypotheses are: (H2) inertia (Slide) reduces selection time, (H3) Push and Slide offer better control and result in less overshoots than Rate, and (H4) inertia (Slide) reduces physical effort in long selections.

Method

The experiment was conducted on a 13 00 Apple MacBook Pro Retina running OS X 10.9.4, with display resolution set to Best (1280⇥800 px, 116 ppi). Input was provided through the laptop's trackpad and a Logitech M90 optical mouse on a plywood desk. The pointing transfer function was set to the fourth tick of the Tracking Speed slider in the respective touchpad and mouse preference panels. Raw device input was used for Push and Slide. Experimental software was written in Objective-C with the Cocoa API. Raw device input was monitored using the I/O Kit and Apple's private multitouch API. The Quartz Event Services API was used to lock the pointer. The software implemented push-edge and slide-edge scrolling, and it used Apple's default rate-based edge-scroll method for the Rate condition. It is defined in the autoscroll: method of the NSClipView class -which scrolls the viewport p pixels every 100 ms, where p is the pointer-to-edge distance.

Procedure, task, and design

Participants were instructed to perform a sequence of top-tobottom text selections as quickly and accurately as possible. The window was displayed at the centre of the screen, 400 pixels tall, containing 21,480 lines of text typeset with Ubuntu Mono Regular 13px (Figure 2). For each trial, participants had to select a section of text framed and coloured blue -starting eight lines from the bottom of the viewport, and varying in size. Supposing that edge-scroll is mainly target-directed, we overlaid a gradient on the scrollbar as a hint of the number of lines participants had to select. As the task was one-dimensional, selecting anywhere on a line selected the entire line. Every selection required concurrently dragging and scrolling using the requested device and edge-scroll technique. Other scrolling methods and devices were disabled.

have things; petticoats, starched ones, shoes, too, real jaunty ones to show off her foot when she has to step over a puddle. Do you understand, sir, do you understand what all that smartness means? And here I, her own father, here I took thirty copecks of that money for a drink! And I am drinking it! And I have already drunk it! Come, who will have pity on a man like me, eh? Are you sorry for me, sir, or not? Tell me, sir, are you sorry or not? He-he-he!" He would have filled his glass, but there was no drink left. The pot was empty. "What are you to be pitied for?" shouted the tavern-keeper who was again near them. Shouts of laughter and even oaths followed. The laughter and the oaths came from those who were listening and also from those who had heard nothing but were simply looking at the figure of the discharged government clerk. To perform a selection, participants had to move the pointer to the top of the blue text, press the pointing device button to start the selection (and enter a dragging mode), and scroll downwards until they reached the bottom of the blue text. Once the pointer was positioned over the last line, they could release the button to complete the selection. If their selection did not exactly capture the blue text, an error was recorded, and the trial was repeated.

Select all the text in the frame

The experiment used a 2 ⇥ 3 ⇥ 3 ⇥ 4 within-subjects design for the factors: device (mouse and trackpad), technique (Push, Slide, and Rate), block (1-3, with the first block serving as opportunity for learning the new methods), and size (shortest: 15; short: 45; long: 135; longest: 405 lines of text to be selected). The order of technique and device was counterbalanced. Sizes were presented from the shortest to the longest, one for each block, with 5 consecutive repetitions for a given size within a block -for a total of 360 trials per participant. Participants completed NASA-TLX worksheets after each technique. The experiment lasted ~45 minutes.

Participants

Twenty-four university staff and students (five female) participated in the experiment (mean age 30.4, s = 7.0). Most used a trackpad as their primary pointing device (58%).

Results

Results for each of the dependent variables (trial time, overshoot distance, and number of clutches) are presented below. Error rate was not significantly different (p=.313) between Push (6.0%), Slide (5.7%), and Rate (6.9%). We excluded error trials from analyses. We also excluded the first repetition of each size, as changes were not immediately noticed.

Trial time

Trial time is the total text selection time, from the first mouse movement after the button was pressed, to the button release.

Repeated-measures ANOVA1 revealed a significant effect of block on trial time (F 2,46 = 28.0, p < .001, h 2 p = .55; block 1: 4.0s, 2: 3.7s, 3: 3.6s), with a significant technique⇥block interaction (F 4,92 = 3.0, p=.023, h 2 p = .11) due to a learning of the new technique behaviour during the first block. As we are concerned with user performance after familiarization, the remaining analysis discards the first block.

There were significant main effects of device (F 1,23 = 28.7, p < .001, h 2 p = .56), technique (F 2,46 = 6.9, p=.002, h 2 p = .23), and size (F 1.1,25.6 = 360.8, p < .001, h 2 p = .94), and a significant technique⇥size interaction (F 2.4,55.3 = 4.2, p=.014, h 2 p = .16) on trial time. Post-hoc tests revealed that trial time was significantly lower for Slide than for Push (p=.009) or Rate (p=.006) with the shortest size, and than for Rate with the short size (p < .004). Trial time was also found significantly higher with Rate than with both Push and Slide for the long (respectively, p=.014 and p=.004) and longest (respectively, p=.029 and p=.003) sizes.

Importantly, Push and Slide consistently outperformed Rate by up to 13% in the longest size (Figure 3, left), supporting H1. There was no significant device⇥technique interaction (p=.405), suggesting that the benefits of Push are consistent across devices. Finally, inertia had no effect on performance as post-hoc tests revealed no significant difference between Push and Slide (p=.405), leading to rejection of H2.

Overshoot distance

Overshoot distance was measured as the maximum scroll distance beyond the target bounds that was reached during the trial. As expected, ANOVA confirmed that overshoot distance increases with size (F 1.4,31.0 = 106.1, p < .001, h 2 p = .82). There was also a significant effect of Technique (F 2,46 = 66.7, p < .001, h 2 p = .74), with Rate having the largest overshoot distance (369 px) and Push the least (262 px; Slide: 279 px); supporting H3. The increased overshooting with Rate likely stems from the need for the user to return the cursor to a location inside the window before scrolling stops after they notice the target has come into view.

Physical effort: clutching and perceived workload

We analysed the number of clutches used on the trackpad, assuming that frequent clutching indicates high physical workload. ANOVA showed significant main effects of technique (F 2,46 = 47.3, p < .001, h 2 p = .67) and size (F 1.3,30.7 = 172.4, p < .001, h 2 p = .88), as well as a significant technique⇥size interaction (F 2.2,50.2 = 28.2, p < .001, h 2 p = .55). As Rate allows continuous scrolling to be maintained without movement, it is unsurprising that it has significantly less clutching than Push. Pairwise comparisons between Slide and Push confirmed that the addition of inertia decreases clutching (Figure 3, right). NASA-TLX responses showed physical demand to be lower with Rate than Push or Slide, for both the mouse (Friedman c 2 (2) = 17.2, p < .001) and trackpad (c 2 (2) = 11.9, p=.003). We also found that for the mouse, Slide had significantly higher perceived effort than Rate (c 2 (2) = 7.7, p=.022), which probably stems from the difficulty of lifting the mouse while maintaining a high tangential velocity. Therefore, clutching and subjective data only partially support H4.

DISCUSSION

Our key finding is that push-edge and slide-edge scrolling both improved text selection time over a rate-based edge-scroll technique. The performance benefits were consistent for both trackpad and mouse pointing devices. Workload measures suggest that push-edge scrolling is best combined with mouse input, and slide-edge scrolling best with a trackpad.

It is likely that the actual performance benefits for push-edge and slide-edge scrolling over rate-based methods are understated with our results. Our experiment intentionally used a relatively small window centered on the screen, which is a 'best case' scenario for traditional techniques because of the large control area outside the window. When a window extends to the edge of the display (e.g., when maximised) the control area is compromised or unavailable, substantially impairing traditional methods -but there is no reason to anticipate any performance detriment in this case with our techniques.

The target of an edge-scrolling action while selecting text is consequently located in the scrolled viewport. However, there are situations where the scope of a dragging action may span multiple viewports (e.g. dragging a file between two windows) where moving the pointer beyond the viewport edge could be interpreted either as the continuation of the dragging action or as a request to activate edge-scroll. Current systems use various methods to distinguish these actions, such as waiting for a pointer to dwell over an interior control area in order to start scrolling. The relative merits of each approach, as well as their applicability to our techniques are left for future work.

Finally, one potential usability issue of our techniques is that they temporarily remove the user's direct control of the cursor as it is 'captured' by the window edge. The immediate feedback of seeing the viewport scroll should mitigate the user's surprise, but not if it reaches its terminus. To help visually communicate the captured cursor state, the techniques could apply a visual resistance effect to the viewport.

CONCLUSION

Edge-scroll is extensively used to enable scrolling while a concurrent dragging state is maintained. Most existing designs use rate control, where cursor position into a control area around the window edge controls scrolling velocity. We showed the potential performance limitations of current techniques and designed push-edge and slide-edge scrolling to work around them by allowing to scroll by "pushing against the viewport edge". Experimental results demonstrate that push-edge and slide-edge scrolling significantly improve performance.

CHAPTER 2 H 1 1 Figure 1 .

 2111 Figure 1. When selecting text with a touchpad, downward movements after crossing the viewport edge will (a) change the rate of automatic scrolling with existing techniques (rate control); or (b) manually scroll the document, stopping the pointer at the edge with push-edge scrolling (position control).

Technique 1 -

 1 165 selections left Technique 1 -165 selections left

Figure 2 .

 2 Figure 2. Experimental software for the text selection task.

Figure 3 .

 3 Figure 3. Mean trial time (left) and number of trackpad clutches (right) by size and technique (error bars denote 95% CI).

Greenhouse-Geisser corrections to the degrees of freedom were applied when sphericity was violated. Pairwise comparisons used Bonferroni correction.

ACKNOWLEDGMENTS

This work was supported by ANR (TurboTouch, ANR-14-CE24-0009). Thanks to UCLIC members, F. Chevalier, G. Bailly, T. Pietrzak and A. Goguey for their helpful comments.