On the reconstruction of the ancestral bacterial genomes in genus Mycobacterium and Brucella

Abstract

Background

To reconstruct the evolution history of DNA sequences, novel models of increasing complexity regarding the number of free parameters taken into account in the sequence evolution, as well as faster and more accurate algorithms, and statistical and computational methods, are needed. More particularly, as the principal forces that have led to major structural changes are genome rearrangements (such as translocations, fusions, and so on), understanding their underlying mechanisms, among other things via the ancestral genome reconstruction, are essential. In this problem, since finding the ancestral genomes that minimize the number of rearrangements in a phylogenetic tree is known to be NP-hard for three or more genomes, heuristics are commonly chosen to obtain approximations of the exact solution. The aim of this work is to show that another path is possible.

Results

Various algorithms and software already deal with the difficult nature of the problem of reconstruction of the ancestral genome, but they do not function with precision, in particular when indels or single nucleotide polymorphisms fall into repeated sequences. In this article, and despite the theoretical NP-hardness of the ancestral reconstruction problem, we show that an exact solution can be found in practice in various cases, encompassing organelles and some bacteria. A practical example proves that an accurate reconstruction, which also allows to highlight homoplasic events, can be obtained. This is illustrated by the reconstruction of ancestral genomes of two bacterial pathogens, belonging in Mycobacterium and Brucella genera.

Conclusions

By putting together automatically reconstructed ancestral regions with handmade ones for problematic cases, we show that an accurate reconstruction of ancestors of the Brucella genus and of the Mycobacterium tuberculosis complex is possible. By doing so, we are able to investigate the evolutionary history of each pathogen by computing their common ancestors. They can be investigated extensively, by studying the gene content evolution over time, the resistance acquisition, and the impacts of mobile elements on genome plasticity.

Background

Mycobacterium tuberculosis (MTB) is the etiologic agent of human tuberculosis (TB), that is one of the oldest recorded human afflictions which is still among the main worldwide death causes. In 2015, more than 10 million people became ill with TB and approximately 2 millions died from the disease, almost exclusively in low and middle income countries. Moreover, it induces a major global health problem, since about one-third of the world’s population has latent TB. Hence this is the first infectious disease declared by the World Health Organization (WHO) as a global emergency. More precisely, tuberculosis is caused by pathogens belonging to the Mycobacterium tuberculosis complex (MTBC) which consists of different species that are typical human pathogens (Micobacterium canettii, africanum, and tuberculosis), rodent ones (M. microti), or even Mycobacteria with a large host spectrum like bovis [1, 2]. Even if these organisms are genetically similar, they exhibit large differences with regard to epidemiology, pathogenicity, and host spectrum. Mycobacterium tuberculosis spreads throughout the human population since thousands of years, as the TB form that attacks bone and causes skeletal deformities can be still identified on individuals who died from it several thousands years ago, like ancient Egyptian mummies with apparent tubercular deformities.

The MTBC species are classified in 6 phylogenetic lineages which can be further divided into sublineages showing phenotypic differences reflecting for example their virulence (pathogenicity). The species members of the Mycobacterium tuberculosis complex have a clonal structure with large genome similarity (more than 99.9 percent of DNA sequences in common [3]). Compared to more ancient species, this complex has more virulent chromosomes [4]. As they have the same ancestor [5], the fact that we can find rodent and human pathogens, and other with a larger spectrum, is indeed surprising. To study M. tuberculosis DNA sequence, its virulent laboratory strain M. tuberculosis H37Rv is commonly used. This strain consists of a single circular chromosome composed by 4,411,532 nucleotides and 3906 protein genes. DNA homology studies and comparison of 16S rRNA coding regions have permitted to establish how they are linked, showing a 95−100% DNA relatedness. For example, there is only one difference between the 16S rRNA gene sequence of M. tuberculosis and the one of M. bovis.

The long-term coevolution of Mycobacterium tuberculosis with humans [6] has led to a more or less large geographic spread of the different phylogenetic lineages of MTBC. Moreover, some of the lineages appear to have a large geographic distribution, while others seem to be restricted to a smaller group of human host populations. Over time, MTBC genomes have evolved through genomic repetition or replacement (insertion sequences, etc.) and genomic modification at different scales of complexity. In this latter case, modifications range from small-scale ones resulting from mutation or indels to larger ones occurring on DNA strands (inversion, duplication, or deletion).

