The evolution of bacterial genome assemblies - where do we need to go next?

CURRENT SEQUENCING TECHNOLOGIES

Our ability to decipher DNA, the building blocks that define the state and evolution of life, has fundamentally changed our understanding and interaction with living beings.

It is becoming commonplace in many laboratories to perform whole-genome sequencing and create a draft genome assembly as the first step of working with a new organism. It was not very long ago where that idea would have been considered cost-prohibitive and unrealistic. Genome sequencing technologies have come a long way in a surprisingly short time. Advances in sequencing technologies and computational pipelines for genomic data reconstruction have made these experiments tenable for non-experts and have allowed scientists to pose questions that were not even imaginable a few years ago. Now in the era of “big data”, each new advancement is pushed along two fronts. The first one is the sequencing technologies themselves in a never-ending battle for more data with increased accuracy using poorer samples. The second is the computational pipelines to handle this information to get a better understanding of what the data is telling us.

Although Sanger receives much of the credit for modern DNA-sequencing, Wu and Taylor^[1] were one of the first to report the determination of a DNA sequence; however, many of the techniques used were adapted from previous work from Sanger^[2-4] and Holley et al.^[5]. Further DNA-sequencing approaches were developed^[6-8]; however, these tended to be rather laborious until Sanger published the now famous “chain-terminating” inhibitors approach^[9]. This technique has gone through many iterations and has been greatly improved, and the core technology remains widely known as “Sanger Sequencing”.

Throughout the 1980s, shot-gun sequencing approaches accelerated the pace of DNA sequencing^[10]. It was realized that fragments of a genome could be randomly digested and ligated into bacterial vectors. These vectors could now be readily sequenced thanks to Sanger-based sequencing approaches, and regions of overlap were used to place the DNA sequences for allowing the sequence of whole genomes to be determined. The idea that whole genomes could now be sequenced spurred The Human Genome Project (HGP).

The impact Sanger sequencing has had on the field cannot be understated, but, as with any scientific discovery, there is always a push for more. Sanger sequencing builds from a single primer present in the sample. Therefore, if more than one primer is present, or if the sample contains another DNA fragment that the primer can bind to, the correct sequence of interest cannot be determined. Also, only a single strand can be sequenced at one time. With many advances in technology and sample automation strategies^[10-12], this technique soon reached its peak performance. However, it should be noted that Sanger sequencing is one of the most widely used tools in laboratories today in terms of sequence verification and molecular biology applications but is less practical in the era of omics data.

It is estimated that the cost associated with the HGP was $2.7 billion USD. Scientists recognized the value of this technology, and companies understood that revenue could be generated by decreasing the costs to make it more accessible. Figure 1 highlights the average cost of sequencing (per MB) as compared to the amount of sequencing data being deposited into NCBI. As we see new technologies emerging, the cost is significantly decreasing. These new technologies are now collectively referred to as next-generation sequencing (NGS) technologies and can be further characterized into short-read (pyrosequencing, Ion Torrent, Illumina), long-read (SMRT: PacBio and Nanopore), and next-gen (next step in sequence evolution).

Figure 1. Cost versus sequences deposited. Sources: genomes.gov & NCBI. Amount of sequence deposited doubles roughly every 18 months and this is likely to accelerate as long-read technologies further mature. Red: Cost to sequence a Mbp of DNA; Black: individual sequences directly deposited into Genbank; Grey: WGS sequence data deposited (April 2002); WGS: whole-genome sequencing.

THE ASSEMBLY OF METAGENOME ASSEMBLED GENOMES - CONCEPTS, STRATEGIES AND CHALLENGES

A very promising aspect of high throughput sequencing in the field of microbiology was the possibility to study microbial communities directly on the genome level. The thousands to millions of generated sequence reads from just one NGS run could cover hundreds of microbial genomes. This increase in sequencing coverage compared to Sanger-based sequencing made it possible to obtain direct insights into the DNA composition of samples containing a multitude of organisms. A sequence data set obtained from such samples is called a metagenome. Metagenomes were studied both directly on sequence read level, or an assembly of the short sequence reads into longer contigs was attempted to obtain complete gene sequences or operons from the data. Metagenome analysis became a powerful strategy to study not (yet) culturable microorganisms^[24]. It was now possible to obtain information about so far unknown taxa and the potential functional properties of microbial genomes directly recovered from natural habitats^[25,26].

