The genome atlas: navigating a new era of reference genomes
Authors: Alice M. Kaye, Wyeth W. Wasserman
Journal: Trends in Genetics (2021)

Current approaches to generating a new reference genome are misguided because they often conflate the two distinct roles of the reference: (i) enabling knowledge exchange; and (ii) increasing computational efficiency. This dichotomy needs to be reflected in the underlying design of a new reference.There is widespread acceptance that the current linear reference genome needs to be expanded to incorporate the common variation observed across different populations for both equity and performance considerations.Switching to methods that index DNA reads rather than reference genomes offers advantages by preventing compatibility issues between reference versions. Indexed read sets have the longevity to adapt as sequencing technology advances, read length increases, and the volume of data being generated surges. The reference genome serves two distinct purposes within the field of genomics. First, it provides a persistent structure against which findings can be reported, allowing for universal knowledge exchange between users. Second, it reduces the computational costs and time required to process genomic data by creating a scaffold that can be relied upon by analysis software. Here, we posit that current efforts to extend the linear reference to a graph-based structure while trying to fulfil both of these purposes concurrently will face a trade-off between comprehensiveness and computational efficiency. In this article, we explore how the reference genome is used and suggest an alternative structure, The Genome Atlas (TGA), to fulfil the bipartite role of the reference genome. The reference genome serves two distinct purposes within the field of genomics. First, it provides a persistent structure against which findings can be reported, allowing for universal knowledge exchange between users. Second, it reduces the computational costs and time required to process genomic data by creating a scaffold that can be relied upon by analysis software. Here, we posit that current efforts to extend the linear reference to a graph-based structure while trying to fulfil both of these purposes concurrently will face a trade-off between comprehensiveness and computational efficiency. In this article, we explore how the reference genome is used and suggest an alternative structure, The Genome Atlas (TGA), to fulfil the bipartite role of the reference genome. The human reference genome has provided the foundation for years of genetic discovery and research, but recently, multiple review papers have highlighted the deficiencies of the linear reference, leading to a growing consensus that a richer reference structure is necessary for continued improvements in the era of widespread whole-genome sequencing (WGS) [1.Ballouz S. et al.Is it time to change the reference genome?.Genome Biol. 2019; 20: 159Crossref PubMed Scopus (51) Google Scholar, 2.Computational Pan-Genomics Consortium Computational pan-genomics: status, promises and challenges.Brief. Bioinform. 2018; 19: 118-135PubMed Google Scholar, 3.Yang X. et al.One reference genome is not enough.Genome Biol. 2019; 20: 104Crossref PubMed Scopus (22) Google Scholar]. The Human Genome Project published its first draft in 2001, which, as the only publicly available human genome at nucleotide resolution, became integrated as the backbone of genome analysis. The emergence of new techniques combined with the falling cost of WGS allowed each new release to resolve gaps, improving its overall stability and accuracy. The increased reliability of the reference genome, coupled with the desire for faster and less computationally intensive genome analysis to keep pace with rapidly evolving sequencing methods, shifted the focus of bioinformatics from genome assembly-based (see Glossary) to genome alignment-based algorithms. The availability of a persistent frame of reference enabled genetic variants and their functional roles to be cataloged using the base pair (bp) coordinate system in large databases, such as gnomAD [4.Karczewski K.J. et al.The mutational constraint spectrum quantified from variation in 141,456 humans.Nature. 2020; 581: 434-443Crossref PubMed Scopus (1713) Google Scholar] and ClinVar [5.Landrum M.J. et al.ClinVar: improving access to variant interpretations and supporting evidence.Nucleic Acids Res. 2018; 46: D1062-D1067Crossref PubMed Scopus (1028) Google Scholar]. Further improvements in computational methods, combined with long read sequencing, expanded such databases to include structural variations (SVs). Despite advances in computing, the size of WGS data forces the continued use of a reference scaffold against which to perform genome assembly, read alignment, and variant calling [6.Pereira R. et al.Bioinformatics and computational tools for next-generation sequencing analysis in clinical genetics.J. Clin. Med. Res. 2020; 9: 132Google Scholar,7.Lightbody G. et al.Review of applications of high-throughput sequencing in personalized medicine: barriers and facilitators of future progress in research and clinical application.Brief. Bioinform. 2019; 20: 1795-1811Crossref PubMed Scopus (39) Google Scholar]. Consequently, the reference genome has matured into an integral component of genome analysis pipelines. Over the past decade, the sharp decrease in sequencing costs has led to a flood of large-scale and population-specific sequencing projects, revealing unique sequences missing from the current reference genome [8.Sherman R.M. et al.Assembly of a pan-genome from deep sequencing of 910 humans of African descent.Nat. Genet. 2019; 51: 30-35Crossref PubMed Scopus (115) Google Scholar,9.Lee Y.-G. et al.Insertion variants missing in the human reference genome are widespread among human populations.BMC Biol. 2020; 18: 167Crossref PubMed Scopus (2) Google Scholar], population-specific differences in common genetic variants [10.Maretty L. et al.Sequencing and de novo assembly of 150 genomes from Denmark as a population reference.Nature. 2017; 548: 87-91Crossref PubMed Scopus (65) Google Scholar,11.Wong K.H.Y. et al.De novo human genome assemblies reveal spectrum of alternative haplotypes in diverse populations.Nat. Commun. 2018; 9: 3040Crossref PubMed Scopus (35) Google Scholar], and increasing recognition that some populations, particularly Indigenous, are at risk of being left behind [12.Garrison N.A. et al.Genomic research through an Indigenous lens: understanding the expectations.Annu. Rev. Genomics Hum. Genet. 2019; 20: 495-517Crossref PubMed Scopus (54) Google Scholar]. The GRCh38 reference assembly has been expanded to include alternate loci (ALT loci), most of which are similar to the primary assembly but contain many small variants that commonly occur together. Consequently, naively aligning to a concatenation of the reference genome and ALT loci results in reads that map to multiple ambiguous locations [13.Church D.M. et al.Extending reference assembly models.Genome Biol. 2015; 16: 13Crossref PubMed Scopus (88) Google Scholar]. Additionally, this approach is not an accurate representation of the underlying biological knowledge. It ignores the location information of ALT loci, their prevalence in different populations, and the link between reference sequences and ALT loci. This jury-rigged solution creates hurdles for traditional mapping algorithms [14.Church D.M. et al.Modernizing reference genome assemblies.PLoS Biol. 2011; 9e1001091Crossref PubMed Scopus (244) Google Scholar] and, while these can be partially overcome with newer methods [15.Li H. Aligning sequence reads, clone sequences and assembly contigs with BWA-MEM.arXiv. 2013; 2013 (arXiv:1303.3997v2)Google Scholar], these methods often provide only minimal gains at the cost of space and time complexity. To better reflect the structure of the genome, focus has shifted to the concept of a pan-genome to capture the complete set of variations present in a species. For an extensive review of human pan-genomics, see [16.Sherman R.M. Salzberg S.L. Pan-genomics in the human genome era.Nat. Rev. Genet. 2020; 21: 243-254Crossref PubMed Scopus (54) Google Scholar]. Sherman and Salzberg noted that, while the concept of a pan-genome is appealing due to its completeness and may be suitable for storing the wealth of genetic information, pan-genomes are likely to be cumbersome for fast genome analysis. This suggests that it is more suitable to use methods that create tailored low-complexity reference genomes for individual analyses [17.Grytten I. et al.Assessing graph-based read mappers against a baseline approach highlights strengths and weaknesses of current methods.BMC Genomics. 2020; 21: 282Crossref PubMed Scopus (7) Google Scholar]. Mathematical graphs are a natural choice to describe pan-genomes. Common graph structures, such as de Bruijn graphs (DBGs) and directed, acyclic graphs (DAGs), have formed the basis of assembly tools, including the Gig Assembler used for the first draft of the reference [18.Kent W.J. Haussler D. Assembly of the working draft of the human genome with GigAssembler.Genome Res. 2001; 11: 1541-1548Crossref PubMed Scopus (93) Google Scholar]. The plethora of computationally efficient assembly algorithms optimized for node-labeled graphs made DBGs and DAGs an obvious starting point in the design of richer reference structures. Graphs built from the reference sequence and commonly observed alleles have been used successfully to improve the accuracy of read alignment in complex regions [19.Dilthey A. et al.Improved genome inference in the MHC using a population reference graph.Nat. Genet. 2015; 47: 682-688Crossref PubMed Scopus (107) Google Scholar] and genotyping algorithms [20.Kim D. et al.Graph-based genome alignment and genotyping with HISAT2 and HISAT-genotype.Nat. Biotechnol. 2019; 37: 907-915Crossref PubMed Scopus (1271) Google Scholar, 21.Eggertsson H.P. et al.GraphTyper2 enables population-scale genotyping of structural variation using pangenome graphs.Nat. Commun. 2019; 10: 5402Crossref PubMed Scopus (28) Google Scholar, 22.Valenzuela D. et al.Towards pan-genome read alignment to improve variation calling.BMC Genomics. 2018; 19: 87Crossref PubMed Scopus (19) Google Scholar, 23.Sibbesen J.A. The Danish Pan-Genome Consortium, Maretty L, Krogh A. Accurate genotyping across variant classes and lengths using variant graphs.Nat. Genet. 2018; 50: 1054-1059Crossref PubMed Scopus (21) Google Scholar]. A global effort, predominantly led by the Global Alliance for Genomics and Health (GA4GH), has produced a variety of graph-based structures suitable for genome analysis [24.Kehr B. et al.Genome alignment with graph data structures: a comparison.BMC Bioinformatics. 2014; 15: 99Crossref PubMed Scopus (22) Google Scholar]. However, many face common stumbling blocks: the need for prohibitively high computer resources or poor performance compared with existing linear methods. The task of producing a computationally tractable reference genome that outperforms state-of-the-art linear methods, while still allowing clear visualization using a human-readable coordinate system, appears to be insurmountable. In this article, we encourage the reader to take a step back from the technical task of adapting existing approaches and consider the role of the reference genome in the oncoming era of pan-genomics. We propose a novel way of dividing the tasks a reference genome performs to circumvent the trade-off between computational requirements and analysis performance. We introduce the concepts of The Genome Atlas (TGA), and blueprint genomes, which we believe will work seamlessly together to satisfy the design brief of a new reference genome. Most current approaches to graph-based references focus on extending the linear reference into a graph where nodes are labeled with contiguous sequences and haplotypes are embedded as paths. One of the more impactful reviews of graph-based reference genomes, written by members of the GA4GH Reference Variation Task Group, coined the term 'genome graph' as an encompassing definition for a graph representation of a pan-genome [25.Paten B. et al.Genome graphs and the evolution of genome inference.Genome Res. 2017; 27: 665-676Crossref PubMed Scopus (125) Google Scholar]. They identified four main types of genome graphs: DBG, DAGs, Sequence Graphs, and Bi-Edged Sequence Graphs (Figure 1). Here, we outline many of the common problems facing genome graphs: how to uniquely define positions and instances in a nonlinear graph structure; how to utilize information such as allele frequency and linkage measures; and how to embed haplotypes, while maintaining desirable properties, such as readability [2.Computational Pan-Genomics Consortium Computational pan-genomics: status, promises and challenges.Brief. Bioinform. 2018; 19: 118-135PubMed Google Scholar]. The Variation graph (Vg) [26.Garrison E. et al.Variation graph toolkit improves read mapping by representing genetic variation in the reference.Nat. Biotechnol. 2018; 36: 875-879Crossref PubMed Scopus (150) Google Scholar] was highlighted by Paten et al. as a promising implementation of a genome graph. Vg has been shown to improve read mapping when reads contain variants, but crucially, it performs worse than the widely used BWA-MEM linear mapper on nonvariant reads. This means that, for a full read set, despite performing better than other graph methods [20.Kim D. et al.Graph-based genome alignment and genotyping with HISAT2 and HISAT-genotype.Nat. Biotechnol. 2019; 37: 907-915Crossref PubMed Scopus (1271) Google Scholar,27.Rakocevic G. et al.Fast and accurate genomic analyses using genome graphs.Nat. Genet. 2019; 51: 354-362Crossref PubMed Scopus (73) Google Scholar], Vg achieves lower accuracy compared with linear-based mappers [17.Grytten I. et al.Assessing graph-based read mappers against a baseline approach highlights strengths and weaknesses of current methods.BMC Genomics. 2020; 21: 282Crossref PubMed Scopus (7) Google Scholar]. Improvements gained from the inclusion of known variants become a hindrance when the percentage of known single nucleotide variants (SNVs) increases or insertions and deletions (indels) are included due to the increased number of equally plausible alignment positions [28.Pritt J. et al.FORGe: prioritizing variants for graph genomes.Genome Biol. 2018; 19: 220Crossref PubMed Scopus (24) Google Scholar]. This strongly suggests that using a full pan-genome reference graph during read mapping is not the best approach. A further challenge facing current graph-based representations arises when replicating the index structure of a linear, text-based reference that is central to computational efficiency. The primary goal of indexing is to enable fast string matching within a large set of sequences. Many compressed data structures have been developed and optimized for such tasks focusing on adapting existing k-mer-based approaches (reviewed in [29.Marchet C. et al.Data structures based on k-mers for querying large collections of sequencing datasets.bioRxiv. 2019; (Published online December 6, 2019. https://doi.org/10.1101/866756)Google Scholar]), by extracting k-mers from paths in the graph [16.Sherman R.M. Salzberg S.L. Pan-genomics in the human genome era.Nat. Rev. Genet. 2020; 21: 243-254Crossref PubMed Scopus (54) Google Scholar,21.Eggertsson H.P. et al.GraphTyper2 enables population-scale genotyping of structural variation using pangenome graphs.Nat. Commun. 2019; 10: 5402Crossref PubMed Scopus (28) Google Scholar,30.Lipman D.J. Altschul S.F. Kececioglu JD. A tool for multiple sequence alignment.Proc. Natl. Acad. Sci. U. S. A. 1989; 86: 4412-4415Crossref PubMed Scopus (359) Google Scholar]. Unfortunately, these methods explode in size and construction time for large collections of genomes [31.Shimbel A. Structural parameters of communication networks.Bull. Math. Biol. 1953; 15: 501-507Google Scholar], forcing compromises to be made, such as ignoring small variants [32.Li H. et al.The design and construction of reference pangenome graphs with minigraph.Genome Biol. 2020; 21: 265Crossref PubMed Scopus (28) Google Scholar]. Until a single graph format is adopted as standard, the task of developing efficient computational methods is secondary to the primary focus of designing the best representation for the reference. However, once a single format is widely adopted, then innovative approaches for optimization will emerge, as has been the case for the linear reference. The tremendous effort already invested in the current linear reference has created a level of inertia towards foundational changes. While an overhaul from the status quo may be perceived as daunting, the approach of incremental change is short sighted. A Genome Graph reference that is only a mild extension to the current linear reference is unlikely to harness the full potential of graph-based representations and will ultimately lead to the need for further redesigns in the future. To ensure the longevity of a new reference, we need to consider its role in the era of pan-genomics and determine a full design brief to enable the production of a suitable framework with the capability to adapt as technology advances while minimizing disruption for scientists and clinicians who depend upon it. The current reference genome has two primary, but very distinct roles: to enable knowledge exchange and to increase the efficiency and accuracy of analysis algorithms (Figure 2). To highlight the distinction between the two applications, we can use the analogy of a library containing our collective wealth of knowledge regarding features of the genome, represented by the individual books (Box 1). The analogy exemplifies the difference between microtasks, which need to be performed efficiently (e.g., searching for a word in a book or collection of books), and interpretative tasks (e.g., searching for connections between books using meta-information, such as author or genre). These tasks can be undertaken within the broader framework of the library, but require distinct approaches.Box 1The library analogyConsider a library of books that represents our collective knowledge regarding the human genome. In place of the words of a story, a book contains a DNA sequence that is an allele or a stretch of immutable nucleotides (i.e., one where this is no variation). The complete diploid genome of an individual can be constructed through two series of books read in a specific order. Each book has two distinct sets of information associated with it: (i) meta-information (e.g., author, publisher, date published, or genre); and (ii) the text within. The separation between meta-information and content also appears in genomic features, with meta-information, such as allelic frequency, disease associations, and functional roles.This divide between meta-information and the text within extends to the types of activity we perform in the library, from searching for similarities between books, or looking within a specific book for individual words. Within a library, we are able to search and filter for books using all of the meta-information (e.g., identifying other books by the same author or that were published within the same year, which is the equivalent of finding all alleles associated with a disease or that have a high prevalence in a population). These searches are done without using the words inside the book, because this would make the search excruciatingly slow. The utilization of meta-information for feature identification and identifying the functional implications of the differences between genomes or populations of genomes form the basis of what we have termed 'knowledge exchange' with regards to a reference genome.In stark contrast to knowledge exchange, performing text-based searches of the words within a book or set of books does not require any meta-information (i.e., we do not need to know the name of the author to know whether the word 'cat' appears in a book). The meta-information also does not specify whether two books contain the same words. Performing searches and comparative actions across every book in the library would be computationally intractable, especially if the library is continually expanding. Therefore, the way we organize the content is critical for our capacity to get the greatest benefit from the library. Consider a library of books that represents our collective knowledge regarding the human genome. In place of the words of a story, a book contains a DNA sequence that is an allele or a stretch of immutable nucleotides (i.e., one where this is no variation). The complete diploid genome of an individual can be constructed through two series of books read in a specific order. Each book has two distinct sets of information associated with it: (i) meta-information (e.g., author, publisher, date published, or genre); and (ii) the text within. The separation between meta-information and content also appears in genomic features, with meta-information, such as allelic frequency, disease associations, and functional roles. This divide between meta-information and the text within extends to the types of activity we perform in the library, from searching for similarities between books, or looking within a specific book for individual words. Within a library, we are able to search and filter for books using all of the meta-information (e.g., identifying other books by the same author or that were published within the same year, which is the equivalent of finding all alleles associated with a disease or that have a high prevalence in a population). These searches are done without using the words inside the book, because this would make the search excruciatingly slow. The utilization of meta-information for feature identification and identifying the functional implications of the differences between genomes or populations of genomes form the basis of what we have termed 'knowledge exchange' with regards to a reference genome. In stark contrast to knowledge exchange, performing text-based searches of the words within a book or set of books does not require any meta-information (i.e., we do not need to know the name of the author to know whether the word 'cat' appears in a book). The meta-information also does not specify whether two books contain the same words. Performing searches and comparative actions across every book in the library would be computationally intractable, especially if the library is continually expanding. Therefore, the way we organize the content is critical for our capacity to get the greatest benefit from the library. Considering further how the reference aids computation, we observe that there are two microtasks that are the fundamental building blocks of both assembly and alignment algorithms: (i) global sequence search (e.g., whether a book contains a specific word); and (ii) local sequence comparison (e.g., whether sentences from two different books match). A global sequence search is looking for occurrences of a single, short sequence within a much larger sequence, or set of sequences. Local sequence comparison occurs between two sequences of roughly equal size, for example between k-mers within a reference genome and a read, or between two reads during assembly. The speed at which these tasks can be performed and the number of times they must be performed are the dominant factors when determining the efficiency of genome analysis algorithms. Global sequence searches rely upon compressed data structures to achieve sublinear search times. The most suitable data sets to use as an index are those that change infrequently, but that will be queried a large number of times. Traditionally, an index of the genome is generated, and short sequences extracted from read sets are searched for in the index to act as 'seeds' during read alignment. As noted earlier, there are diverse indexing approaches to enable efficient global searches, with some specifically optimized to minimize index size [33.Kryukov K. et al.Nucleotide Archival Format (NAF) enables efficient lossless reference-free compression of DNA sequences.Bioinformatics. 2019; 35: 3826-3828Crossref PubMed Scopus (11) Google Scholar], or query time [34.Liu Y. et al.Fast detection of maximal exact matches via fixed sampling of query K-mers and Bloom filtering of index K-mers.Bioinformatics. 2019; 35: 4560-4567Crossref PubMed Scopus (7) Google Scholar]. Most methods rely on the reference genome being a static structure (i.e., unchanged during analysis). It is possible to add to an index structure without regenerating it entirely, but for larger sequences it is often quicker to regenerate the index (W. Gerlach, MSc thesis, Bielefeld University, 2007). As we move towards the application of personalized medicine and population-specific sequencing studies, the reference genome is no longer sufficiently static (Box 2), with personalized reference genomes providing more accurate analysis results [17.Grytten I. et al.Assessing graph-based read mappers against a baseline approach highlights strengths and weaknesses of current methods.BMC Genomics. 2020; 21: 282Crossref PubMed Scopus (7) Google Scholar]. The rate at which new genomic information is discovered is dramatically increasing. We posit that, to obtain the most accurate analysis results, we wish to use the newest reference genome and databases available and, given that these are being continually updated, the reference is no longer a suitable data set to construct an index from. Thus, we must reconsider the data sets involved to identify those that are static and more appropriate for indexing.Box 2Indexing a set of readsOn the technical side, there is an additional consideration about the relationship between a reference genome and DNA-sequencing data from an individual. The space and time complexity for creating an indexed reference genome is high. Naturally, the most suitable data sets to use as an index are those that change infrequently, but will be queried a large number of times. In this sense, the data set to be indexed could be labeled as static, because it is rarely, if ever, updated. While the current reference genome is not updated frequently, the databases of genetic variations are growing rapidly due to widespread, low-cost genome sequencing. Therefore, if we define the reference genome as including our knowledge of variants, it can no longer be considered a static data set and, thus, is an unsuitable data set to use for an index.However, if we take a broader look at the data sets involved in WGS analysis, we can see that a read set generated for a genome is unchanged during analysis, with the exception of preprocessing and error correction. Reads are reported sequentially and, thus, it is entirely possible to design an indexing algorithm that is built incrementally in real time as the reads are outputted by the sequencing machine [35.Rathee S. Kashyap A. StreamAligner: a streaming based sequence aligner on Apache Spark.J. Big Data. 2018; 5: 8Crossref Scopus (6) Google Scholar].Recent advances in induced sorting algorithms to generate the suffix array of a string [42.Kärkkäinen J. et al.Engineering external memory induced suffix sorting.in: Fekete S. Ramachandran V. 2017 Proceedings of the Meeting on Algorithm Engineering and Experiments (ALENEX). Society for Industrial and Applied Mathematics, 2017: 98-108Crossref Scopus (21) Google Scholar] have the potential to be optimized specifically for DNA sequences by recognizing that the suffix of a string comprising a concatenated set of reads is unique after the sentinel character that separates the individual reads, provided that duplicate reads are removed through preprocessing. The removal of duplicate reads is not an unreasonable step because increasing read lengths reduces the possibility of identical reads to a negligible level, and allows for the same sentinel character to be used as the divider between all reads. This uniqueness could lead to the development of algorithms that treat each read individually with a small memory footprint, rather than as a concatenated string, which is often too large to fit into computer memory. On the technical side, there is an additional consideration about the relationship between a reference genome and DNA-sequencing data from an individual. The space and time complexity for creating an indexed reference genome is high. Naturally, the most suitable data sets to use as an index are those that change infrequently, but will be queried a large number of times. In this sense, the data set to be indexed could be labeled as static, because it is rarely, if ever, updated. While the current reference genome is not updated frequently, the databases of genetic variations are growing rapidly due to widespread, low-cost genome sequencing. Therefore, if we define the reference genome as including our knowledge of variants, it can no longer be considered a static data set and, thus, is an unsuitable data set to use for an index. However, if we take a broader look at the data sets involved in WGS analysis, we can see that a read set generated for a genome is unchanged during analysis, with the exception of preprocessing and error correction. Reads are reported sequentially and, thus, it is entirely possible to design an indexing algorithm that is built incrementally in real time as the reads are outputted by the sequencing machine [35.Rathee S. Kashyap A. StreamAligner: a streaming based sequence aligner on Apache Spark.J. Big Data. 2018; 5: 8Crossref Scopus (6) Google Scholar]. Recent advances in induced sorting algorithms to generate the suffix array of a string [42.Kärkkäinen J. et al.Engineering external memory induced suffix sorting.in: Fekete S. Ramachandran V. 2017 Proceedings of the Meeting on Algorithm Engineering and Experiments (ALENEX). Society for Industrial and Applied Mathematics, 2017: 98-108Crossref Scopus (21) Google Scholar] have the poten