Obviously, understanding the past and future evolution of the MTBC would be of great interest, leading to the ability to study the ancestors and to understand the evolution history of species, and finally to an improved knowledge of the mechanisms of resistance and virulence acquisition in human tuberculosis. Fortunately, the relatively short time-frame during which the MTBC emerged (this bacteria is quite recent [7]), the relatively low genomes lengths and the recombination scarcity, together with an easier access to ancient and current DNA sequences, are favourable factors to address this question. Therefore, it should be possible to design a model of evolution for this set of genomes, in order to recover their evolution history and to predict their future evolution.

Another interesting group of pathogenic bacteria to be investigated is the genus Brucella which causes Brucellosis, a disease that primarily affects animals, especially domesticated livestock, producing abortion and other reproductive disorders. Human can also be infected, mainly through animal-to-person spread, in which case long-lasting flu-like symptoms are observed. Like tuberculosis, brucellosis is a global problem, since it is the most common bacterial infection spread from animals to humans worldwide. After the recent identification of the species B. vulpis, a total of eleven species have been identified within the genus Brucella according to their pathogenicity and preferential animal host [8, 9], among which the six classically recognized species are: B. melitensis, B. abortus, B. suis, B. ovis, B. canis, and B. neotomae. B. abortus and B. melitensis are the most important species regarding prevalence and morbidity in humans and domestic animals.

Clearly, a detailed knowledge of the Brucella phylogeny would also be of great interest. First, the phylogenetic reconstruction can lead to an enhanced understanding of the ecology, evolutionary history, and host relationships of this genus. Second, it can be used to discover suitable genotyping methods for rapid detection and diagnostic measures, used for example in epidemiological studies to facilitate human disease research. Moreover, as the Brucella genus is highly conserved and has low genetic variation, the phylogenetic reconstruction is still a challenge, even if the Brucella genus is probably easier to tackle than the MTBC.

This requires the development of new algorithms for the detection and evolution of genomic changes. Researchers studying this question focus mainly on the nucleotidic mutations prediction, and take specific forms for the matrix of mutations that seem not in accordance with recent experimental evaluations, see [10]. These evolutionary models must be constructed in a different manner, to better reflect what really occurred. Moreover, the important effects of other genome changes (such as nucleotide insertions and deletions, large-scale recombination, or repeated sequence changes) have to be considered more deeply, and an effective ancestral reconstruction of ancient bacteria should be carried out.

This research work is an extension of an article presented to the 5th International Work-Conference on Bioinformatics and Biomedical Engineering (IWBBIO 2017, [11]). Its main objective is to show that, if we focus on strongly related bacterial chromosomes, the reconstruction of their most recent common ancestors is possible in practice. In order to do so, we propose a pragmatic approach that mix already published reconstruction algorithms with new original scripts and a human cross-validation. As an illustrative example, we provide the ancestral reconstruction of 65 genomes of the Mycobacterium tuberculosis complex, and of the 47 Brucella genomes that are available on the NCBI database.

The dynamics of the evolution process in DNA sequences results from local evolutionary events that consist in SNPs or indels. Genomic rearrangements, which are larger alterations of the genetic organization, can take the form of inversions and transpositions, or occur by chromosome fusion and fission. Obviously, over time such large-scale mutations have affected gene order and content, therefore they have a prominent role in speciation [12]. A key problem when studying evolutionary change at the level of a DNA sequence, which is investigated by the research work presented in this article, is the problem of ancestral sequence reconstruction. This one is as follows: given an evolutionary tree relating organisms and the DNA genomic sequences of the leaf species, predict the DNA sequence of all ancestral species in the tree. Many biological studies have addressed this problem and thus various methods have been proposed for inferring ancestral sequences. Apart from ancestral genome reconstruction problem, biomolecular evolution is usually devised through the evolution of core and pan-genome. Below is a brief overview on ancestral genome reconstruction.

Similarities in sequences or in the gene order (genome composition) are usually considered in up-to-date ancestral reconstruction methods. The first case, based on sequence similarity, can be considered as resolved now, at least when indels are not considered [13–21]. Indeed, considering a phylogenetic tree and its associated DNA alignment, Bayesian inference or maximum likelihood approaches can be applied to estimate ancestral states of nucleotides [22, 23]. The main problem is the insertion-deletion case, which is usually disregarded [24]. The small number of models that consider indels focus on the parcimony approach, or consider the evolutionary model called Thorne-Kishino-Felsenstein [25]. Combinatorics investigations are applied in the case of larger modifications, by modeling these recombinations as permutations of homologous sequences. This reformulation leads to the well-known genome rearrangement problem [26], in which the shortest edit operations that can map one chromosome to another are searched. Note that this NP-hard problem [12, 27] is directly related to the sequence length and the number of mutations, while genomes considered in this article are quite small and have faced only a low amount of recombination: the difficulty can be circumvented for such genomes.