The sequencing coverage per NGS run continued to increase, and assembly tools were adjusted to fit metagenomic data, which eventually enabled the recovery of entire microbial genomes from sequence data of complex microbial communities. As shotgun sequencing evolved from elaborate and resource-consuming procedures to be much more time and cost-effective, so did the recovery of metagenome-assembled genomes (MAGs) (metagenome-assembled genomes) from the generated sequence data.

The first published MAGs were obtained from shotgun Sanger sequencing of a low complexity metagenome^[27], yet as NGS technologies advanced to be widely available, metagenomes were then mainly generated via Solexa/Illumina, 454 or Ion Torrent sequencing. However, the assembly of MAGs from these relatively short read data was more challenging than the assembly of longer reads. In particular, the short reads made it difficult or impossible to assemble repetitive sequences or regions encoding highly conserved genes with several copies present in the genome, such as ribosomal genes^[28,29]. Sequencing strategies such as mate pair sequencing and artificial long reads mitigated some of these drawbacks, but it made the sample processing and sequence generation, as well as the post-sequencing data processing steps, more elaborate and time-consuming.

Therefore, initially, a “MAG” was not a complete closed genome, but rather an assembly of several (up to hundreds of) contigs that were assigned to a taxon, often described as an operational taxonomic unit (OTU). Assigning contigs to an OTU or MAG usually involved a binning step, based on factors such as k-mer frequency, the abundance of assembled genes, alignments to reference databases, or a combination of several strategies^[30-32]. The quality and completeness of such MAGs were then determined by the presence of genes encoding RNAs and universal single-copy genes^[33,34].

In order to obtain high-quality MAGs or even closed genomes from metagenomic data, long read sequencing technology platforms described above lead to substantially improved assembly results. Parallel to this, software tools optimized for metagenome assembly and binning were developed, making it possible to assemble single contig MAGs from complex metagenomes^[35,36].

FUTURE CHALLENGES

Matured short- and long-read sequencing technologies have opened the doors for deep shotgun metagenome sequencing and, more recently, enabled the assembly of MAGs at larger scales. New scientific horizons have begun to emerge. With increasing sequencing depth and read length, sufficient coverage of rare members of microbial communities will become possible, enabling the creation of high-quality MAGs that will help to ascertain the role and contribution of these rare members to the overall community fitness. As described above, the quality of MAGs is ever-increasing, and we anticipate high-quality draft and closed MAGs to emerge routinely over the next five years. As these (complete) MAGs become increasingly available not only for major clades within a microbial community but also for increasingly rare members, we will also be able to begin capturing strain variants within a population, including their respective variant frequencies. It is a long established fact that even in homogeneous bacterial cultures originating from a single cell, sub-populations with genetic variations emerge over time. The same is true in natural microbiomes where ongoing and sometimes rapid evolution creates a reservoir of sub-populations with genetic and possibly phenotypic variations that are not negatively affecting their respective survival. These silent sub-populations enable rapid responses and adaptation to change environmental pressures and are a key driver to the survival of a given microbial species and strains. Without these sub-populations, the microbiome’s ability to quickly respond and adapt to any environmental pressures would likely be severely compromised. But these sub-populations on strain level have so far eluded recognition and (full) characterization, despite their importance for ongoing bacterial genomic evolution and lifestyle adaptations, impact on personalized medicines, and biotechnological and pharmaceutical fermentation^[55,56].

Maturing long-read and emerging next-gen^[57] sequencing technologies, that will again deliver longer read lengths at reduced error rates compared to current long-read technologies, will deliver on the promise of deep MAGs that will detail genomic variation in unprecedented detail. However, these new frontiers come with formidable and interlinked challenges.

