Detection and characterization of horizontal transfers in prokaryotes using genomic signature

Sequence signature

The sequence signature is defined as the frequencies of the whole set of short oligonucleotides observed in a sequence (26,31). It can be easily obtained thanks to a very fast algorithm derived from the Chaos Game Representation (CGR) (33), which allows coping with a 1 Mb sequence in a few seconds on a laptop computer. Signatures may be visualized as square images where the color (or gray level) of each pixel represents the frequency of a given oligonucleotide (called word thereafter) (31) (for examples of signatures, see Supplementary Materials 2, 4 and 6).

DNA sequences

DNA sequences are gathered from GenBank. The genomes of 22 prokaryotes are scanned for HTs, B.subtilis, E.coli and H.influenzae genomes being given a special attention to illustrate our approach. In particular, B.subtilis and E.coli provide valuable benchmark thanks to the set of previous works addressing that very issue (12,14,16,34–37). Signatures of about 12 000 species are obtained from genomic sequences longer than 1.5 kb. Sequences derived from the same species are concatenated for accuracy purposes. Species from the three domains of life, archaea (∼260 species), bacteria (∼3950 species) and eukarya (∼6750 species) as well as viruses (∼1300 species), are represented for a total amount of 1.0 Gb.

DNA sampling

The detection of atypical regions is based on the observation of deviation of local signatures (i.e. signature of small fragments of DNA) from the genomic signature of the recipient species. Genomes are consequently sampled by means of a sliding window with an appropriate size. In fact, it would be interesting to have windows the smallest as possible for highest sampling accuracy. However, intra-genomic variability of signature increases for small windows. In addition, variability depends on species and word length. Base composition (1-letter word), 2- and 3-letter words are poorly species-specific: they do not allow a good discrimination between species (25,27). As a general rule, the longer the words (up to 9-letter long), the higher the specificity of the signature (25,27,31). However, counts of long words in small windows are too low to allow a reliable estimation of the parameters. In our hands, the analysis of 4-letter words in a sliding window of 5 kb (with a 0.5 kb step) offers a good trade-off between reliability of count, file size and computational charge, whatever the species. In addition, a double-strand signature (called local signature thereafter) is computed for each window to get rid of variations induced by strand asymmetry (38–42).

For illustration purposes, local signatures are developed as vertical vectors and stacked together in genome order to give an overall picture of word usage variations along each genome. In such plots, horizontal lines show the variation in frequency of words along the genome, whereas local changes in word usage appear as vertical breaks (Figure 1).

Recipient species signature

Considering that the greatest part of the genome is species-typical, the signature of the recipient species might have been estimated from the analysis of the whole sequence. Although the vast majority of local signatures look mostly the same (believed to be instances of the recipient species signature), some of them may greatly differ. In order to avoid potential biases linked to these outliers, it has been subsequently decided to select typical local signatures on the basis of their similarities, observed after clustering. The underlying idea is that typical local signatures aggregate in few large groups, whereas outliers are found in small complementary groups at a great distance from the recipient genome signature. Groups were consequently determined with the K-means clustering tool, using every scheme of clusters between 3 and 8 for each species. Finally, the best scheme of clusters was obtained by a decision tree-based partition [CART algorithm (43)]. The purpose of the CART algorithm is to predict values of a categorical dependent variable (clusters of local signatures in this work, each signature being characterized by its distance to the estimated genomic signature) from one or more continuous and/or categorical predictor variables [the different clustering schemes (3–8 clusters) in this work]. The CART algorithm thus provides an optimal split between groups collecting signatures close to the estimated recipient genome signature and the others groups. For each species, a clustering scheme is selected (e.g. the 5-group clustering) and a partition offered (continued example: group 2 and 3 on one side; 1, 4 and 5 on the other). The recipient species signature is subsequently calculated as the mean of the signatures of the groups belonging to the partition with the smallest distance to the estimated genomic signature.

Atypical regions