The remainder of this article is organized as follows. The methodology proposed for ancestral reconstruction is detailed in the next section. Results of the application of this approach on the Mycobacterium tuberculosis complex (specifically on two of its species, namely M. tuberculosis and M. canettii) and on the Brucella genus case (focusing specifically on the B. abortus and B. melitensis species) are investigated in the third section. Finally, this research article ends with a discussion and a conclusion with future work.

Methods

Let us now detail our concrete ancestral reconstruction for bacterial genomes, illustrated through a first set of strains detailed hereafter.

Data acquisition and processing

A python script has firstly been written to automatically download all the complete genomes of Mycobacterium genus available on the NCBI database, encompassing 2 africanum, 15 bovis, 5 canettii, 1 microti, and 42 tuberculosis. Note that canettii and tuberculosis are well represented in this dataset, which is helpful to study how virulence has appeared in the first species, and if the second one is at the origin of the MTBC complex 40,000 years ago. Details about these 65 genomes are provided in Table 1.

Table 1 The considered Mycobacterium strains

After the data acquisition stage, the next step is to align the downloaded sequences [28, 29]. Prior to the Multiple Sequence Alignment (MSA), genomes must be operated such that each sequence starts to the same location and is read in the same direction: we deal with circular genomes. This is why a sequence of reference (200 bp from M. tuberculosis H37Rv) and its reverse complement have been blasted locally. Then, a circular shift and/or a reverse complement of the whole sequence have been applied when required.

Most of the well-known alignment tools have failed to align these genomes, due to their size, while we do not want to split the sequences, to reduce the complexity of the alignment, as this multiplies the intermediate steps, increasing by doing so the risks of errors. It was not the case of AlignSeqs, available in the R module called decipher [30]. This latter achieved to perform the MSA in an accurate and rapid way. With this tool, multiple sequence alignments are done by aligning 2 genomes first, and then adds a third genome, etc., until all the sequences are aligned [31].

Phylogeny

The alignment of multiple genomes of Mycobacterium leads to the visualization of synteny blocks, emphasizing the location of large inversions.

A manual reverse of these inversions were possible, leading to an improvement of the alignment of the 65 genomes. This is beneficial for the next stage of the pipeline, namely the phylogenetic investigation. This stage has been performed using RAxML, in which the phylogenetic tree is reconstructed according to a maximum likelihood approach [32]. Note that, thanks to the manual reverse of inversions, the obtained tree has been computed using almost all the complete genomes (only columns with indels are ignored), while without this manual operation, all columns inside the inversion are disregarded. Being based on almost all the genomes, and being strongly supported according to bootstrap values, the obtained tree is trustworthy, and we can reasonably consider it as a backbone to reconstruct ancestral states of MTBC nucleotides.

The proposed ancestral reconstruction is in two parts: 1-length modifications (SNPs and indels) are first considered, before investigating larger modifications (insertion, deletion, or duplication of large scale subsequences). These two case are detailed below.

Ancestral reconstruction: the mononucleotidic variants case

The treatment is divided in two sub-parts: insertion-deletions on the one hand, and single nucleotide polymorphisms on the other hand. The second case is simple, and its difficulty is only in the separation between real SNPs and polymorphism induced by an indel recombination. The first case is more complicated, as indels may be related to mobile elements or tandem repeats. These two cases are detailed below.

Ancestral reconstruction of SNPs is realized as follows. We first compute the marginal probability distributions in each nucleotide of internal vertices in the phylogeny obtained previously. Assuming a site independence, we have applied the sum-product message passing method [33] to calculate these distributions. This method has been applied by using PHAST [34], which is able to reconstructs ancestral indels too (parsimony approach).

Ancestral reconstruction: the case of larger variants

In the case of mid-size modifications over time, a string algorithm has been first designed to detect sequence inversions (even in the case of small and noisy ones). However, and due to the fact that MTBC complex is reputed to evolve in a clonal manner, only artifacts have been detected by applying this algorithm on supercomputer facilities. This will not be the case if this pipeline is applied to more recombinating bacteria like the Pseudomonas or Yersinia genus. Note that, up to now, duplications have not yet been regarded, as the synteny block analysis performed previously has shown that large scale duplications have not occurred in the MTBC case.