The growing volume of sequence information originating from the deep community and single-cell sequencing will necessitate reliable and automatic genome assemblies paired with high throughput detailed genome annotations. Although advances in long-read technologies are predicted to generate reads that could span > 100 kb, making it easier to assemble closed genomes, this will be offset by the amount of sequence generated with the maturation of single-cell genome sequencing. Reliable and fully automated structurally correct assembly of genomes and genome closure will require more advanced use of artificial intelligence (AI) and machine learning. AI systems will need to be able to self-train structurally correct assemblies based on existing reference assemblies and deep machine learning algorithms capable of recognizing and reconciling structural genomic variations. As new assembly and annotation algorithms emerge, AI systems will autonomously test, evaluate and, if deemed suitable, incorporate these into their existing repertoire of tools. Ongoing improvement and optimization will likely require cloud-based systems capable of providing required computational resources. For end-users or confidential data, static snapshots will be spawned on a regular basis that can be installed and executed on less powerful, local workstations. It is hard to predict how quickly such systems will evolve and whether they will be deployed on conventional silicon or make use of quantum computers capable of solving more complex networks. While the prospect of quantum computing is tantalizing, its use is currently limited to very specific types of problems and programming languages, and concepts have barely begun to emerge, allowing for broader fields of applications^[58]. However, forecasts predict a much wider use of quantum computers over the next decade^[59], and their application in next-generation genome spaces, defined as the collective of (accessible) draft and complete genomes, is highly likely.

The second challenge is directly enabled by the emerging ability to investigate populations in detail on the genome level and thus capture genetic variances on a strain level, initially as a static snapshot and eventually including longitudinal changes in response to environmental changes. How can this variation be captured and recognized in a (human) readable format?

Historically and at present, genome assemblies of microbial populations are assembled based on determining a single consensus sequence. Sub-populations on the strain level that feature genomic variations are discarded. However, with the emergence of single-cell sequencing and the availability of high-quality draft and closed genomes assembled from communities and microbiomes, opportunities arise to investigate the contribution of such sub-populations to the rate of evolution, the fitness of a population, and their ability to respond to environmental challenges and pressures. Such a deep understanding of genomic variation below the strain level would also enhance our ability to predict the impact on health and well-being of humans and animals as we begin to better understand how beneficial and pathogenic microbes evolve, how genomic hotspots with rapid evolution emerge and contribute to phenotypic change, and at what frequencies such sub-populations are maintained or begin to disappear.

To begin addressing this emerging opportunity, we first need to consider thresholds upon which base call errors can be distinguished from the presence of subpopulations. This will likely be sequence technology-specific and, from today’s perspective, comprises a statistical approach combining known error rates with sequence depth coverage for at least one long-read and one short-read sequencing chemistry. This implies a minimum coverage needed for each genome position that enables not only the establishment of an error rate threshold but also allows for sufficient information to detect sub-populations.

Identifying sequence reads from metagenomes or (assembled) genomes from single-cell sequencing that belong within the same strain genome space is the next step in the evolution of genome assemblers. Here, the definition of a strain genome space requires some consideration and may need to be flexible enough to accommodate different evolutionary strategies deployed by different microbes. Perhaps a first model could consider the presence and absence of SNPs. In such a model, only minor changes between individual sequence reads and/or assembled genomes are permitted. The capture of such variations can then be readily realized by representing each genome position by five possible states: A, T, G, C, and Gap. Each state can be attributed a normalized distribution value, ranging from 0 to 100% of the total sequence coverage at this genome position. Consequently, the consensus genome for a strain genome space would appear to be slightly larger than individual sub-population specific genomes to accurately reflect insertion events present in individual sub-population genomes.

Current methods to record sequence information and annotation cannot effectively capture this genomic variation. Data formats such as FASTQ^[60] capture the quality of a base call at a given position with a sequence read, while FASTA represents consensus sequences without any quality qualifiers. Similarly, file formats such as Genbank can provide information on the first higher abstraction level - gene models and annotation - but cannot represent genomic sequence variation. Sub-population variation on strain levels could be represented by a variation of the FASTA format, as described above. A potential new FASTA file format, tentatively named “FASTV”, is shown in Figure 2 (FASTV). Each genome position is represented by either a one-dimensional matrix (FASTV A) or a more human-readable notation (FASTV B). The FASTV format captures positional variations but does not differentiate between individual subpopulations. This FASTV format could then be seamlessly integrated into a Genbank file format analog. Here, the way gene models (and hence annotation) are described also needs to be modified. One approach could be to allow multiple ORFs per locus tag qualifier, again assigning a normalized frequency to each locus tag associated with a given subpopulation.

