BACKGROUND
investigating the molecular footprint of recombination in viral gene sequences is a multifaceted discipline.
such studies encompass statistical testing for the occurrence of recombination , inferring the rate of recombination in a population, identifying parental and recombinant sequences , and mapping breakpoints in mosaic genomes .
various recombination detection methods have been developed that generally focus on one or a subset of these tasks  <cit>  or attempt to address increasing levels of complexity .
most popular methods to investigate viral recombination patterns in a sequence alignment use a graphical sliding window approach to either scan the similarity of a query sequence against a set of parental  sequences, e.g.
rip  <cit>  and simplot  <cit> , or to compare the phylogenetic branching patterns of the query sequence relative to the parental sequences along the genome, e.g.
bootscan  <cit> .
however, martin et al.
noted that the proper use of these methods often crucially depends on the prior identification of a suitable set of non-recombinant parental sequences , and, indeed, inclusion of recombinant reference sequences can result in flawed interpretation.
in hiv- <dig>  for example, subtype g has been classified as a pure subtype for several years, but detailed analyses have recently shown that this strain most likely evolved following recombination between one lineage classifiable as the so-called circulating recombinant form crf02_ag and another belonging to the 'pure'  lineage subtype j  <cit> .
crf02_ag on the other hand, was found to be a 'pure' subtype and its classification as a circulating recombinant form  could be attributed to the inclusion of subtype g as a reference sequence in the analysis.
in this case, the sampling history of subtypes and crfs has caused a misinterpretation of the evolutionary history of hiv- <dig> group m: the parental subtype j complete genome was not available at the time subtype g was completely sequenced, so a recombinant pattern for subtype g could not be clearly established   <cit> .
therefore, subtype g was assigned as a reference strain in the subsequent analysis of crf02_ag, resulting in a misleading bootscan profile that suggested a mosaic genome for crf02_ag  <cit> .
the crf02_ag example demonstrates that the 'query versus reference' approach will only be suitable to identify recombinants if valid reference sequences can be assigned.
such prior specification is essential for many methods, including dedicated hiv- <dig> methods and web servers .
an additional shortcoming of graphical sliding window based methods, like bootscan, is that they commonly lack a formal test for recombination.
although bootstrap support values are well studied in phylogenetics, there exists no clear statistical basis to conclude significance evidence for recombination based on bootstrap variability across genome regions.
in fact, by assessing significance conditional on reference sequences and crossover points that maximize the same test statistic, the bootscan approach falls into a sequential testing trap  <cit> .
to avoid the problem of assigning appropriate reference sequences, exploratory methods have been proposed that consider every sequence as a potential recombinant or parent in a scanning procedure.
recscan, for example, is a modified bootscan algorithm that checks all combinations of three sequences  for changing nearest neighbor relationships.
if different relationships are observed with bootstrap support above a user-specified value, a test of recombination is applied to the triplet  <cit> .
with appropriate multiple testing correction, a wide array of test approaches can be applied to all triplets to identify those combinations that show evidence for recombination, including also substitution distribution methods like maximum chi <dig> and geneconv  <cit> .
none of these triplet exploration procedures, however, identify the mosaic sequence within a triplet that shows evidence for recombination.
to tease apart recombinant from parental descendent sequences, heath et al.
proposed a series of additional heuristic tests that examine which sequence relationships change the most across a recombination breakpoint.
a weighted consensus of these tests is eventually responsible for calling the recombinant sequence.
unfortunately, this weighting scheme is an arbitrary choice that is currently based on how accurately hiv- <dig> crfs are identified as mosaic genomes.
given the problems concerning crfs and subtype reference sequences discussed above, approaches that do not require empirical decisions would be a useful alternative for recombinant identification.
other methods have been proposed for recombinant identification  <cit> , but accurate prior classification of sequences remains essential for their performance.
moreover, while the performance of different methods to detect the presence of recombination has been extensively evaluated  <cit> , no such efforts have been undertaken to investigate how well methods perform in teasing out recombinant sequences.
here, we present an alternative approach to detect recombinants without prior identification of non-recombinant reference strains.
we employ quartet-trees to rapidly scan for phylogenetic inhomogeneity along a sequence alignment and demonstrate how this information can also be employed to identify those sequences responsible for detectable recombination signals.
using simulated data sets we evaluate different quartet incongruence measures, different approaches to assess the significance of various recombination detection statistics, and compare their performance at identifying recombinant sequences with that of triplet-based methods.
finally, the usefulness of these methods in investigating recombination is demonstrated at different scales of primate immunodeficiency virus  evolution, including both hiv and siv.
methods
we describe the method developments in several subsections.
the procedure we present here is an extension of a previously developed visual recombination detection  method.
we propose a scanning method that evaluates whether all combinations of four sequences  in a sequence alignment jointly provide evidence for recombination.
when this is the case, we rank taxa or groups of taxa according to their contribution to this recombination signal.
the method can be employed to sequentially prune putative recombinants until no significant recombination evidence can be found in the alignment.
we start by briefly explaining the previously developed visual recombination detection method based on quartet scanning.
following the description of a novel quartet mapping approach based on a distance-based method, we present measures for phylogenetic inhomogeneity of quartets, ranking measures for taxa and groups of taxa, and a global test statistic for recombination.
we then describe how null distributions can be obtained for this test statistic.
in addition to the visrd method, we briefly explain alternative triplet approaches used in our comparisons.
we conclude the methods section by providing details on the simulated and empirical data used in this study.
visual recombination detection
the visrd method is designed to visually inspect sequence alignments for recombination events  <cit> .
visrd works by computing quartets, or unrooted phylogenetic trees on four leaves, for each possible four-taxa set.
more specifically, at each window  in an alignment, a support si for each of the three possible quartet topologies ti is computed using a statistical geometry approach  <cit> .
essentially, the support value si is computed by summing the number of site patterns that support ti.
following the quartet mapping approach  <cit>  , a relative support is then derived that is defined as si/.
subsequently, the relative supports for each quartet topology are summarized in a quartet-mapping triangle; in this triangle, each quartet  is represented as a point , whose co-ordinates reflect the support for each of the quartet topologies.
the corners of the triangle represent the three fully resolved topologies for four taxa.
a point located close to such a corner reflects a high support for this particular topology; the three relative support values are represented by the lengths of the perpendiculars from the point to the triangle sides  <cit> .
by considering all of the windows along the alignment, a trajectory in the triangle is then computed between the points representing consecutive alignment partitions, the rationale being that if recombination has occurred then points will follow highly variable trajectories.
these trajectories are then filtered to reduce noise, and depicted in a highway plot , which can be visually inspected for recombination signals.
in the highway plot, the horizontal axis represents the sites in the alignment, whereas the three lanes represent the regions supporting different topologies .
mapping points from a quartet triangle onto the vertical axis of the highway plot is achieved by the polar coordinates of the points; the most variable trajectories between all the points are filtered out to make up the highway plot  <cit> .
distance-based quartet mapping
as indicated above, the quartet mapping approach used in visrd is based on statistical geometry  <cit> .
although this eases computation , it does not take into account the complexities of nucleotide or amino acid substitution patterns.
to address this shortcoming in the new version of visrd, but still resort to fast quartet computation, we have implemented a distance-based minimum evolution  approach to quartet-mapping.
briefly, branch lengths for each quartet topology are obtained using least-squares solutions  <cit> , and the relative support for topology ti is given by , where me is the sum of all positive branch lengths.
for nucleotide sequences, the jukes-cantor, kimura-two-parameter, felsenstein ' <dig>  and tamura-nei substitution models were implemented.
for protein sequences, we implemented a protein distance estimator that uses a logarithmic correction of observed divergence based on the alignment score derived from empirical transition probability matrices  <cit> .
this estimator provides accurate and robust protein distances and can be easily adjusted to use any empirical log odds substitution matrix.
in addition to the blosum <dig> model, we also included hiv-specific matrices  <cit>  and constructed a matrix specifically for analyses of more divergent piv sequences.
for the latter, we downloaded the available amino acid alignments for the gag, pol, env and nef proteins from the hiv database , including  <dig>   <dig>   <dig> and  <dig> sequences respectively.
ambiguously aligned regions were deleted and phylogenetic trees were inferred using phyml  <cit> .
the stochastic amino acid substitution model was inferred using a maximum likelihood phylogenetic approach that estimates the  <dig> evolutionary rates, defining the general time reversible model of amino-acid substitution, jointly from a set of sequence alignments.
finally, a similarity matrix with an expected sequence dissimilarity of 38% was generated as previously described using hyphy  <cit> .
nucleotide sequence alignments were analyzed using the f <dig> model and a sliding window with a window size of  <dig> bp and a step size of  <dig> bp.
amino acid sequence data was analyzed using the piv log odds matrix, with a window size of  <dig> aa and a step size of  <dig> aa.
since nucleotide sequence alignments can include a number of ambiguous characters, particularly for a population of virus sequences, we also employed a random resolution of ambiguous sites to one of their possible states based on the nature of the ambiguity.
for example, a site denoted with the iubmb symbol 'm' could be either an 'a' or 'c' nucleotide; each is picked with equal probability when distances are computed using one of the above models.
taxon ranking, group ranking and test statistics
we now introduce methods for ranking taxa according to their phylogenetic variability, and we derive a test statistic for the presence of recombination in the global data set.
the purpose of the ranking is to identify the most likely recombinants in case the sequence alignment set contains significant evidence for recombination.
the ranking of taxa is achieved using measures of variability for quartets that contain these taxa.
the first measure of quartet trajectory variability computes a score for each quartet q by determining the average position of the corresponding point pi in the triangle across all alignment windows, , and then computing the maximum distance of any point from this position across all windows  = max|pi - |, figure 1).
this is repeated for all quartets, and, for each taxon t, a rank  is then computed as the average dmax value taken over all quartets that include t. in particular, defining qt as being the set of quartets that involve taxon t, and letting |qt| denote the size of qt, the rank rt is given by
we also introduce a grouping model as an extension to the four taxa model in visrd, allowing any number of groups, perhaps representing well-established virus subtypes or genotypes.
in particular, for a pre-specified group of taxa g, a group rank  is computed by taking the average dmax value over all taxa in group g, and over all quartets qt that include t, that is,
using the group rank , we can focus on recombinant events between pre-specified groups while alleviating some of the computational burden of the quartet scanning procedure.
the mean taxon ranking  and mean group ranking , which quantify the phylogenetic incongruence in the alignment, are then taken to be the average of these values taken over all taxa and all groups, respectively.
alternatively, a rank for each quartet can be calculated based on the total distance covered by the trajectory for each quartet ) or the total distance from each point to the average position of all points in the trajectory ).
by again computing the average dtot or dav value taken over all n quartets qt that include t, and averaging over a group, we can obtain a taxon ranking and group ranking  and  respectively).
the taxon and group ranking values  are computed using all quartet trajectories .
the ranking values can be used to assess the contribution to phylogenetic inhomogeneity of individual taxa and groups respectively , whereas the mean taxon and group ranking values  can be individually used as statistics for detecting recombination in a sequence alignment .
to identify all plausible recombinants in a sequence alignment using the quartet scanning procedure, we follow a sequential deletion procedure.
when significant evidence of recombination is detected using  or , the top-ranked sequence is considered as a putative recombinant and removed from the alignment.
we repeat this until no significant evidence of recombination remains in the data set.
because the absolute values of dmax, and hence rt and rg, are sensitive to the evolutionary models and topology evaluation criteria used, we report the ranking values as percentages of the highest ranking taxon or group.
significance
we evaluate three different ways to generate a null distribution  for the observed test statistics  and  derived from the taxon-ranking procedure:  monte carlo simulation ,  permutation, and  'redistribution' of alignment columns.
for the mc-simulation procedure, maximum likelihood phylogenetic trees are reconstructed using phyml  <cit> , employing the general time-reversible substitution model  and gamma-distributed rate variation among sites.
because a null distribution needs to be generated under the hypothesis of 'no recombination', we simulate replicate data sets using seq-gen  <cit>  down a single phylogeny reconstructed from the complete alignment.
the same models and parameters were used as obtained by the maximum likelihood tree inference.
since taxon ranking is based on topological incongruence, we also use a permutation procedure that homogenizes any phylogenetic incongruence along the genome .
this is achieved by random shuffling of the original alignment sites to produce replicate data sets.
to accommodate varying phylogenetic signal along the genome, which is often observed in real data due to different degrees of conservation within the sequence, we also use a permutation procedure that mimics the variability distribution of the alignment columns .
essentially, this procedure reshuffles alignment columns according to their shannon entropy score  <cit> .
replicate data sets are constructed by randomly drawing for each alignment position an alignment column with the same entropy score.
we use these three procedures to generate  <dig> replicate data sets each, so as to compute null distributions for our test statistics .
a monte carlo p-value is then estimated by counting the proportion of times the  or  statistic on a replicate data set is larger than the value of  or  observed for the original alignment.
if the null hypothesis is true, this proportion is expected to be  <dig> .
such a p-value represents the probability of obtaining a result at least as extreme as observed for the original alignment in the absence of recombination.
rdp <dig> recombinant detection
the rdp <dig> software provides access to multiple recombination signal detection methods that can be used in conjunction with one another  <cit> .
the methods  in rdp <dig> search for recombination signals within sequence triplets sampled from an input alignment.
when a recombination signal is detected, the member of the sequence triplet that is the recombinant is identified using a weighted consensus of a series of fifteen heuristic tests .
these tests fall into four main categories:  phylogenetic profile based ;  branch pruning and re-grafting based ;  multiple recombination signal pattern analysis based ; and  simple distance-based metrics  that, for example, compare pair-wise genetic, and phylogenetic tree distances of different subsets of sequences within the alignment on opposite sides of the detected recombination breakpoint position.
these latter tests collectively identify the recombinant by determining which sequence within an alignment have relationships that change the most across a recombination breakpoint.
the weighting scheme used in rdp <dig> to combine the results of all  <dig> tests into a single consensus score is based on how accurately the individual tests identified hiv-1m inter-subtype recombinants  from amongst a background of supposedly non-recombinant hiv-1m subtype sequences.
simulation studies
two different simulation studies were undertaken:  to evaluate how accurately the methods detected specific recombinant sequences, and  to evaluate the overall power of methods at detecting recombination signals within sequence alignments .
for , sequence data sets were simulated using seq-gen  <cit>  according to different trees for different data partitions .
data sets were simulated using both a series of symmetric and asymmetric trees.
in the symmetric setup, taxon  <dig> is a relatively recent recombinant with a single breakpoint, while taxa  <dig> and  <dig> share a recombination event deeper in the tree resulting in two breakpoints for each sequence.
in the asymmetric setup, taxa  <dig> and  <dig> share a relatively recent recombination event with two breakpoints, while taxon  <dig> represents an older recombinant with a single breakpoint.
a thousand nucleotides were evolved on each of these trees using the hasegawa-kishino-yano nucleotide substitution model with gamma-distributed rate variation , resulting in a total sequence length of  <dig> nucleotides.
to simulate data sets with different degrees of diversity, we employed total tree depths of  <dig> ,  <dig>  and  <dig>  substitutions per site .
in addition to rate variation among sites, we also investigated the impact of rate variation among different alignment regions by employing different relative rates for each  <dig> nucleotide partition  and different relative rates for  <dig> nucleotide partitions .
for each treatment ,  <dig> data sets were simulated.
for , to evaluate power and false positive rates in detection of recombination signals, we examined simulated data sets previously generated for comparing the performance of  <dig> recombination detection methods  <cit> .
these data sets are now often used to evaluate recombination detection power and false positives of new recombination detection programs  <cit> .
empirical data
the hiv- <dig> group m alignments analyzed here were derived from a full genome alignment provided by the hiv database  and manually edited using se-al .
the accession numbers of the included sequences are listed in additional file  <dig>  two data sets were used in our analyses: a crf03_ab data set  and a larger crf02_ag data set .
the hiv- <dig> group m alignments had nucleotide diversities ranging from  <dig>  substitutions/site to  <dig>  substitutions/site, which roughly corresponds to θ =  <dig> in the coalescent-based simulations.
the siv data set, previously analyzed for recombination by bailes et al.
, includes protein sequences representative of eight major siv clades.
their accession numbers are also listed in additional file  <dig>
bootstrapped phylogenetic trees were reconstructed using phyml  <cit> , as described above  and network analyses were performed using the neighbor-net method  <cit>  implemented in splitstree  <cit> .
to represent robust network-like relationships, the neighbor-nets only display the splits that were present in 75% of  <dig> bootstrapped replicates.
intrasubtype recombination was analyzed using gard method, available at  <cit> .
tree topology tests were performed using tree-puzzle  <cit>  and consel  <cit> .
RESULTS
we present the results of three main analysis sections:  simulations to evaluate the performance of different statistics and measures of significance for the quartet scanning method, and simulation studies to compare quartet scanning with other recombination detection methods in terms of both detecting recombination signals and identifying recombinant sequences;  an analysis of hiv- <dig> group m inter- and intra-subtype recombination; and  recombinant identification using amino acid data from more divergent piv lineages.
simulated data analysis
to establish the most powerful test statistic and significance assessment for the distance-based quartet scanning procedure, we analyzed data sets simulated on both symmetric and asymmetric trees, which generates three recombinant sequences per data set .
to restrict the amount of computation in these simulation analyses, we focus on the first sequence ranked on top by the quartet scanning approach and perform sequential deletion whilst significant evidence for recombination can be detected.
the results are summarized in figure  <dig>  dashed lines represent the frequency by which a recombinant is correctly ranked on top, independent of whether significant evidence for recombination can be detected in the replicate set.
the full lines represent the same frequency, but now supported with significant evidence for recombination in the replicate set.
the upper graphs summarize the results for symmetric trees , while the lower graphs summarize the results for asymmetric trees ; the graphs on the right present the results for data sets simulated with different relative rates correlated with the alignment partitions .
overall, recombination detection was more successful in symmetric trees and a correctly identified recombinant was usually a 'recent' recombinant .
the total trajectory distance  appeared to be a poor test statistic for recombination detection; the total distance from the mean position in the quartet triangles  was most powerful in symmetric trees but the maximum distance from the mean position  appeared to be the best overall statistic.
we thus chose dmax as test statistic for further simulation analyses.
increasing sequence diversity resulted in increased power of recombinant detection in the symmetric trees .
somewhat surprisingly, this was usually not the case for analysis of asymmetric trees.
in general, simulation led to stronger statistical support than permutation.
permutation, in turn, was more powerful than redistribution.
redistribution also resulted in a markedly lower number of correctly identified recombinants, with statistical support for recombination in the replicate data set, in the presence of rate variation among alignment partitions.
the fact that the partitions in between recombination breakpoints coincide with the relative rate partitions in our simulation may be at least partly responsible for this, as sites may not be efficiently permutated across breakpoints in the redistribution procedure.
to assess the performance of the distance-based quartet scanning method as a statistical test for recombination signal, rather than a method for recombinant identification, we also reanalyzed nucleotide data simulated under varying levels of divergence and recombination originally presented in  <cit>  and compared this with other methods .
in general, these analyses revealed that the performance of our test is comparable to similar methods, like recscan  <cit> , but these are less powerful than the maximum chi <dig> method  <cit> , which has been reported as one of the most powerful nonparametric recombination detection methods.
finally, we compared the quartet scanning with triplet-based approaches and their associated heuristic tests in their ability to identify the mosaic sequences.
because this requires a sequential deletion procedure in visrd, we restricted the computation by focusing on the data sets simulated using symmetric trees  and a tree depth of  <dig>  substitution per site, similar to the hiv- <dig> group m data sets we analyze below.
figure  <dig> summarizes the frequencies with which taxa  <dig>   <dig> and  <dig> were ranked on top at some stage in the sequential deletion procedure, and hence identified as recombinants, by the different approaches, as well as the total number false positives  in the absence and presence of rate variation among sequence regions .
the performance of the different methods is evaluated using the adjusted rand index .
the adjusted rand index is the rand index corrected for chance events, which is commonly used to measure the performance of data clustering.
in our case, the adjusted rand index reflects the percentage of correctly identified recombinant and non-recombinant sequences corrected for chance.
therefore, higher values of this measure indicate higher classification accuracy.
the methods mainly differ in the frequencies with which taxa  <dig> and  <dig> – the sequences descended from a common recombinant ancestor – are identified as recombinants.
without rate variation among partitions, quartet scanning based on simulation achieves the best rate of detection of these recombinants and obtains the highest performance score.
using permutation and redistribution, top-ranked recombinants were still relatively frequently associated with significant evidence for recombination in the alignment, but the lower detection rate of taxon  <dig> indicates some sensitivity decrease.
the presence of rate variation among partitions increases the false positive rate in many methods ; quartet scanning with simulation is particularly sensitive to this.
a permutation procedure, however, still achieves good detection rates without an excessive number of false positives.
as also suggested by the first set of simulation analyses, redistribution becomes less powerful in detecting recombination signal with the simulated rate variation scheme.
to investigate whether this results from inefficient permutation across breakpoints correlated with rate variation, we also performed additional simulations using a rate variation scheme that does not coincide with recombination breakpoints .
also in this case, redistribution did not prove to be a very powerful approach to assess significance.
we therefore restrict significance assessment to simulation and permutation in subsequent analyses of real sequence data.
hiv- <dig> intersubtype recombination
we applied the visrd based recombination detection methods to the hiv- <dig> group m full coding genome data set which included representative sequences for subtypes a-d, and f and the circulating recombinant form crf03_ab .
the crf03_ab sequences have been clearly identified as having descended from a common ab recombinant form that circulated amongst injecting drug users in russian and ukrainian cities  <cit> .
in the quartet scanning procedure, the sequences were analyzed by grouping them according to their subtype or crf assignment.
this group model only considers quartets composed of sequences from four different pre-specified groups, which allows the analysis to focus exclusively on inter-subtype recombination events while reducing the computational burden.
in figure  <dig>  a highway plot provides a visual summary of the most variable trajectories obtained during the quartet scanning procedure.
the three lanes in the plot represent the different topology-supporting regions in a quartet triangle , and the change in support is calculated as a function of the alignment position for each quartet.
the highway plot clearly indicates two breakpoints at approximately nucleotide positions  <dig> and  <dig> .
table  <dig> lists the results of the taxon-ranking procedure for the different groups in this data set, indicating that quartets containing a crf03_ab sequence show the highest variability in phylogenetic clustering along the genome based on both dmax and dav.
evidence for recombination in the data set is statistically supported by both simulation and permutation .
the parental subtype a and b sequences are ranked just below the recombinant group due to the high trajectory variability of quartets including one or both parental subtypes and the recombinants.
even after deletion of the crf03_ab sequences, there was still significant evidence for recombination; a more in depth investigation of hiv- <dig> group m recombination is provided below.
the group ranking values are listed as percentages of the highest ranked group; these values are shown for the maximum distance  and total distance  from the mean position in the quartet-mapping triangle.
the number of sequences per group is indicated between brackets
the crf03_ab example indicates that the  values allow us to test for the presence of recombination while taxon ranking itself indicates the most likely recombinants in the data set without assigning query and parental sequences a priori.
together, this represents an ideal approach to evaluate recombination patterns in cases where prior identification of non-recombinant reference strains is problematic, like crf02_ag  <cit> .
to investigate the crf02_ag issue and other putative recombination events in hiv- <dig> group m, we compiled a full genome data set that is more representative of hiv- <dig> group m diversity.
also for the analysis of this data set, we employed a grouping model for the different subtypes and crfs.
although the  value for the data set including crf02_ag indicates the presence of recombination , the taxon ranking does not suggest crf02_ag sequences as most plausible mosaic genomes .
in agreement with the recent analysis of abecasis et al.
, the 'pure' subtype g sequences are proposed as the most likely recombinants while crf02_ag and subtype j are suggested as parental lineages because they were ranked right below the putative recombinant.
when the subtype g sequences were removed, recombination in the remaining data could not be excluded based on the permutation statistics .
this time, subtype k was at the top of the taxon ranking .
it is interesting to note that a recombinant origin for subtype k has been proposed  <cit> , but its designation was not changed since it was already an established subtype in the literature  <cit> .
after excluding both subtype g and k, there was still significant evidence for recombination in the remaining sequences  with subtype a heading the taxon ranking .
to our knowledge, subtype a has not been proposed as a intersubtype recombinant before and the recombination signal we detected might result from intra-subtype recombination or intersubtype recombination for particular lineages within subtype a. we investigated this in more detail by removing the grouping model for the different a sub-subtypes and by investigating the evolutionary relationships using neighbor-nets  <cit>  .
when no grouping for subtype a was used, sub-subtypes a <dig>  a <dig> and a <dig> were ranked on top but not sub-subtype a <dig> .
the neighbor-net indeed demonstrates network-like behaviour among these lineages within subtype a .
having removed sub-subtype a <dig>  a <dig> and a <dig>  our recombination test was still significant  and subtype j was now at the top of the ranking .
although recombination has not been clearly suggested for this subtype, it does indeed cluster differently in each major gene tree of the hiv genome .
further removing of subtype j, there was still significant recombination signal , this time with subtype d heading the ranking .
without the grouping model, three subtype d sequences appeared to be the descendents of a recombinant lineage , which was also confirmed by the neighbor-net analysis .
the final recombination signal was detected with subtype f heading the ranking, despite a stable clustering of this subtype in the major hiv genes encoding structural proteins .
recombination remained detectable after deletion of subtype f using the mc-simulation approach, but given the more likely confounding effect of rate variation for mc-simulation and the increasing difficulty of correctly identifying 'deeper' recombinants , we did not explore further taxa deletion.
a different order in putative recombinant deletion was suggested after identification of recombination in subtype g, k and subsubtypes a using dav as a test statistic, and subtype j was not suggested as a recombinant by this procedure .
the intrasubtype recombination detected in subtype a and d in both procedures was confirmed by gard analysis   <cit> .
we also employed our method to investigate a highly divergent hiv- <dig> variant from the democratic republic of congo   <cit> .
although the drc strain was positioned differently in phylogenetic trees of different genome regions , it did not cluster with any of the know subtypes and hence it was suggested as a candidate for a new subtype  <cit> .
as reference strains, we included the sequence set for which no significant recombination could be detected .
our taxon-ranking procedure indicated that quartets containing this unclassified variant show the highest phylogenetic variability along the genome , and this was associated with significant evidence for recombination in the alignment .
to confirm alternative clustering patterns for this lineage and to demonstrate significant incongruence despite convincing bootstrap support values, we reconstructed phylogenies for the major genes  and performed tree topology tests .
the trees indicate different clustering for the drc lineage in the three major genes, and despite low bootstrap support for the alternative clustering, tree topology tests consistently indicated significant incongruence.
finally, we also analyzed crf01_ae sequences, for which only the subtype a lineage is available as closely related parental descendent.
the taxon-ranking procedure still revealed significant recombinant signal and implicated this crf as the lineage responsible for the incongruence .
recombination in primate immunodeficiency virus evolution
recombination in more divergent simian immunodeficiency viruses  was investigated for the amino acid data set of bailes et al.
, which includes representatives of eight major siv lineages for which full-length sequences were available at the that time  <cit> .
because there are no equidistant siv clades similar to hiv- <dig> group m subtypes, we did not apply the grouping model in this case.
however, when sister taxa representing one of the eight major lineages were ranked as the top two sequences, we removed both taxa in the subsequent taxon-ranking test.
the protein sequences were analyzed using a corrected distance measure based on an amino acid substitution matrix inferred from a larger piv data set , and significance was assessed using permutation.
additional file  <dig> summarizes the results using the dmax statistic, while the results using dav are provided in additional file  <dig>
the ranking procedure provided evidence for recombination in the siv sequences .
moreover, in agreement with the findings of a hybrid origin of siv in chimpanzees  <cit> , sivcpz was suggested as contributing the most to this recombination signal .
however, recombination was also detected during subsequent deletion of several siv lineages , leaving only a single quartet that did not deviate significantly from clonality .
there were some differences in the order by which sivs were identified as recombinants using dav , but the final result was the same.
discussion
in this study, we have developed an extension of the visual recombination detection method to rank taxa or groups of taxa according to their contribution to phylogenetic inhomogeneity, and to test whether the overall phylogenetic incongruence has been significantly shaped by recombination.
during a sequence alignment scan, phylogenetic inhomogeneity is quantified using a quartet mapping approach.
taxa or groups of taxa are subsequently scored using the average inhomogeneity of all the quartets they belong to.
to assess whether significantly more inhomogeneity is observed than expected under a recombination-free scenario, we tested mc-simulation, permutation, and a redistribution approach.
simulation analyses suggested that identifying recombinant sequences without prior specifications is a challenging problem.
only relatively recent recombinants were readily identified in symmetric trees in which recombination was rare.
in asymmetric trees, correct recombination identification proved to be more difficult.
exploratory analyses using rdp <dig> indicated that this was also the case for the triplet-based methods .
however, the recent recombination event in the asymmetric trees occurred between somewhat more closely related parental strains, with a different mosaic structure, and involved two descendants of the same recombinant .
hence, symmetric and asymmetric tree simulations cannot be easily compared.
asymmetric topologies are generally more problematic in phylogenetic inference; more sequence data is required to recover such topologies and they are much more susceptible to substitution saturation than symmetrical ones  <cit> .
the latter may provide an explanation as to why increased sequence divergence did not result in increased power for recombination detection, as is the case for symmetric trees.
we found marked performance differences between different measures of quartet variability.
the poor performance of the total path distance  may be explained by the fact that the sliding window procedure follows the path between two regions supporting alternative topologies in the quartet-mapping triangle in a stepwise manner, with each step represented by different alignment windows.
steps covering similar distances could also be taken within a region supporting a particular topology, and thus lead to similar total distances in the absence of real topology changes.
no matter how many steps need to be taken to the point most distant from the average position, the dmax measure only considers how remote this point is from the average positions.
the power of the maximum distance or the average distance from the mean position will depend on the mosaic pattern of the recombinant.
dmax may be more sensitive to a small region with a different evolutionary history, while dav may be more sensitive in case of multiple breakpoints.
the simulation setup was necessarily simplistic and we have not explored the variety of recombination scenarios that occur in reality.
when recombination is frequent, for example, all the methods would be expected to suffer from an elevated failure rate  <cit> .
also the sampling rates of recombinants relative to the parental descendent sequences will have a profound impact on the ability to identify mosaic sequences without prior reference specifications.
when few representatives of the recombinant lineage are sampled and relatively closely related descendents of both parents are present, as is the case in our phylogenetic simulations, than the quartet variability should be attributed relatively easily to the recombinant sequences by averaging the dmax measure over the quartets that contain these sequences.
when, however, many sequences are included that share the same recombination history, more quartets will be composed of multiple recombinant sequences of the same type.
these quartets may not generate strong phylogenetic inhomogeneity, and in the case of all four sequences representing the same recombinant, there should be no recombination signal detectable.
although sampling rate inequalities for pre-defined lineages may be controlled to some extent by the grouping constraint we proposed, it will no doubt adversely affect taxa ranking.
this may also at least partly explain the poor performance for the asymmetric tree simulations.
here, the more recent recombinant event had two descendents, and frequently, a parental descendent sequence  was top ranked instead of one of the recombinants.
this again highlights that the current simulations only offer a starting point for comparison of different methods and their relative ability to identify recombinant sequences, and more complexity needs to be explored in the future.
for the quartet scanning method, generating null distributions using mc-simulation was the most powerful procedure for identifying evidence of recombination, but it also suffered from a higher false positive rate in the presence of rate variation among genome regions.
the simulated patterns of rate variation are within the range of rate variation that can be expected in hiv genomes.
in this case, permutation still achieved good detection rates with far fewer false positives.
this difference illustrates that the mc-simulation technique performs very well when the simulated evolution meets the assumptions in the parametric bootstrapping procedure.
under model misspecification, however, it may still be powerful but could also be particularly misleading.
although rate variation can be accommodated in the simulation procedure using a discrete gamma distribution, it assumes that rate variation is independent and identically distributed along the genome.
this assumption will rarely be met for real data situations and it can have an important impact on comparative analyses of data partitions across the genome.
the redistribution approach, which was designed to homogenize phylogenetic incongruence, but still mimic the pattern of rate variation present in the real data sets, was the least sensitive.
even when rate variation does not correlate with breakpoints, this approach does not easily lead to rejection of clonality.
this indicates that permuting sites only with the same entropy score might be too stringent to efficiently homogenize recombination signal.
in this respect, it is important to realize that a sequence alignment can include sites with many different entropy values, and the number of sites per entropy value will naturally show a large variation.
for a typical simulated data set comprising  <dig> taxa , site patterns can break down in about  <dig> different site entropy values with, for example, a few patterns having unique entropy values and about 13% of patterns being single nucleotide polymorphisms that are represented by a single entropy value.
therefore, more efficient shuffling under more relaxed rate variation constraints may be achieved by binning site patterns based on their entropy value.
taxon ranking based on quartet scanning is geared towards the detection of recombinant sequences.
as a test for recombination signal, more powerful methods have been devised , and we do not provide a formal test for breakpoint locations.
also for the latter, we refer to other developments .
however, automated detection of recombinant lineages on its own is already a difficult problem to tackle.
the quartet scanning approach performed well relative to various triplet-based methods in identifying recombinants simulated using a phylogenetic approach.
however, more comprehensive simulations are needed for a better assessment of the powers and pitfalls in different situations.
the essential difference between the quartet scanning approach and the various rdp <dig> methods does not lie in the use of quartets versus triplets, but how their information is further employed to identify the recombinants.
while the quartet approach deduces a ranking for taxa or groups of taxa from the phylogenetic variability in all these quartets, additional tests are used by rdp <dig> to distinguish recombinants from parental strains in triplets that resulted in significant evidence for recombination.
it is interesting to note that the weighting scheme used by rdp <dig> to combine the results of these tests into a single consensus score is currently based on how accurately the individual tests identified hiv-1m inter-subtype recombinants  <cit> .
because of the crf02_ag-subtype g classification artefact, a new weighting optimization might be required.
moreover, quartet scanning might make a useful addition to the tests employed by rdp <dig>
the quartet approach presented here has particular practical limitations.
the number of quartets increases as a fourth order polynomial in the number of taxa, and although the group constraint can alleviate some of the computational burden, both memory requirements and speed of computation will become impractical for very large data sets.
this will be a particular burden when many sequential deletions of recombinants are required and significance needs to be assessed each time using replicate data sets.
although these procedures can be automated and quartet information for a permutation alignment can be stored and reused, it would be very useful to be able to test every taxon or taxon group without having to rescan pruned data sets.
our experience with the current implementation leads us to conclude that  <dig> taxa is roughly the upper limit for this method in practice.
there are obviously less triplets than quartets for any data set, and rdp <dig> does not evaluate evidence for recombination in triplets by analyzing replicate data sets, which makes them faster in computational practice.
in sliding window analyses, it is important to determine how much phylogenetic discordance can be expected from the stochastic nature of the substitution process in the absence of recombination.
in bootscan analysis  <cit> , a query sequence is generally inferred as a recombinant when a high bootstrap support is observed for clustering with different reference sequences in different gene regions.
however, this criterion will not be adequate to assess variable clustering deeper in the phylogenetic trees.
this is illustrated by the highly divergent hiv- <dig> variant from the drc.
although phylogenetic inhomogeneity was originally noted for this strain, the absence of strong bootstrap support for the different clustering deeper in the tree led to the suggestion that it could be a new subtype  <cit> .
however, our test indicates that there is more phylogenetic inhomogeneity than expected by chance and thus classifies the drc strain as a probable recombinant lineage.
significant phylogenetic incongruence between major genes was also confirmed by various tree topology tests.
both subtype g in hiv- <dig> group m and sivcpz in piv are confirmatory examples, but the application of quartet scanning demonstrates how recombinants can be easily identified for cases that required in depth analysis to distinguish parental from recombinant lineages  <cit> .
the pervasive recombination in piv lineages has also been suggested by exploratory bootscanning and network analyses  <cit> .
in addition to piv recombination and hiv- <dig> intersubtype recombination, our test was also sensitive to intrasubtype recombination.
both in subtype a and d, we identified recombination signals that were supported by network and gard analyses, and in both cases, the recombination signal did not stem from a single mosaic sequence but from a recombination event that was shared among several sequences.
although we have identified recombinant genomes at different evolutionary scales in immunodeficiency virus evolution, it may prove a challenging endeavor to clearly delineate breakpoints for many of these recombinants.
similar to other sliding window analyses, breakpoints will be readily identified if viruses from two well-supported phylogenetic clusters recombined relatively recently.
however, our method may also pick up recombinants in trees without well-supported groups , ancestral recombinants that were generated before well-supported clades were established , or ancestral recombinants of lineages that have now diverged considerably .
CONCLUSIONS
in conclusion, we have developed a recombination detection method that does not require a priori assignment of reference and query sequences, statistically evaluates the recombination signal in nucleotide and amino acid alignments, and identifies the most plausible recombinants.
taxon ranking and the associated statistical test complement the visual recombination detection method.
collectively, these developments comprise a powerful and versatile methodology that will make a useful addition to the array of recombination detection methods.
abbreviations
crf: circulating recombinant form; piv: primate immunodeficiency virus; siv: simian immunodeficiency virus; visrd: visual recombination detectionl; rdp: recombination detection program; me: minimum evolution; f84: felsenstein ' <dig> model of evolution; gtr: general time-reversible substitution model
competing interests
the authors declare that they have no competing interests.
authors' contributions
pl conceived of the study, designed the analysis strategy, analyzed the data and drafted the paper.
ml carried out the visrd software implementation and supported data analyses.
dm wrote a modified rdp <dig> version to analyze simulated data and supported data analysis using triplet-based methods and its interpretation.
vm conceived of the study, participated in its design and coordination, and helped to draft the manuscript.
all authors read and approved the manuscript.
supplementary material
additional file 1
table a1: hiv and siv sequence data.
a table listing the accession numbers and subtype/crf assignment or piv lineage for the sequences analyzed in this study.
click here for file
additional file 2
coalescent simulations.
this file describes the coalescent simulations and presents the results for different visrd settings and for a comparison of visrd to other approaches.
click here for file
additional file 3
additional table a <dig>  this table lists the group ranking for quartet scanning of the crf02_ag data set using dmax.
click here for file
additional file 4
additional table a <dig>  this table lists the group ranking for quartet scanning of the crf02_ag data set using dav.
click here for file
additional file 5
additional table a <dig>  this table lists the group ranking based on quartet scanning of hiv- <dig> group m including a divergent drc lineage.
click here for file
additional file 6
tree reconstructions and topology tests for gag, pol and env data sets including the drc lineage.
this file includes a figure presenting the results of a phylogenetic analysis for gag, pol and env data sets including the drc lineage and a table listing tree topology test results for the same data.
click here for file
additional file 7
additional table a <dig>  this table lists the group ranking based on quartet scanning of hiv- <dig> group m including crf01_ae.
click here for file
additional file 8
additional table a <dig>  this table lists the taxon ranking for quartet scanning of the primate immunodeficiency viruses using dmax.
click here for file
additional file 9
additional table a <dig>  this table lists the taxon ranking for quartet scanning of the primate immunodeficiency viruses using dav.
click here for file
acknowledgements
we thank liz bailes for providing the siv sequence data and david posada for providing the simulated data sets and previous results for other recombination detection methods.
we thank joe parker for programming assistance.
pl was supported by a marie-curie and fwo  postdoctoral fellowship.
pl, ml and vm thank the isaac newton institute for mathematical sciences, cambridge, uk, for hosting them in the context of its phylogenetics programme where part of the research presented in this paper was carried out.

 investigating the molecular footprint of recombination in viral gene sequences is a multifaceted discipline. such studies encompass statistical testing for the occurrence of recombination, inferring the rate of recombination in a population, identifying parental and recombinant sequences, and mapping breakpoints in mosaic genomes. various recombination detection methods have been developed that generally focus on one or a subset of these tasks or attempt to address increasing levels of complexity. most popular methods to investigate viral recombination patterns in a sequence alignment use a graphical sliding window approach to either scan the similarity of a query sequence against a set of parental sequences, e.g. rip and simplot, or to compare the phylogenetic branching patterns of the query sequence relative to the parental sequences along the genome, e.g. bootscan. however, martin et al. noted that the proper use of these methods often crucially depends on the prior identification of a suitable set of non-recombinant parental sequences, and, indeed, inclusion of recombinant reference sequences can result in flawed interpretation. in hiv- for example, subtype g has been classified as a pure subtype for several years, but detailed analyses have recently shown that this strain most likely evolved following recombination between one lineage classifiable as the so-called circulating recombinant form crf02_ag and another belonging to the'pure' lineage subtype j. crf02_ag on the other hand, was found to be a'pure' subtype and its classification as a circulating recombinant form could be attributed to the inclusion of subtype g as a reference sequence in the analysis. in this case, the sampling history of subtypes and crfs has caused a misinterpretation of the evolutionary history of hiv- group m: the parental subtype j complete genome was not available at the time subtype g was completely sequenced, so a recombinant pattern for subtype g could not be clearly established. therefore, subtype g was assigned as a reference strain in the subsequent analysis of crf02_ag, resulting in a misleading bootscan profile that suggested a mosaic genome for crf02_ag. the crf02_ag example demonstrates that the'query versus reference' approach will only be suitable to identify recombinants if valid reference sequences can be assigned. such prior specification is essential for many methods, including dedicated hiv- methods and web servers. an additional shortcoming of graphical sliding window based methods, like bootscan, is that they commonly lack a formal test for recombination. although bootstrap support values are well studied in phylogenetics, there exists no clear statistical basis to conclude significance evidence for recombination based on bootstrap variability across genome regions. in fact, by assessing significance conditional on reference sequences and crossover points that maximize the same test statistic, the bootscan approach falls into a sequential testing trap. to avoid the problem of assigning appropriate reference sequences, exploratory methods have been proposed that consider every sequence as a potential recombinant or parent in a scanning procedure. recscan, for example, is a modified bootscan algorithm that checks all combinations of three sequences for changing nearest neighbor relationships. if different relationships are observed with bootstrap support above a user-specified value, a test of recombination is applied to the triplet. with appropriate multiple testing correction, a wide array of test approaches can be applied to all triplets to identify those combinations that show evidence for recombination, including also substitution distribution methods like maximum chi and geneconv. none of these triplet exploration procedures, however, identify the mosaic sequence within a triplet that shows evidence for recombination. to tease apart recombinant from parental descendent sequences, heath et al. proposed a series of additional heuristic tests that examine which sequence relationships change the most across a recombination breakpoint. a weighted consensus of these tests is eventually responsible for calling the recombinant sequence. unfortunately, this weighting scheme is an arbitrary choice that is currently based on how accurately hiv- crfs are identified as mosaic genomes. given the problems concerning crfs and subtype reference sequences discussed above, approaches that do not require empirical decisions would be a useful alternative for recombinant identification. other methods have been proposed for recombinant identification, but accurate prior classification of sequences remains essential for their performance. moreover, while the performance of different methods to detect the presence of recombination has been extensively evaluated, no such efforts have been undertaken to investigate how well methods perform in teasing out recombinant sequences. here, we present an alternative approach to detect recombinants without prior identification of non-recombinant reference strains. we employ quartet-trees to rapidly scan for phylogenetic inhomogeneity along a sequence alignment and demonstrate how this information can also be employed to identify those sequences responsible for detectable recombination signals. using simulated data sets we evaluate different quartet incongruence measures, different approaches to assess the significance of various recombination detection statistics, and compare their performance at identifying recombinant sequences with that of triplet-based methods. finally, the usefulness of these methods in investigating recombination is demonstrated at different scales of primate immunodeficiency virus evolution, including both hiv and siv. methods we describe the method developments in several subsections. the procedure we present here is an extension of a previously developed visual recombination detection method. we propose a scanning method that evaluates whether all combinations of four sequences in a sequence alignment jointly provide evidence for recombination. when this is the case, we rank taxa or groups of taxa according to their contribution to this recombination signal. the method can be employed to sequentially prune putative recombinants until no significant recombination evidence can be found in the alignment. we start by briefly explaining the previously developed visual recombination detection method based on quartet scanning. following the description of a novel quartet mapping approach based on a distance-based method, we present measures for phylogenetic inhomogeneity of quartets, ranking measures for taxa and groups of taxa, and a global test statistic for recombination. we then describe how null distributions can be obtained for this test statistic. in addition to the visrd method, we briefly explain alternative triplet approaches used in our comparisons. we conclude the methods section by providing details on the simulated and empirical data used in this study. visual recombination detection the visrd method is designed to visually inspect sequence alignments for recombination events. visrd works by computing quartets, or unrooted phylogenetic trees on four leaves, for each possible four-taxa set. more specifically, at each window in an alignment, a support si for each of the three possible quartet topologies ti is computed using a statistical geometry approach. essentially, the support value si is computed by summing the number of site patterns that support ti. following the quartet mapping approach, a relative support is then derived that is defined as si/. subsequently, the relative supports for each quartet topology are summarized in a quartet-mapping triangle; in this triangle, each quartet is represented as a point, whose co-ordinates reflect the support for each of the quartet topologies. the corners of the triangle represent the three fully resolved topologies for four taxa. a point located close to such a corner reflects a high support for this particular topology; the three relative support values are represented by the lengths of the perpendiculars from the point to the triangle sides. by considering all of the windows along the alignment, a trajectory in the triangle is then computed between the points representing consecutive alignment partitions, the rationale being that if recombination has occurred then points will follow highly variable trajectories. these trajectories are then filtered to reduce noise, and depicted in a highway plot, which can be visually inspected for recombination signals. in the highway plot, the horizontal axis represents the sites in the alignment, whereas the three lanes represent the regions supporting different topologies. mapping points from a quartet triangle onto the vertical axis of the highway plot is achieved by the polar coordinates of the points; the most variable trajectories between all the points are filtered out to make up the highway plot. distance-based quartet mapping as indicated above, the quartet mapping approach used in visrd is based on statistical geometry. although this eases computation, it does not take into account the complexities of nucleotide or amino acid substitution patterns. to address this shortcoming in the new version of visrd, but still resort to fast quartet computation, we have implemented a distance-based minimum evolution approach to quartet-mapping. briefly, branch lengths for each quartet topology are obtained using least-squares solutions, and the relative support for topology ti is given by, where me is the sum of all positive branch lengths. for nucleotide sequences, the jukes-cantor, kimura-two-parameter, felsenstein' and tamura-nei substitution models were implemented. for protein sequences, we implemented a protein distance estimator that uses a logarithmic correction of observed divergence based on the alignment score derived from empirical transition probability matrices. this estimator provides accurate and robust protein distances and can be easily adjusted to use any empirical log odds substitution matrix. in addition to the blosum model, we also included hiv-specific matrices and constructed a matrix specifically for analyses of more divergent piv sequences. for the latter, we downloaded the available amino acid alignments for the gag, pol, env and nef proteins from the hiv database, including and sequences respectively. ambiguously aligned regions were deleted and phylogenetic trees were inferred using phyml. the stochastic amino acid substitution model was inferred using a maximum likelihood phylogenetic approach that estimates the evolutionary rates, defining the general time reversible model of amino-acid substitution, jointly from a set of sequence alignments. finally, a similarity matrix with an expected sequence dissimilarity of 38% was generated as previously described using hyphy. nucleotide sequence alignments were analyzed using the f model and a sliding window with a window size of bp and a step size of bp. amino acid sequence data was analyzed using the piv log odds matrix, with a window size of aa and a step size of aa. since nucleotide sequence alignments can include a number of ambiguous characters, particularly for a population of virus sequences, we also employed a random resolution of ambiguous sites to one of their possible states based on the nature of the ambiguity. for example, a site denoted with the iubmb symbol'm' could be either an'a' or'c' nucleotide; each is picked with equal probability when distances are computed using one of the above models. taxon ranking, group ranking and test statistics we now introduce methods for ranking taxa according to their phylogenetic variability, and we derive a test statistic for the presence of recombination in the global data set. the purpose of the ranking is to identify the most likely recombinants in case the sequence alignment set contains significant evidence for recombination. the ranking of taxa is achieved using measures of variability for quartets that contain these taxa. the first measure of quartet trajectory variability computes a score for each quartet q by determining the average position of the corresponding point pi in the triangle across all alignment windows,, and then computing the maximum distance of any point from this position across all windows = max|pi - |, figure 1). this is repeated for all quartets, and, for each taxon t, a rank is then computed as the average dmax value taken over all quartets that include t. in particular, defining qt as being the set of quartets that involve taxon t, and letting |qt| denote the size of qt, the rank rt is given by we also introduce a grouping model as an extension to the four taxa model in visrd, allowing any number of groups, perhaps representing well-established virus subtypes or genotypes. in particular, for a pre-specified group of taxa g, a group rank is computed by taking the average dmax value over all taxa in group g, and over all quartets qt that include t, that is, using the group rank, we can focus on recombinant events between pre-specified groups while alleviating some of the computational burden of the quartet scanning procedure. the mean taxon ranking and mean group ranking, which quantify the phylogenetic incongruence in the alignment, are then taken to be the average of these values taken over all taxa and all groups, respectively. alternatively, a rank for each quartet can be calculated based on the total distance covered by the trajectory for each quartet ) or the total distance from each point to the average position of all points in the trajectory ). by again computing the average dtot or dav value taken over all n quartets qt that include t, and averaging over a group, we can obtain a taxon ranking and group ranking and respectively). the taxon and group ranking values are computed using all quartet trajectories. the ranking values can be used to assess the contribution to phylogenetic inhomogeneity of individual taxa and groups respectively, whereas the mean taxon and group ranking values can be individually used as statistics for detecting recombination in a sequence alignment. to identify all plausible recombinants in a sequence alignment using the quartet scanning procedure, we follow a sequential deletion procedure. when significant evidence of recombination is detected using or, the top-ranked sequence is considered as a putative recombinant and removed from the alignment. we repeat this until no significant evidence of recombination remains in the data set. because the absolute values of dmax, and hence rt and rg, are sensitive to the evolutionary models and topology evaluation criteria used, we report the ranking values as percentages of the highest ranking taxon or group. significance we evaluate three different ways to generate a null distribution for the observed test statistics and derived from the taxon-ranking procedure: monte carlo simulation, permutation, and'redistribution' of alignment columns. for the mc-simulation procedure, maximum likelihood phylogenetic trees are reconstructed using phyml, employing the general time-reversible substitution model and gamma-distributed rate variation among sites. because a null distribution needs to be generated under the hypothesis of'no recombination', we simulate replicate data sets using seq-gen down a single phylogeny reconstructed from the complete alignment. the same models and parameters were used as obtained by the maximum likelihood tree inference. since taxon ranking is based on topological incongruence, we also use a permutation procedure that homogenizes any phylogenetic incongruence along the genome. this is achieved by random shuffling of the original alignment sites to produce replicate data sets. to accommodate varying phylogenetic signal along the genome, which is often observed in real data due to different degrees of conservation within the sequence, we also use a permutation procedure that mimics the variability distribution of the alignment columns. essentially, this procedure reshuffles alignment columns according to their shannon entropy score. replicate data sets are constructed by randomly drawing for each alignment position an alignment column with the same entropy score. we use these three procedures to generate replicate data sets each, so as to compute null distributions for our test statistics. a monte carlo p-value is then estimated by counting the proportion of times the or statistic on a replicate data set is larger than the value of or observed for the original alignment. if the null hypothesis is true, this proportion is expected to be. such a p-value represents the probability of obtaining a result at least as extreme as observed for the original alignment in the absence of recombination. rdp recombinant detection the rdp software provides access to multiple recombination signal detection methods that can be used in conjunction with one another. the methods in rdp search for recombination signals within sequence triplets sampled from an input alignment. when a recombination signal is detected, the member of the sequence triplet that is the recombinant is identified using a weighted consensus of a series of fifteen heuristic tests. these tests fall into four main categories: phylogenetic profile based; branch pruning and re-grafting based; multiple recombination signal pattern analysis based; and simple distance-based metrics that, for example, compare pair-wise genetic, and phylogenetic tree distances of different subsets of sequences within the alignment on opposite sides of the detected recombination breakpoint position. these latter tests collectively identify the recombinant by determining which sequence within an alignment have relationships that change the most across a recombination breakpoint. the weighting scheme used in rdp to combine the results of all tests into a single consensus score is based on how accurately the individual tests identified hiv-1m inter-subtype recombinants from amongst a background of supposedly non-recombinant hiv-1m subtype sequences. simulation studies two different simulation studies were undertaken: to evaluate how accurately the methods detected specific recombinant sequences, and to evaluate the overall power of methods at detecting recombination signals within sequence alignments. for, sequence data sets were simulated using seq-gen according to different trees for different data partitions. data sets were simulated using both a series of symmetric and asymmetric trees. in the symmetric setup, taxon is a relatively recent recombinant with a single breakpoint, while taxa and share a recombination event deeper in the tree resulting in two breakpoints for each sequence. in the asymmetric setup, taxa and share a relatively recent recombination event with two breakpoints, while taxon represents an older recombinant with a single breakpoint. a thousand nucleotides were evolved on each of these trees using the hasegawa-kishino-yano nucleotide substitution model with gamma-distributed rate variation, resulting in a total sequence length of nucleotides. to simulate data sets with different degrees of diversity, we employed total tree depths of, and substitutions per site. in addition to rate variation among sites, we also investigated the impact of rate variation among different alignment regions by employing different relative rates for each nucleotide partition and different relative rates for nucleotide partitions. for each treatment, data sets were simulated. for, to evaluate power and false positive rates in detection of recombination signals, we examined simulated data sets previously generated for comparing the performance of recombination detection methods. these data sets are now often used to evaluate recombination detection power and false positives of new recombination detection programs. empirical data the hiv- group m alignments analyzed here were derived from a full genome alignment provided by the hiv database and manually edited using se-al. the accession numbers of the included sequences are listed in additional file two data sets were used in our analyses: a crf03_ab data set and a larger crf02_ag data set. the hiv- group m alignments had nucleotide diversities ranging from substitutions/site to substitutions/site, which roughly corresponds to θ = in the coalescent-based simulations. the siv data set, previously analyzed for recombination by bailes et al., includes protein sequences representative of eight major siv clades. their accession numbers are also listed in additional file bootstrapped phylogenetic trees were reconstructed using phyml, as described above and network analyses were performed using the neighbor-net method implemented in splitstree. to represent robust network-like relationships, the neighbor-nets only display the splits that were present in 75% of bootstrapped replicates. intrasubtype recombination was analyzed using gard method, available at. tree topology tests were performed using tree-puzzle and consel.  we present the results of three main analysis sections: simulations to evaluate the performance of different statistics and measures of significance for the quartet scanning method, and simulation studies to compare quartet scanning with other recombination detection methods in terms of both detecting recombination signals and identifying recombinant sequences; an analysis of hiv- group m inter- and intra-subtype recombination; and recombinant identification using amino acid data from more divergent piv lineages. simulated data analysis to establish the most powerful test statistic and significance assessment for the distance-based quartet scanning procedure, we analyzed data sets simulated on both symmetric and asymmetric trees, which generates three recombinant sequences per data set. to restrict the amount of computation in these simulation analyses, we focus on the first sequence ranked on top by the quartet scanning approach and perform sequential deletion whilst significant evidence for recombination can be detected. the results are summarized in figure dashed lines represent the frequency by which a recombinant is correctly ranked on top, independent of whether significant evidence for recombination can be detected in the replicate set. the full lines represent the same frequency, but now supported with significant evidence for recombination in the replicate set. the upper graphs summarize the results for symmetric trees, while the lower graphs summarize the results for asymmetric trees; the graphs on the right present the results for data sets simulated with different relative rates correlated with the alignment partitions. overall, recombination detection was more successful in symmetric trees and a correctly identified recombinant was usually a'recent' recombinant. the total trajectory distance appeared to be a poor test statistic for recombination detection; the total distance from the mean position in the quartet triangles was most powerful in symmetric trees but the maximum distance from the mean position appeared to be the best overall statistic. we thus chose dmax as test statistic for further simulation analyses. increasing sequence diversity resulted in increased power of recombinant detection in the symmetric trees. somewhat surprisingly, this was usually not the case for analysis of asymmetric trees. in general, simulation led to stronger statistical support than permutation. permutation, in turn, was more powerful than redistribution. redistribution also resulted in a markedly lower number of correctly identified recombinants, with statistical support for recombination in the replicate data set, in the presence of rate variation among alignment partitions. the fact that the partitions in between recombination breakpoints coincide with the relative rate partitions in our simulation may be at least partly responsible for this, as sites may not be efficiently permutated across breakpoints in the redistribution procedure. to assess the performance of the distance-based quartet scanning method as a statistical test for recombination signal, rather than a method for recombinant identification, we also reanalyzed nucleotide data simulated under varying levels of divergence and recombination originally presented in and compared this with other methods. in general, these analyses revealed that the performance of our test is comparable to similar methods, like recscan, but these are less powerful than the maximum chi method, which has been reported as one of the most powerful nonparametric recombination detection methods. finally, we compared the quartet scanning with triplet-based approaches and their associated heuristic tests in their ability to identify the mosaic sequences. because this requires a sequential deletion procedure in visrd, we restricted the computation by focusing on the data sets simulated using symmetric trees and a tree depth of substitution per site, similar to the hiv- group m data sets we analyze below. figure summarizes the frequencies with which taxa and were ranked on top at some stage in the sequential deletion procedure, and hence identified as recombinants, by the different approaches, as well as the total number false positives in the absence and presence of rate variation among sequence regions. the performance of the different methods is evaluated using the adjusted rand index. the adjusted rand index is the rand index corrected for chance events, which is commonly used to measure the performance of data clustering. in our case, the adjusted rand index reflects the percentage of correctly identified recombinant and non-recombinant sequences corrected for chance. therefore, higher values of this measure indicate higher classification accuracy. the methods mainly differ in the frequencies with which taxa and – the sequences descended from a common recombinant ancestor – are identified as recombinants. without rate variation among partitions, quartet scanning based on simulation achieves the best rate of detection of these recombinants and obtains the highest performance score. using permutation and redistribution, top-ranked recombinants were still relatively frequently associated with significant evidence for recombination in the alignment, but the lower detection rate of taxon indicates some sensitivity decrease. the presence of rate variation among partitions increases the false positive rate in many methods; quartet scanning with simulation is particularly sensitive to this. a permutation procedure, however, still achieves good detection rates without an excessive number of false positives. as also suggested by the first set of simulation analyses, redistribution becomes less powerful in detecting recombination signal with the simulated rate variation scheme. to investigate whether this results from inefficient permutation across breakpoints correlated with rate variation, we also performed additional simulations using a rate variation scheme that does not coincide with recombination breakpoints. also in this case, redistribution did not prove to be a very powerful approach to assess significance. we therefore restrict significance assessment to simulation and permutation in subsequent analyses of real sequence data. hiv- intersubtype recombination we applied the visrd based recombination detection methods to the hiv- group m full coding genome data set which included representative sequences for subtypes a-d, and f and the circulating recombinant form crf03_ab. the crf03_ab sequences have been clearly identified as having descended from a common ab recombinant form that circulated amongst injecting drug users in russian and ukrainian cities. in the quartet scanning procedure, the sequences were analyzed by grouping them according to their subtype or crf assignment. this group model only considers quartets composed of sequences from four different pre-specified groups, which allows the analysis to focus exclusively on inter-subtype recombination events while reducing the computational burden. in figure a highway plot provides a visual summary of the most variable trajectories obtained during the quartet scanning procedure. the three lanes in the plot represent the different topology-supporting regions in a quartet triangle, and the change in support is calculated as a function of the alignment position for each quartet. the highway plot clearly indicates two breakpoints at approximately nucleotide positions and. table lists the results of the taxon-ranking procedure for the different groups in this data set, indicating that quartets containing a crf03_ab sequence show the highest variability in phylogenetic clustering along the genome based on both dmax and dav. evidence for recombination in the data set is statistically supported by both simulation and permutation. the parental subtype a and b sequences are ranked just below the recombinant group due to the high trajectory variability of quartets including one or both parental subtypes and the recombinants. even after deletion of the crf03_ab sequences, there was still significant evidence for recombination; a more in depth investigation of hiv- group m recombination is provided below. the group ranking values are listed as percentages of the highest ranked group; these values are shown for the maximum distance and total distance from the mean position in the quartet-mapping triangle. the number of sequences per group is indicated between brackets the crf03_ab example indicates that the values allow us to test for the presence of recombination while taxon ranking itself indicates the most likely recombinants in the data set without assigning query and parental sequences a priori. together, this represents an ideal approach to evaluate recombination patterns in cases where prior identification of non-recombinant reference strains is problematic, like crf02_ag. to investigate the crf02_ag issue and other putative recombination events in hiv- group m, we compiled a full genome data set that is more representative of hiv- group m diversity. also for the analysis of this data set, we employed a grouping model for the different subtypes and crfs. although the value for the data set including crf02_ag indicates the presence of recombination, the taxon ranking does not suggest crf02_ag sequences as most plausible mosaic genomes. in agreement with the recent analysis of abecasis et al., the'pure' subtype g sequences are proposed as the most likely recombinants while crf02_ag and subtype j are suggested as parental lineages because they were ranked right below the putative recombinant. when the subtype g sequences were removed, recombination in the remaining data could not be excluded based on the permutation statistics. this time, subtype k was at the top of the taxon ranking. it is interesting to note that a recombinant origin for subtype k has been proposed, but its designation was not changed since it was already an established subtype in the literature. after excluding both subtype g and k, there was still significant evidence for recombination in the remaining sequences with subtype a heading the taxon ranking. to our knowledge, subtype a has not been proposed as a intersubtype recombinant before and the recombination signal we detected might result from intra-subtype recombination or intersubtype recombination for particular lineages within subtype a. we investigated this in more detail by removing the grouping model for the different a sub-subtypes and by investigating the evolutionary relationships using neighbor-nets. when no grouping for subtype a was used, sub-subtypes a a and a were ranked on top but not sub-subtype a. the neighbor-net indeed demonstrates network-like behaviour among these lineages within subtype a. having removed sub-subtype a a and a our recombination test was still significant and subtype j was now at the top of the ranking. although recombination has not been clearly suggested for this subtype, it does indeed cluster differently in each major gene tree of the hiv genome. further removing of subtype j, there was still significant recombination signal, this time with subtype d heading the ranking. without the grouping model, three subtype d sequences appeared to be the descendents of a recombinant lineage, which was also confirmed by the neighbor-net analysis. the final recombination signal was detected with subtype f heading the ranking, despite a stable clustering of this subtype in the major hiv genes encoding structural proteins. recombination remained detectable after deletion of subtype f using the mc-simulation approach, but given the more likely confounding effect of rate variation for mc-simulation and the increasing difficulty of correctly identifying'deeper' recombinants, we did not explore further taxa deletion. a different order in putative recombinant deletion was suggested after identification of recombination in subtype g, k and subsubtypes a using dav as a test statistic, and subtype j was not suggested as a recombinant by this procedure. the intrasubtype recombination detected in subtype a and d in both procedures was confirmed by gard analysis. we also employed our method to investigate a highly divergent hiv- variant from the democratic republic of congo. although the drc strain was positioned differently in phylogenetic trees of different genome regions, it did not cluster with any of the know subtypes and hence it was suggested as a candidate for a new subtype. as reference strains, we included the sequence set for which no significant recombination could be detected. our taxon-ranking procedure indicated that quartets containing this unclassified variant show the highest phylogenetic variability along the genome, and this was associated with significant evidence for recombination in the alignment. to confirm alternative clustering patterns for this lineage and to demonstrate significant incongruence despite convincing bootstrap support values, we reconstructed phylogenies for the major genes and performed tree topology tests. the trees indicate different clustering for the drc lineage in the three major genes, and despite low bootstrap support for the alternative clustering, tree topology tests consistently indicated significant incongruence. finally, we also analyzed crf01_ae sequences, for which only the subtype a lineage is available as closely related parental descendent. the taxon-ranking procedure still revealed significant recombinant signal and implicated this crf as the lineage responsible for the incongruence. recombination in primate immunodeficiency virus evolution recombination in more divergent simian immunodeficiency viruses was investigated for the amino acid data set of bailes et al., which includes representatives of eight major siv lineages for which full-length sequences were available at the that time. because there are no equidistant siv clades similar to hiv- group m subtypes, we did not apply the grouping model in this case. however, when sister taxa representing one of the eight major lineages were ranked as the top two sequences, we removed both taxa in the subsequent taxon-ranking test. the protein sequences were analyzed using a corrected distance measure based on an amino acid substitution matrix inferred from a larger piv data set, and significance was assessed using permutation. additional file summarizes the results using the dmax statistic, while the results using dav are provided in additional file the ranking procedure provided evidence for recombination in the siv sequences. moreover, in agreement with the findings of a hybrid origin of siv in chimpanzees, sivcpz was suggested as contributing the most to this recombination signal. however, recombination was also detected during subsequent deletion of several siv lineages, leaving only a single quartet that did not deviate significantly from clonality. there were some differences in the order by which sivs were identified as recombinants using dav, but the final result was the same. discussion in this study, we have developed an extension of the visual recombination detection method to rank taxa or groups of taxa according to their contribution to phylogenetic inhomogeneity, and to test whether the overall phylogenetic incongruence has been significantly shaped by recombination. during a sequence alignment scan, phylogenetic inhomogeneity is quantified using a quartet mapping approach. taxa or groups of taxa are subsequently scored using the average inhomogeneity of all the quartets they belong to. to assess whether significantly more inhomogeneity is observed than expected under a recombination-free scenario, we tested mc-simulation, permutation, and a redistribution approach. simulation analyses suggested that identifying recombinant sequences without prior specifications is a challenging problem. only relatively recent recombinants were readily identified in symmetric trees in which recombination was rare. in asymmetric trees, correct recombination identification proved to be more difficult. exploratory analyses using rdp indicated that this was also the case for the triplet-based methods. however, the recent recombination event in the asymmetric trees occurred between somewhat more closely related parental strains, with a different mosaic structure, and involved two descendants of the same recombinant. hence, symmetric and asymmetric tree simulations cannot be easily compared. asymmetric topologies are generally more problematic in phylogenetic inference; more sequence data is required to recover such topologies and they are much more susceptible to substitution saturation than symmetrical ones. the latter may provide an explanation as to why increased sequence divergence did not result in increased power for recombination detection, as is the case for symmetric trees. we found marked performance differences between different measures of quartet variability. the poor performance of the total path distance may be explained by the fact that the sliding window procedure follows the path between two regions supporting alternative topologies in the quartet-mapping triangle in a stepwise manner, with each step represented by different alignment windows. steps covering similar distances could also be taken within a region supporting a particular topology, and thus lead to similar total distances in the absence of real topology changes. no matter how many steps need to be taken to the point most distant from the average position, the dmax measure only considers how remote this point is from the average positions. the power of the maximum distance or the average distance from the mean position will depend on the mosaic pattern of the recombinant. dmax may be more sensitive to a small region with a different evolutionary history, while dav may be more sensitive in case of multiple breakpoints. the simulation setup was necessarily simplistic and we have not explored the variety of recombination scenarios that occur in reality. when recombination is frequent, for example, all the methods would be expected to suffer from an elevated failure rate. also the sampling rates of recombinants relative to the parental descendent sequences will have a profound impact on the ability to identify mosaic sequences without prior reference specifications. when few representatives of the recombinant lineage are sampled and relatively closely related descendents of both parents are present, as is the case in our phylogenetic simulations, than the quartet variability should be attributed relatively easily to the recombinant sequences by averaging the dmax measure over the quartets that contain these sequences. when, however, many sequences are included that share the same recombination history, more quartets will be composed of multiple recombinant sequences of the same type. these quartets may not generate strong phylogenetic inhomogeneity, and in the case of all four sequences representing the same recombinant, there should be no recombination signal detectable. although sampling rate inequalities for pre-defined lineages may be controlled to some extent by the grouping constraint we proposed, it will no doubt adversely affect taxa ranking. this may also at least partly explain the poor performance for the asymmetric tree simulations. here, the more recent recombinant event had two descendents, and frequently, a parental descendent sequence was top ranked instead of one of the recombinants. this again highlights that the current simulations only offer a starting point for comparison of different methods and their relative ability to identify recombinant sequences, and more complexity needs to be explored in the future. for the quartet scanning method, generating null distributions using mc-simulation was the most powerful procedure for identifying evidence of recombination, but it also suffered from a higher false positive rate in the presence of rate variation among genome regions. the simulated patterns of rate variation are within the range of rate variation that can be expected in hiv genomes. in this case, permutation still achieved good detection rates with far fewer false positives. this difference illustrates that the mc-simulation technique performs very well when the simulated evolution meets the assumptions in the parametric bootstrapping procedure. under model misspecification, however, it may still be powerful but could also be particularly misleading. although rate variation can be accommodated in the simulation procedure using a discrete gamma distribution, it assumes that rate variation is independent and identically distributed along the genome. this assumption will rarely be met for real data situations and it can have an important impact on comparative analyses of data partitions across the genome. the redistribution approach, which was designed to homogenize phylogenetic incongruence, but still mimic the pattern of rate variation present in the real data sets, was the least sensitive. even when rate variation does not correlate with breakpoints, this approach does not easily lead to rejection of clonality. this indicates that permuting sites only with the same entropy score might be too stringent to efficiently homogenize recombination signal. in this respect, it is important to realize that a sequence alignment can include sites with many different entropy values, and the number of sites per entropy value will naturally show a large variation. for a typical simulated data set comprising taxa, site patterns can break down in about different site entropy values with, for example, a few patterns having unique entropy values and about 13% of patterns being single nucleotide polymorphisms that are represented by a single entropy value. therefore, more efficient shuffling under more relaxed rate variation constraints may be achieved by binning site patterns based on their entropy value. taxon ranking based on quartet scanning is geared towards the detection of recombinant sequences. as a test for recombination signal, more powerful methods have been devised, and we do not provide a formal test for breakpoint locations. also for the latter, we refer to other developments. however, automated detection of recombinant lineages on its own is already a difficult problem to tackle. the quartet scanning approach performed well relative to various triplet-based methods in identifying recombinants simulated using a phylogenetic approach. however, more comprehensive simulations are needed for a better assessment of the powers and pitfalls in different situations. the essential difference between the quartet scanning approach and the various rdp methods does not lie in the use of quartets versus triplets, but how their information is further employed to identify the recombinants. while the quartet approach deduces a ranking for taxa or groups of taxa from the phylogenetic variability in all these quartets, additional tests are used by rdp to distinguish recombinants from parental strains in triplets that resulted in significant evidence for recombination. it is interesting to note that the weighting scheme used by rdp to combine the results of these tests into a single consensus score is currently based on how accurately the individual tests identified hiv-1m inter-subtype recombinants. because of the crf02_ag-subtype g classification artefact, a new weighting optimization might be required. moreover, quartet scanning might make a useful addition to the tests employed by rdp the quartet approach presented here has particular practical limitations. the number of quartets increases as a fourth order polynomial in the number of taxa, and although the group constraint can alleviate some of the computational burden, both memory requirements and speed of computation will become impractical for very large data sets. this will be a particular burden when many sequential deletions of recombinants are required and significance needs to be assessed each time using replicate data sets. although these procedures can be automated and quartet information for a permutation alignment can be stored and reused, it would be very useful to be able to test every taxon or taxon group without having to rescan pruned data sets. our experience with the current implementation leads us to conclude that taxa is roughly the upper limit for this method in practice. there are obviously less triplets than quartets for any data set, and rdp does not evaluate evidence for recombination in triplets by analyzing replicate data sets, which makes them faster in computational practice. in sliding window analyses, it is important to determine how much phylogenetic discordance can be expected from the stochastic nature of the substitution process in the absence of recombination. in bootscan analysis, a query sequence is generally inferred as a recombinant when a high bootstrap support is observed for clustering with different reference sequences in different gene regions. however, this criterion will not be adequate to assess variable clustering deeper in the phylogenetic trees. this is illustrated by the highly divergent hiv- variant from the drc. although phylogenetic inhomogeneity was originally noted for this strain, the absence of strong bootstrap support for the different clustering deeper in the tree led to the suggestion that it could be a new subtype. however, our test indicates that there is more phylogenetic inhomogeneity than expected by chance and thus classifies the drc strain as a probable recombinant lineage. significant phylogenetic incongruence between major genes was also confirmed by various tree topology tests. both subtype g in hiv- group m and sivcpz in piv are confirmatory examples, but the application of quartet scanning demonstrates how recombinants can be easily identified for cases that required in depth analysis to distinguish parental from recombinant lineages. the pervasive recombination in piv lineages has also been suggested by exploratory bootscanning and network analyses. in addition to piv recombination and hiv- intersubtype recombination, our test was also sensitive to intrasubtype recombination. both in subtype a and d, we identified recombination signals that were supported by network and gard analyses, and in both cases, the recombination signal did not stem from a single mosaic sequence but from a recombination event that was shared among several sequences. although we have identified recombinant genomes at different evolutionary scales in immunodeficiency virus evolution, it may prove a challenging endeavor to clearly delineate breakpoints for many of these recombinants. similar to other sliding window analyses, breakpoints will be readily identified if viruses from two well-supported phylogenetic clusters recombined relatively recently. however, our method may also pick up recombinants in trees without well-supported groups, ancestral recombinants that were generated before well-supported clades were established, or ancestral recombinants of lineages that have now diverged considerably.  in conclusion, we have developed a recombination detection method that does not require a priori assignment of reference and query sequences, statistically evaluates the recombination signal in nucleotide and amino acid alignments, and identifies the most plausible recombinants. taxon ranking and the associated statistical test complement the visual recombination detection method. collectively, these developments comprise a powerful and versatile methodology that will make a useful addition to the array of recombination detection methods. abbreviations crf: circulating recombinant form; piv: primate immunodeficiency virus; siv: simian immunodeficiency virus; visrd: visual recombination detectionl; rdp: recombination detection program; me: minimum evolution; f84: felsenstein' model of evolution; gtr: general time-reversible substitution model competing interests the authors declare that they have no competing interests. authors' contributions pl conceived of the study, designed the analysis strategy, analyzed the data and drafted the paper. ml carried out the visrd software implementation and supported data analyses. dm wrote a modified rdp version to analyze simulated data and supported data analysis using triplet-based methods and its interpretation. vm conceived of the study, participated in its design and coordination, and helped to draft the manuscript. all authors read and approved the manuscript. supplementary material additional file 1 table a1: hiv and siv sequence data. a table listing the accession numbers and subtype/crf assignment or piv lineage for the sequences analyzed in this study. click here for file additional file 2 coalescent simulations. this file describes the coalescent simulations and presents the results for different visrd settings and for a comparison of visrd to other approaches. click here for file additional file 3 additional table a this table lists the group ranking for quartet scanning of the crf02_ag data set using dmax. click here for file additional file 4 additional table a this table lists the group ranking for quartet scanning of the crf02_ag data set using dav. click here for file additional file 5 additional table a this table lists the group ranking based on quartet scanning of hiv- group m including a divergent drc lineage. click here for file additional file 6 tree reconstructions and topology tests for gag, pol and env data sets including the drc lineage. this file includes a figure presenting the results of a phylogenetic analysis for gag, pol and env data sets including the drc lineage and a table listing tree topology test results for the same data. click here for file additional file 7 additional table a this table lists the group ranking based on quartet scanning of hiv- group m including crf01_ae. click here for file additional file 8 additional table a this table lists the taxon ranking for quartet scanning of the primate immunodeficiency viruses using dmax. click here for file additional file 9 additional table a this table lists the taxon ranking for quartet scanning of the primate immunodeficiency viruses using dav. click here for file acknowledgements we thank liz bailes for providing the siv sequence data and david posada for providing the simulated data sets and previous results for other recombination detection methods. we thank joe parker for programming assistance. pl was supported by a marie-curie and fwo postdoctoral fellowship. pl, ml and vm thank the isaac newton institute for mathematical sciences, cambridge, uk, for hosting them in the context of its phylogenetics programme where part of the research presented in this paper was carried out.