Conversely, midsize indels and SNPs have been investigated in details by using PHAST. This investigation has allowed us to notice that: (1) In most of the cases, the situation is obvious, leading either to a deletion or an insertion at a well specific location inside the phylogenetic tree, like in Figs. 1 and 2. (2) These larger variants events are rare in various lineages (e.g., tuberculosis), as illustrated in Table 2. (3) In the case of indels of size ≥ 2, the parsimony approach of PHAST produces frequently a wrong ancestral state deduction, which must be modified by hand. Note that its competitors have been tested too, and they all presented worse reconstructions on our specific dataset. (4) The inserted sequence has, in general, not faced additional mutations over times.

Fig. 1

Indels on internal nodes of the tree of some M. canettii species

Fig. 2

Ancestral reconstruction of one problematic indel in the alignment

Table 2 Single nucleotide polymorphism between species (100.X is the name of an ancestral species, cf. the phylogeny)

This semi-automatic pipeline for ancestral genomes has finally succeeded to reconstruct the genomes at each internal node of the tree, which can be done because the number of recombination of more than one nucleotide is low. These recombinations have mainly been deduced manually, while state-of-the-art tools have not been able to reach an acceptable level of accuracy.

Figure 3 summarizes all the ancestral reconstruction process, in which the gray boxes are operated manually, while the other stages are automatic. Indeed, obtained results on mononucleotidic variants have been carefully checked by naked eye, as the number of such variants is lower than one hundred, while ad hoc algorithms were designed to deal with variants of larger size, see Fig. 4.

Fig. 3

Flowchart of the proposed approach

Fig. 4

Ancestral reconstruction of a M. canettii SNP

CRISPR investigation

Another particular DNA pattern that can evolve through Evolution is the so-called CRISPR one. CRISPR refers to repeated DNA sequences that help to preserve organisms from noticeable threats like viruses. These sequences are a fundamental component of some immune systems, which helps to protect their organism’s health. Such repeated DNA sequences are found in archaeal and bacterial genomes. These sequences range in size from 23 to 47 base pairs.

The name of CRISPR refers to an acronym which stands for Clustered Regularly Interspaced Short Palindromic Repeat [35, 36]. The CRISPR system was initially found as part of an immune system of sorts in some bacteria, used for cutting apart foreign DNA. It consists of two parts of the protein itself, which is the workhorse of the CRISPR system: a bacterial enzyme named Cas9, and a small RNA, called the guide RNA, that matches the DNA sequence to be nicked [37].

Results

The Case of Mycobacterium Tuberculosis Complex

All the 65 Mycobacterium genomes have been aligned thanks to the AlignSeqs function described previously. We thus obtained a first representation of synteny of all of them, see Fig. 5. As can be seen, genomes are very similar in the MTBC case, and only a low number of recombinations have occurred within these genomes.

Fig. 5

Synteny blocks of Mycobacterium strains available online

As an illustrative example of the phylogenetic study depicted in Sec. 3, the phylogeny of M. canettii is represented in Fig. 6 (outgroup: M. tuberculosis). We selected the GTR Gamma model of nucleotide substitution as recommended by JModelTest 2.0, and the tree has been computed by RAxML. Note that the obtained tree is well-supported, as well as in the M. tuberculosis cased, whose supports are larger than 98% (cf. Fig. 7). Indeed, with these bacteria, we have not to find the most supported tree based on the largest subset of core genes, as aligning the whole complete genomes leads to a well supported tree: it is not possible to improve the results, which is nice as the core genome is many times greater than in the chloroplast case (and so, it is not sure that the heuristic approach presented in our previous articles [32, 38, 39] can succeed to find the optima).

Fig. 6

M. canettii phylogeny (outgroup: M. tuberculosis)

Fig. 7

M. tuberculosis phylogeny (GTR Gamma model and outgroup:M. africanum)

The obtained results on mononucleotidic variants have been humanly verified, which has been possible due to a low number of variants (cf., for instance, to Tables 3 and 4).

Table 3 Number of columns of the MSA with SPNs or indels for M. canettii (large deletions are counted character by character)

Table 4 Variations in the alignment of the M. tuberculosis clade under consideration

166 indels and 2956 SNPs have finally been detected, when considering the 5 M. canettii (see Fig. 8). Figure 9, for its part, collects the positions of the 25 indels and 394 SNPs that have been detected in the clade of the 8 M. tuberculosis.