Figure 2. Potential FASTA successor formats FASTV and FASTM. For FASTV, each genome position is depicted either by a matrix “()” (FASTV A) or a human readable notation (FASTV B), where a single genome location comprises the normalised distribution value of one of the four bases or a gap, ignoring bases that do not occur at this position. In FASTV, the frequency of nucleotides per position is captured but not their respective occurrences across subpopulations. For FASTM, a two-dimensional matrix enables capturing not only the relative distribution of nucleotide variations for every genome position, but the matrix further enables the tracking of individual subpopulations.

The main drawback of such a simplistic expansion of available file formats is that the sub-population structure is lost, and genetic changes within individual sub-population genomes cannot be traced.

Possible evolution of the FASTV file format that enables accurate sub-genome tracking could be a two-dimensional matrix file format where one dimension is reserved for genome position and the second dimension recognizes each identified sub-population. A tentative file format, dubbed FASTM, is suggested in Figure 2 (FASTM), where each matrix presents a single genome position. With each matrix, the first column represents identified sub-populations while the remaining columns depict the normalized frequency of each possible state in the order A, C, T, G, and - (gap). Therefore, the FASTM format is a representation of a sequence data set, with the sequence reads already aligned and assembled. Sequencing errors have already been checked and removed using the appropriate tools/statistics etc., assuming sufficient sequence coverage. The FASTM matrix then captures, with the required significance, different variants within a single strain or different strains within one species in their respective abundances. Similar to the FASTV format, FASTM matrices could be incorporated into Genbank-like file formats and sub-population-specific gene models indicated through locus-tag qualifiers [Figure 3].

Figure 3. Possible translation of sub-population sequence variations into a Genbank-like format. Sequence variations that do not change the gene model retain the same locus tag (e.g. LT_0234) but change the appendix to refer to the respective subpopulation. Where gene model changes occurs, e.g., sub-populations 2 and 4, new locus tags are issued along with the respective sup-population appendix.

These new concepts for recording and analyzing genetic variations on sub-strain populations are far from perfect, and we suggest a concerted effort of the science community to further refine these file format suggestions to be prepared for these upcoming challenges.

We also acknowledge the limitations of these new concepts and formats. While we believe that, in principle, they will be capable of accurately reflecting nucleotide changes, insertion/deletion events, and horizontal gene transfer events of any size, they are not suitable to capture genome rearrangements of any size. It may be a worthwhile discussion to define at which point observed genome rearrangements (or the loss and acquisition of loci) still represent a sub-population of the same strain and when these changes warrant the definition of a new reference strain. If genome rearrangements are allowable features within a strain sub-population, the suggested FASTV and FASTM format will no longer be suitable and should be replaced. Another limitation is the relevance these new formats will carry with respect to the quality of meta-genome assemblies. Ideally, these formats are applied to closed genomes or high-quality draft assemblies comprising long scaffolds or contigs. For assemblies of lower quality and thereby decreasing scaffold/contig lengths, it will be more difficult to confidently assign metagenome reads and detected sequence variants to individual assemblies and contigs therein and scaffolds/contigs to genomes. With ongoing advances in sequencing technologies and assembly tools, this limitation may only be a temporary one until high-quality assemblies and closed genomes have become the norm.

CONCLUSION

Ongoing advances in both sequencing technologies and assembly tools are paving the way for a deeper and more detailed understanding (and tracing) of evolutionary changes in the genetic blueprint of individual strains. Current data formats are not fully suitable for managing these new levels of detail, and here we have provided two initial formats that may begin addressing this challenge. We hope that these suggested formats may serve as a starting point for a larger science community discussion on how to best capture sub-population dynamics on a genome assembly level and identify which advances will be required in assemblers and annotation pipelines to leverage this new level of detail.

Looking further into the future, such new formats alongside advances in machine-learning driven assemblers may also serve to describe longitudinal changes in populations and their respective assembled MAGs, perhaps for the first time allowing to follow directly emerging genetic adaptation and population changes on a sub-strain level.