Comparison of signatures is made possible, thanks to an Euclidian metric, accounting for differences in word usage. It must be pointed out that distances between signatures are calculated for high dimensional data (256 dimensions corresponding to the 256 different 4-letter words) and are consequently subjected to the so-called ‘concentration of measure phenomenon’ (44). All distances in a high dimension space seem to be comparable since they increase with the square root of the dimension of the space, whereas the variance of their distribution remains unchanged. In fact, the radius of the hyper sphere holding 99% of the signatures of our database is only seven times the nearest neighbor distance (smallest distance between two species). Small differences in distance may consequently be considered highly significant.

For each species, a set of recipient-specific distances is obtained, every local signature belonging to the large clusters being given a distance to the host signature. In order to select outlying signatures, a cut-off distance is chosen on the basis of the distribution of distances observed for each species. It appears that the 99% percentile offered a good trade-off between sensibility and specificity for outlier detection (for impact of the threshold on detection of atypical regions, see Results). Most signatures from minority clusters are detected in this way. Isolated signatures are detected as well, while very few signatures from the recipient species clusters are selected (1%). Outliers together with the flanking regions on the genome are later on reanalyzed with smaller window and step (1/10^th of the original size typically) in order to more accurately determine their limits, when signal-to-noise ratio allows it.

Finally, the gene content of all detected regions is analyzed with the help of species dedicated databases [Genome Information Broker, http://gib.genes.nig.ac.jp/]. A BlastN search (GenBank, default settings) is carried out for each atypical region in order to identify the origin of potential HTs if homology is high enough.

Search for the origin of atypical regions

About 12 000 species (including chromosomal, plasmidic, mitochondrial and chloroplastic DNA) from GenBank are found eligible for a genomic signature. Given the signature of an atypical DNA fragment, species with a close signature might be considered as potential donors. Such a screening is performed for every atypical region of the 22 species under consideration. The first five nearby species are retained when their distance to the outlier was donor-compatible.

B.subtilis genome analysis

A clustering with three classes allows assessing the signature of B.subtilis. The most populated class (collecting 84% of the segments) is chosen to represent B.subtilis. For this sub-population, the mean distance (arbitrary unit) to the recipient (centroid of the class) and the cut-off distance are 126 and 204, respectively (Table 1). Runs of contiguous outlying windows sharing the same cluster are considered as single transfer events. As a consequence, 58 regions (Figure 2a and Supplementary Material 2) fall beyond the cut-off distance and are thus potential candidates for hosting foreign DNA (for a segmentation of the B.Subtilis genome in terms of genes, see Supplementary Material 3). Figure 2b illustrates the accuracy of segmentation of an atypical region obtained by using a sliding window of 0.5 kb with a 50 bp step.

rRNA genes make up ∼1.1% of B.subtilis genome (Table 1). All rRNA genes are found in the outlier population. In addition, all windows containing rRNA are assigned to a specific cluster. In fact, it is known that rRNA has its own signature, which is at variance from the host signature (12). rRNA genes account for 7% of the outliers (tRNAs are not considered in this study, because their size is too small to generate a significant deviation from the host signature if they are isolated).

A total of 86% of the B.subtilis genome should be considered as B.subtilis typical (Table 1). When looking for the origin of B.subtilis segments in the 12 000 signature database, B.subtilis appears in the 10 first potential donors for 84% of the whole set of 5 kb sequences that can be derived from its genome. This result confirms that segments having signatures belonging to the predominant clusters are good representatives of the recipient species signature.

The 49 rRNA-free atypical regions vary in size from 1.5 to 135 kb and make up 13% of the total genome (Table 1). About 50% of atypical regions are less than (or around) 6 kb long. Distances of outlier from first potential donor often fall within the intra-genomic range (Table 1 and Supplementary Material 2): 75% of B.subtilis atypical regions have first donors lying at a distance <234 (mean cut-off distance for the 22 species analyzed in this paper). Potential donors are clearly distinct from B.subtilis (Supplementary Material 2). Although outliers are distributed all over the genome, several types can be distinguished on the basis of their signatures; potential donors seem to come from a few sets of species with similar signatures (Table 1 and Supplementary Material 2). The most important cluster includes bacteriophage SPBc2 and its neighbors, another concerns the Enterococcus genus.

However, in some instances, the outlier-to-donor distance is too great to consider the ‘closest’ species as potential donor. In contrast, unusual small values deserve a specific attention. In particular, the very small distance between bacteriophage SPBc2 and ‘2150751–2285750’ atypical region (d = 2) allows to spot the part of B.subtilis genome where bacteriophage SPBc2 is incorporated (12,47). Other regions in the genome are also found similar (in terms of signature) to bacteriophage SPBc2. Most of them correspond to bacteriophages, imbedded in B.subtilis genome, whose free forms are not sequenced (12,47). Observed similarities with SPBC2 are, however, expected since signatures of phages usually share some characteristics with the species they infect (48). The SPBc2 sequence is the only foreign sequence identified in B.subtilis, using homology as criterion (BlastN, with parameters set to default). In fact, Blast analysis of B.subtilis outliers leads to contrasted results. Besides SPBc2 and 7 out of 9 prophages imbedded in the genome, the only atypical regions identified are those containing the 30 rRNA genes coded in B.subtilis genome. The only few genes that are homologous to parts of atypical regions are found in species belonging to the Bacillus genus. It is interesting to note that no house-keeping genes (except rRNA) are detected in atypical regions. In fact, a great number of genes in atypical regions (except bacteriophage genes and rRNA) have no known function.

H.influenzae genome analysis

A clustering with five classes is required to determine the recipient species signature of H.influenzae. The three most populated classes (collecting 94% of the segments) are chosen to calculate the H.influenzae signature. Mean distance to host and cut-off distance is subsequently found equal to 130 and 239, respectively (Table 1). Similarly to B.subtilis, one cluster (1.5% of H.influenzae genome) is devoted to the 18 rRNA gene copies (Table 1). A total of 91% of rRNA is labeled atypical and account for 29% of the outliers.

Analysis of Table 1 shows that 95% of the H.influenzae genome should be considered as H.influenzae typical. In fact, H.influenzae is one of the 10 first potential donors for 92% of all 5 kb sequences that can be derived from its genome. As already observed for B.subtilis, the concordance of these two percentages corroborates the partition procedure used for the selection of typical/atypical fragments.

The 13 rRNA-free atypical regions vary in size from 1.5 to 19.5 kb and make up 3.3% of the genome (Table 1, Annex 4 and Figure 3, see Annex 5 for a segmentation of the H.influenzae genome in terms of genes). About 50% of atypical regions are less than (or around) 2.5 kb long. Numbers for H.influenzae are clearly at variance with those for B.subtilis: a smaller percentage of the genome qualifies as atypical and the average size of atypical regions is also smaller. This result is examined below in the context of intra-species signature variability (see Discussion).

E.coli genome analysis

A clustering with six classes is required to determine the recipient species signature of E.coli. The main features are summarized in Table 1. The potential donors of the 84 RNA-free atypical regions are given in Annex 6 (for a segmentation of the E.coli genome in terms of genes, see Annex 7). It is worth noticing that 56% of E.coli potential donors belong to the Enterobacteriales family. Segmentation in terms of genes is displayed in Annex 7. The analysis of this genome is particularly useful for the comparison with literature (see below).

Comparison with other methods

Numerous approaches for detecting horizontal gene transfers have been proposed in the last 2 decades. Phylogenetic trees of protein or DNA sequences, unusual distribution of genes, nucleotide composition (including codon biases) are some of the HT features that are considered within the framework of these models (16,34), Hidden Markov Models (HMMs) (12,14,35) and Factorial Correspondence Analysis (FCA) (37) are some criteria that are currently employed. Each of the resulting models has its own advantages and caveats (28–30). As it has been recently pointed out by Ragan (49) and Lawrence and Ochman (50), each approach deals with a particular subset of HTs, being for example more efficient for detecting recent transfers, or more effective for the detection of ancient HTs. Our approach, which is clearly based on oligonucleotide composition, assumes that different species have different signatures but does not rely on any other assumption. It is not surprising, therefore, that the genomic signature approach provides results (in terms of % of DNA transferred) in reasonable agreement with those proposed by Garcia-Vallve (16) and Nakamura et al. (14) for the 22 species that were analyzed in common. Correlations between percentages of HTs found by these three methods are highly significant (R = 0.42, P = 10⁻³ and R = 0.62 and P < 10⁻⁴ between this work and Nakamura et al. and Garcia-Vallve et al., respectively; R = 0.40 and P = 2 10⁻³ between Nakamura et al. and Garcia-Vallve et al.). However, agreement in percentage of atypical DNA does not imply that the same genome regions are detected by different methods.

Two species are extensively studied for HT content: B.subtilis (five methods including ours) and E.coli (six methods including ours). H.influenzae is also analyzed by Garcia-Vallve (16) and Nakamura (14). Comparisons of methods are presented in Tables 2–4 and detailed in Supplementary Materials 3, 5 and 7. A voting procedure (majority rule) has been implemented to determine the status of genes with respect to atypicality. For that task, our initial analysis is converted in terms of genes (Supplementary Materials 3, 5 and 7). Degree of agreement between methods is subsequently observed using the statistical Kappa coefficient (51). Kappa measures the degree of agreement on a scale from minus infinity to 1. A Kappa of one indicates full agreement, a Kappa of zero indicates that there is no more agreement than expected by chance and negative values are observed if agreement is weaker than expected by chance (a very rare situation).

B.subtilis. Garcia-Vallve et al. (16), Nicolas et al. (12), Nakamura et al. (14), Moszer et al. (36) and we are the voters concerned with the analysis of B.subtilis genome (Table 2, B.subtilis). Proportions of horizontally transferred genes are quite similar (14, 13, 11, 13 and 15%, respectively). The number of detected genes per method is close, ranging from 457 for Nakamura (14) to 599 for this work (median 537). Detailed votes are given in Table 2. Among the 4100 genes of B.subtilis genome, 1011 genes are detected by at least one method (about 25% of B.subtilis genes). The number of ‘single vote’ genes ranges from 116 for Garcia-Vallve (16) to 47 for Nicolas (12). A total of 470 genes make up the majority consensus set and we detected 453 of them, which is the best score of the five methods. The best agreement with the majority consensus (in terms of Kappas) is reached by Nicolas (12), followed by our method and Moszer (36) (Table 2). Our method gets the best agreement with Nicolas (12) and the worst with the other HMM method used by Nakamura (14) (pairwise Kappa comparison, Table 2 and Supplementary Material 3). In fact, Nakamura approach is at variance with every other approach (14). It gets the lowest Kappa with the majority consensus or with whatever other methods. From Table 2, the probable number of HT genes in B.subtilis would range from 230 to 1011 with a ‘reasonable’ estimation around 470 corresponding to the majority consensus. It is to be noted that our method is unable to find two genes that are detected by every other methods (Supplementary Material 3). These genes are 338 and 236 nt long, respectively, as compared with 2500 nt, the median size of atypical regions detected by our method (Table 1). Clearly, our method is not appropriate for detecting short isolated atypical genes.

H.influenzae. Garcia-Vallve (16), Nakamura et al. (14) and we are the voters concerned with the analysis of the H.influenzae genome (Supplementary Material 5 and Table 3, H.influenzae). The originality of results obtained by Nakamura (14) is the salient feature of this comparison. The number of detected HT genes is more than twice higher for Nakamura et al., whereas the part belonging to the majority consensus is the smallest (Table 3). Eleven genes are detected both by Garcia-Vallve and Nakamura (14,16) but not by our method; however, the small number of voters precludes any specific comment in this respect. The probable number of HT genes in H.influenzae would range between 11 and 273, with a ‘reasonable’ estimation around 60 (majority consensus of 57) (Table 3).

E.coli. Garcia-Vallve et al. (16), Hayes and Borodovsky (35), Lawrence and Ochman (34), Nakamura et al. (14), Medigue et al. (37) and we are the voters considered for the comparison. Proportions of horizontally transferred genes are 8, 15, 28, 17, 9 and 12%, respectively. Among the 4288 genes of E.coli, 1732 are detected at least once (40% of E.coli genes), but full agreement is only observed for 64 genes (1.5%) (Table 4, E.coli). Obviously, full agreement is expected to be small if voters are numerous. The majority consensus amounts to 342 genes. This number may provide a decent estimation of HT genes in E.coli. The approach of Lawrence (34) gets the greatest number of genes in common with the majority consensus (335 out of 342) at the price of a low specificity (Kappa = 0.43, Table 4), whereas Hayes' approach (35) gets the smallest number (186 out of 342). Our method detects 261 consensual genes, a score very close to the median score of the six methods (252). Four genes (Supplementary Material 7) are detected by all methods except ours. They are isolated and their lengths are ∼750 nt. As already pointed out, our method is not suitable for the detection of short isolated HTs.

The three methods in agreement with the consensus are the FCA of Medigue (37), the base composition and CAI approach of Garcia-Vallve (16) and our method (Kappa values: 0.81, 0.74 and 0.67, respectively). Pair wise comparison (Table 4), as well as correspondence analysis (Table 4, inset), shows that the HMM approach of Hayes (35) and the base composition and CAI method of Lawrence and Ochman (34) provide original results. They exhibit low Kappa values when checked against the majority consensus as well as other methods (Table 4). It is surprising that Hayes and Nakamura methods (14,35) get discordant results though they are based on similar approaches. The situation is similar for Garcia-Vallve and Lawrence methods (16,34) (Table 4).

The results obtained by Hayes and Borodovsky (35) are clearly at variance with the others (Table 4). Although the proportion of claimed outliers is within the range of published numbers for E.coli (14,16,24,34,35,37), 37% of them are method-specific, and the agreement with other methods is weak (Table 4). Hayes and Borodovsky have obviously developed an approach based on HMM dealing with specific outliers. Lawrence and Ochman (34) also get a poor rating especially because they detect about twice as many genes as the other authors do (Table 4).

It is worth noting that if the cut-off distance for our method is lowered, i.e. 95% instead of 99% for instance, some of the ‘single vote’ genes are dug out (for details about the impact of the cut-off distance, see Supplementary Material 7). Meanwhile, the percentage of outliers as reported by our approach rises to 20% and the percentage of ‘single vote’ genes reaches 24%. As expected, a high cut-off distance provides few single vote genes at the risk of missing some potentially transferred genes. Lowering the cut-off increases the proportion of single vote genes with the advantage of detecting most of the potential transfers (Supplementary Material 7). There is obviously a continuous grading in gene ‘atypicality’. It is suggested to first consider most ‘consensual’ genes as potential HTs and then apply amelioration models to explain the grading.

Detection of recent HTs

It is difficult to assess the relevancy of proposed donors, because genes detected as potential HT have generally undergone amelioration (8). The comparison of recently diverged genomes (species or strains) provides the opportunity to find recent HTs, for which corresponding homologous genes in the donor species may be detected (52). Such a study is performed for five E.coli strains (two K12 strains: E.coli MG1655, E.coli W3110, one uropathogenic strain: E.coli CFT073, two enterohaemorrhagic strains: E.coli O157-H7 RIMD 0509952, E.coli O157-H7_EDL933) and two Shigella flexneri strains (S.flexneri 2a 2457T, S.flexneri 2a 301). These seven strains/species have recently diverged, genome sizes are different and the proportion of horizontally transferred genes varies from one strain/species to another (14,52). For instance, only ∼40% of the non-redundant set of proteins is common to E.coli strains CFT073, 0157-H7 EDL 9333 and MG1655 (53). These strains/species can be clustered in four groups with respect to phylogeny (Table 5).

Two criteria are used to searching for ‘recent horizontally transferred genes’: atypical regions (window size 1 kb, step 0.5 kb) (i) must have a signature that differs greatly from that of the host [distance to host must be at least >325, 2.5 times the E.coli intrinsic mean distance (Table 1)] and (ii) must be present in a limited number of strains/species to ascertain their recentness. In fact, outliers meeting the first criterion generally aggregate into several heterogeneous clusters (K-means clustering) that usually include samples from each strain/species. In some instances, however, some strains/species were absent from the cluster. It was subsequently considered that the corresponding regions might have been recently acquired by the relevant strains/species.

Table 5 shows a selection of potential recently transferred genes. Each cluster of atypical regions contains genes present in a specific set of strains. Some atypical genes are strain-specific, some are only absent in the non-pathogenic K12 strains and intermediate situations are also encountered. FASTA and Blast searches confirm that these genes are absent from some of the tested strains as already observed in the analysis complete genomes (53–55). In a large number of cases, we are able to find a well-conserved homologous gene in another species (Table 5). It is interesting to note that some of the suggested donors using our 12 000 signature database are in agreement with the species found by alignment methods. When no homologous gene is found, the proposed donors give credit to the known mechanisms of gene transfer (bacteriophages or plasmids) (Table 5).

It is worth noticing that most of the selected genes that are absent in K12 strains are involved in the pathogenicity of the other strains (52). E.coli 0157-H7 is the strain exhibiting the greatest number of genes absent in K12 strains [about 1400 (54)]. It has the greatest number of genes for which no homolog can be found (Table 5). Moreover, we are unable to propose a donor for a great part of these genes (Table 5). Many selected genes for E.coli 0157-H7 lie in the Ter region of the genome (between positions 2 000 000 and 2 500 000) in agreement with the published results (56).

Putative HTs

We have observed that most genomic regions are typical of the genome they belong to, using the signature as endpoint. Considering that the genomic signature is species-specific, atypicality of a region in terms of oligonucleotide usage has been promoted as a criterion for the detection of HTs. However, atypicality-based methods suffer several caveats that reduce their effectiveness in such a way that only a part of HTs can be detected. In fact, transfers between species with close signatures cannot be detected: significant differences between characteristics of transferred DNA and recipient species DNA are required. For similar reasons, HTs that were drastically ameliorated following their introduction cannot be detected either (8,14). The most stringent constraint, however, results from the size of the screening window. On the one hand, ideally, the best signal-to-noise ratio would be obtained when windows and HTs have a comparable size. On the other hand, the window size must be large enough to provide significant word counts, a requirement that strengthens with the size of the words under consideration and the intrinsic variability of the genomic signature along the genome. All together, the trade-off that has been implemented in this paper allows detecting atypical regions as small as 1 kb. In fact, rRNA regions sharing this characteristic were consistently detected. It must be pointed out that smaller fragments can be eventually detected if their signatures are radically atypical.

G+C% atypicality has often been considered as criterion for detecting HTs (8,24), but this approach suffered several drawbacks (28–30). It is to be noted that our signature-based method detects regions for which the G+C% lies within one standard deviation from the mean G+C% of the species (for instance, regions 2675251–2676250 in B.subtilis or 534751–535250 in H.Influenzae, see also Supplementary Materials 2 and 4).

As already observed by Nicolas et al. (12) for B.subtilis, rRNA has definitely an atypical signature. It is systematically classified as outlier, whatever the species (Table 1). Although transfer of rRNA from one species to another is unlikely (11,57), it cannot be firmly ruled out. However, it is clear that the atypical signature of rRNA does not imply that they are horizontally transferred.

The signature approach has an interesting property (that it shares with HMM) (7,12,28): detection is not bound to any specific function in the genome. In contrast with most other methods, the signature approach not only detects genes, but whole transferred regions as well, in agreement with the described mechanisms of DNA exchange between species. It is to be noticed that the method allows detecting several atypical non-coding regions (Supplementary Materials 3, 5 and 7). One major difference between HMM and signature method lies beyond the time required for the learning process, in the few resources that HMM can mobilize to deal with a short ‘one of its kind’ HT. On the other hand, HTs shorter than 1 kb can hardly be detected by a signature-based approach. An innovative HT detector is likely to result from an adequate fusion of both methods.

Potential donors

Several factors contribute to the efficiency of the search for donors. Of course, distance between putative HT and donor signatures is essential. Accuracy of signatures, linked to the length of available sequences, density of signatures in the ‘vicinity’ of HT, amount of amelioration sustained by HT during its presence in the host are also of importance [P. Deschavanne, S. Lespinats and B. Fertil, unpublished results; (25,27,31)]. Distance between the signature of a putative HT and the closest species varies to a large extent, but usually the shortest ones fall within the intra-genomic range (Table 1, Supplementary Materials 1, 2, 4 and 6). In some cases, the distance between the closest donor signature and the atypical segment signature is so great that no potential donor can be proposed (Supplementary Materials 1, 2, 4 and 6).

When strong similarities between a given DNA sequence and a foreign species are observed, the hypothesis for an underlying transfer is highly strengthen. However, the ‘true’ donor has to be previously sequenced and included in our bank of signatures to allow such a situation to occur. Moreover, we must take into account the intrinsic variability of short DNA segment signature (which is a function of their size, but also species-specific) when compared with the signature of a complete genome or any other large species sample (25,27,31). In the present state, our signature database is in no way representative of the diversity and richness of life. However, it must be noticed that there is already an obvious structure (in terms of distances between signatures) expressing taxonomy relationships between species in our signature database (31,58–61). Related species are often found close to one another. Clusters of potential donors may consequently provide pertinent information about the origin of HTs.

The diversity of signatures of putative HTs that can be observed for most of the species analyzed in this paper reveals the multiplicity of transfer events and donors (Supplementary Materials 2, 4 and 6). However, several outliers, not necessarily neighbors in the genome, are given the same set of potential donors (Table 1, Supplementary Materials 1, 2, 4 and 6). In general, the potential donors belong to few sets of taxonomically close species (Table 1) and share the biotope of the host (Supplementary Materials 1, 2, 4 and 6). For instance, B.subtilis, H.Influenzae and E.coli live in distinct biotopes; their potential donors do so as well. It is particularly encouraging to find that most of the potential donors that our approach has pointed out have had the opportunity to exchange DNA material with the recipient species.

Numerous viruses and plasmids qualify as potential donors (Tables 1 and 5, Supplementary Materials 1, 2, 4 and 6). It is not really surprising since they are known as HT vectors. They are often totally or partially inserted together with transferred genes in the host genome (14).

Some atypical DNA segments are particularly peculiar. They are isolated, have a specific signature (distances from neighbors are great), so that they cannot be given a credible set of donors (Supplementary Materials 1, 2, 4 and 6). Lack of data in the search domain, shift of signature features after a substantial amelioration process, structural constraints serving special functions or roles (14,62) (as it is for rRNA coding regions) are some of the tracks that remain to explore in these circumstances.

It would be interesting to localize the region the transfer may come from when the complete genome of the donor is available. However, homology (at the DNA level) is not a pertinent criterion for the comparison of sequences as soon as amelioration has taken place (8,14). In fact, homology is sometimes weak, e.g. between genes of Escherichia and Salmonella although these species have ‘recently’ diverged (34). It is clear that a more powerful search for the origin of putative HTs would have to embody models of amelioration [such as the one designed by Lawrence and Ochman (8)].

When searching for very recent horizontally transferred genes, in different strains of a species for instance, it was possible to find a great homology between detected genes and some genes from other species (Table 5). In numerous cases, the selection of donors is consistent with FASTA results (Table 5). This confirms the pertinence, beyond the similarity of signature between putative HTs and donors, of the proposed method to retrieve the species of origin of a transferred region. It seems that the search for origin of HTs on the basis of genomic signature is a powerful approach to understand some of the mechanisms of evolution (13,63).