Fig. 8

SNPs location of mononucleotidic variants of M. canettii

Fig. 9

SNPs location of mononucleotidic variants of M. turberculosis

In the considered Mycobacterium strains, only a few important inversions have been detected, such as the inversion present in the last ancestor of 140070010, CIPT 140010059, 140070017, 140060008, and 140070008, as shown in Fig. 10. 99% of DNA sequence identity has been obtained when considering all the blocks of synteny of tuberculosis. We can conclude that these genomes are highly conserved: highly similar regions without any rearrangement, with only small indels and a large inversion.

Fig. 10

Synteny blocks in M. canettii. Each genome is colored according to the position of the corresponding region in the first genome (gray if a region is unshared)

We can conclude from this study that ancestral genome reconstruction is possible when considering close or clonal bacteria, and all the material needed in such a pipeline has been designed. But, for the sake of comparison, it may be interesting to deep investigate the results of this semi-automatic reconstruction method on a quite more stable genus, namely the Brucella, on which human validation of algorithm results is easier (see Tables 5, 6 and 8 for an illustration of their alignment and SNP differences). Such new investigations are conducted in the next section.

Table 5 Differences in the alignment on chromosome 1 of abortus

Table 6 Single nucleotide polymorphism in Brucella melitensis

The Case of Brucella genus

The pipeline presented in the previous section is now applied on another genus, namely the Brucella one, for the sake of comparison and to broader the discussion. Complete sequences of the 47 available genomes have been downloaded from NCBI, namely by species: B. abortus (14 genomes), melitensis (8), sui (16), ovis (1), canis (3), ceti (2), pinnipedialis (2), neotomae (0), microti, inopinata, and vulpis, as described in Table 7.

Table 7 Brucella genus: genome information

Table 8 Single nucleotide polymorphism in Brucella abortus

Note that the genome of Brucella abortus has two circular chromosomes. The first one is 2,124,241 bp long in the Brucella abortus biovar 1 str. 9-941 reference genome, while the second chromosome is of 1,162,204 bp. Other species in the Brucella genus are comparable in genome size. For instance, the Brucella melitensis strain 16M is constituted of 3,294,931 bp disseminated in two circular chromosomes: chr. I has 2,117,144 bp, while chromosome II has 1,177,787 bp. On both of these chromosomes, approximately 3100 ORFs were predicted. In the latter, genes encoding for DNA replication, protein synthesis, core metabolism, and cell-wall biosynthesis can be found on both chromosomes [40, 41].

We operated the sequences so that they share the same orientation (which may need a transconjugate operation) and the same sequence of 200 nucleotides as starting point (which may require a circular shift), if we except local SNPs. This has been achieved using a local blast, with the beginning of Brucella abortus 2308 as an arbitrary reference. After such operations, a syntheny representation of Brucella genomes can be obtained, as shown in Fig. 11. The particular case of B. abortus is depicted in Fig. 12.

Fig. 11

Brucella, chromosome 1: a high sequence similarity with little recombination events

Fig. 12

Synteny map of Brucella abortus (a) chromosome 1 and (b) chromosome 2. Genomes investigation tends to show a high sequence similarity with little recombination events. Each genome is colored according to the position of the corresponding region in the first genome, or gray if a region is unshared

A few inversions have appeared in this representation. For instance, in the B.abortus case, we found a significant inversion at the last common ancestor of strains“biovar_1_str._9-941”, S19, A13334, “strain_BDW”, “bv._2_str._86/8/59”, and 104M. We have manually reversed these inversions, so that an accurate alignment of the whole genomes can be performed. Using this alignment, a very well supported phylogenetic tree has been obtained. For the sake of illustration, a subtree corresponding to the phylogeny of the Brucella abortus species is depicted in Fig. 13, and in Fig. 14 for B. melitensis. It has been obtained using the entire genome sequences with RaxML, GTR Gamma model, and Brucella melitensis as outgroup. As can be shown, all branches exhibit a 100% bootstrap support value.

Fig. 13

Well-supported phylogeny of Brucella abortus species calculated on the entire chromosome 1. The outgroup is melitensis, while RaxML has been launched with the GTR Gamma model

Fig. 14

Well supported phylogeny of Brucella melitensis species