recombination has a profound impact on the evolution of viruses, but characterizing recombination patterns in molecular sequences remains a challenging endeavor. we also demonstrate hiv- intrasubtype recombination, confirm the hybrid origin of siv in chimpanzees and further disentangle the recombinant history of siv lineages in a primate immunodeficiency virus data set. we compare the quartet approach with triplet-based methods that employ additional heuristic tests to identify parental and recombinant sequences. in agreement with recent findings, we provide evidence that the presumed circulating recombinant crf02_ag is a'pure' lineage, whereas the presumed parental lineage subtype g has a recombinant origin. here, we extend a quartet-mapping based recombination detection method to enable identification of recombinant sequences without prior specifications of either query and reference sequences. quartet scanning makes a valuable addition to triplet-based methods for identifying recombinant sequences without prior specifications of either query and reference sequences. analysis of phylogenetic simulations reveal that identifying the descendents of relatively old recombination events is a challenging task for all methods available, and that quartet scanning performs relatively well compared to the triplet based methods.
 recombination has a profound impact on the evolution of viruses, but characterizing recombination patterns in molecular sequences remains a challenging endeavor. despite its importance in molecular evolutionary studies, identifying the sequences that exhibit such patterns has received comparatively less attention in the recombination detection framework. here, we extend a quartet-mapping based recombination detection method to enable identification of recombinant sequences without prior specifications of either query and reference sequences. through simulations we evaluate different recombinant identification statistics and significance tests. we compare the quartet approach with triplet-based methods that employ additional heuristic tests to identify parental and recombinant sequences.  analysis of phylogenetic simulations reveal that identifying the descendents of relatively old recombination events is a challenging task for all methods available, and that quartet scanning performs relatively well compared to the triplet based methods. the use of quartet scanning is further demonstrated by analyzing both well-established and putative hiv- recombinant strains. in agreement with recent findings, we provide evidence that the presumed circulating recombinant crf02_ag is a'pure' lineage, whereas the presumed parental lineage subtype g has a recombinant origin. we also demonstrate hiv- intrasubtype recombination, confirm the hybrid origin of siv in chimpanzees and further disentangle the recombinant history of siv lineages in a primate immunodeficiency virus data set.  quartet scanning makes a valuable addition to triplet-based methods for identifying recombinant sequences without prior specifications of either query and reference sequences. the new method is available in the visrd v. package. 