At this stage, all the material required to attack the ancestral reconstruction of Brucella genomes are on hand. We first have focused on the abortus and metilensis reference species, to investigate the potential origin and the history of the global spread of these Brucellas. We have considered the global alignment of both chromosomes 1 and 2 of the available complete strains, using decipher R package [42], and the tree depicted in Figs. 13 and 14. We firstly achieved a comparative whole-genome single nucleotide polymorphism analysis of these strains collected and downloaded from the NCBI. 32 indels and 373 SNPs have been detected in the clade containing these 6 variants of chromosome 2, and 609 SNPs and 325 indels in chromosomes 1, as shown in Fig. 15. The same has been computed for B. melitensis, leading to 6178 variants and 335 indels, see Fig. 16. This has been achieved using homemade python scripts on aligned sequences.

Fig. 15

SNPs location in Brucella abortus species. (a) Chromosome 1, (b) chromosome 2

Fig. 16

Single nucleotide polymorphism in Brucella melitensis species

At mononucleotidic variant level, the treatments of SNPs and of indels have been separated. Examples of mononucleotidic ancestral reconstructions are provided in Fig. 17. Differences between ancestors and their children are, for their part, provided in Tables 6 (abortus) and 8 (melitensis).

Fig. 17

Nucleotides in the ancestral nodes and their children on Brucella abortus species. a Chromosome 1 b chromosome 2

Figure 12 shows homologous regions among many Brucella abortus genomes, as identified by FindSynteny (R). On the one hand, the similarity and preservation of synteny blocks on Brucella abortus are especially pronounced in chromosome 1, with highly similar regions and without rearrangement of homologous backbone sequences as shown in Fig. 12a. Chromosome 2, on the other hand, is more diverse. There is above all a significant reversal in the Brucella abortus genomes of the clade consisting of abortus 0, 1, 2, 4, 10, and 12 as shown in Fig. 12b. The same information is provided for B. melitensis (chromosome 1) in Fig. 18. These differences most likely represent distinct evolutionary origins, for instance related to the nature of functional genes in the two chromosomes.

Fig. 18

Dotplot of Brucella melitensis species, chromosome 1

We finally analyzed the CRISPR locus sequences of 14 Brucella abortus strains by using CRISPRs web service (http://crispr.i2bc.paris-saclay.fr). The orthologous sequence shared between Brucella abortus genomes and the CRISPR spacer have shown a significant similarity of the spacer sequences. Figure 19, for its part, shows the CRISPR space sequence lengths and their positions inside abortus genomes. For the B. melitensis case, information are provided in Fig. 20.

Fig. 19

Brucella abortus phylogenetic tree: estimation of the CRISPRs length and locations by using the CRISPRFinder web server [36]

Fig. 20

CRISPR investigation in B. melitensis

Discussion

Various algorithms and methods can be found in the literature to resolve, at least partially, the ancestral genome reconstruction problem. We have shown that these existing methods are not accurate and mature enough to be applied on a real case scenario. This is particularly evident when indels or single nucleotide polymorphisms are mixed with repeated sequences. The main drawback of these methods is that they intend to solve all the cases, while some situations are up-to-now too difficult to be resolved automatically. However, in mid-size genomes that have faced a low number of recombinations over time, as for Brucella and Mycobacterium, these problematic situations can be signaled, and a human cross-validation can reinforce the accuracy of the ancestral reconstruction algorithm.

As a proof of concept, all ancestral genomes of all M. canettii available on the NCBI database have been reconstructed, as well as all the ancestors of the available M. tuberculosis complete genomes. At each time, the single nucleotide polymorphism level has first been investigated, before considering the cases of indels and large scale recombination.

Obtained results show that a concrete and accurate reconstruction can be achieved by coupling human decisions on problematic situations with automatic inference of ancestral states in easy to resolve ones, at least for some non recombinant bacteria. With such a reconstruction, it may be possible to deeply investigate the evolution of genomes over time, and possibly to predict their future modifications.

Conclusion

In this article, we presented a semi-automatic pipeline that achieves to completely and accurately reconstruct the ancestral genomes of some clonal bacteria. In this pipeline, the case of SNPs and indels of 1 nucleotide has been resolved using the sum-product message passing algorithm, while larger modifications have been studied by a parsimony approach coupled with a manual deduction.

The obtained ancestors have not yet been investigated in this study, as it was not the objective of this proof of concept. They will be studied with ad hoc algorithms to design, to investigate the evolution of gene content on the one hand, and of mobile elements on the other hand [43, 44]. The rate at which such loss or gain occurs will be examined carefully, and we will study if some particular functionality are more affected by these mutations. To say this differently, we will investigate if modifications have a real impact during the evolution of genomes.