Genomic fluidity: an integrative view of gene diversity within microbial populations
© Kislyuk et al; licensee BioMed Central Ltd. 2011
Received: 13 October 2010
Accepted: 13 January 2011
Published: 13 January 2011
Skip to main content
© Kislyuk et al; licensee BioMed Central Ltd. 2011
Received: 13 October 2010
Accepted: 13 January 2011
Published: 13 January 2011
The dual concepts of pan and core genomes have been widely adopted as means to assess the distribution of gene families within microbial species and genera. The core genome is the set of genes shared by a group of organisms; the pan genome is the set of all genes seen in any of these organisms. A variety of methods have provided drastically different estimates of the sizes of pan and core genomes from sequenced representatives of the same groups of bacteria.
We use a combination of mathematical, statistical and computational methods to show that current predictions of pan and core genome sizes may have no correspondence to true values. Pan and core genome size estimates are problematic because they depend on the estimation of the occurrence of rare genes and genomes, respectively, which are difficult to estimate precisely because they are rare. Instead, we introduce and evaluate a robust metric - genomic fluidity - to categorize the gene-level similarity among groups of sequenced isolates. Genomic fluidity is a measure of the dissimilarity of genomes evaluated at the gene level.
The genomic fluidity of a population can be estimated accurately given a small number of sequenced genomes. Further, the genomic fluidity of groups of organisms can be compared robustly despite variation in algorithms used to identify genes and their homologs. As such, we recommend that genomic fluidity be used in place of pan and core genome size estimates when assessing gene diversity within genomes of a species or a group of closely related organisms.
The advent of technologies to rapidly sequence entire genomes provides a resource of sequenced genomes spanning the entire tree of life [1–4]. Indeed, as the cost and time to sequence genomes have decreased, it has become possible to sequence multiple individuals from within a species. Re-sequencing efforts have led to the following discovery: the representation of gene families in isolates from the same bacterial species is highly variable [5–9]. This variability poses conceptual as well as applied problems. Conceptually, the variability suggests the need to further re-visit species definitions that rely upon comparisons of highly conserved components of the genome, such as 16S rRNA sequences [10–14]. In addition, horizontal gene transfer and other genome rearrangements such as gene deletions and duplications can radically change the phenotype of a bacterium, even within individuals of the same species . For example, the introduction of toxin genes can render a bacterium pathogenic. Hence, from an applied perspective, there is an increasing need to quantify the gene diversity of a species or genus with pathogenic potential [6, 7, 16–19]. The core and pan genome concepts have been proposed as a way to characterize the distribution of gene families within a group of organisms, e.g., within a species or genus [5, 6, 16, 18, 20–22]. The core genome is the set of genes found in every organism within a group (whether sequenced or not). The pan genome is the set of all genes found within organisms of a group (whether sequenced or not), including core genes and genes which appear in a fraction of genomes. Intuitively, the core genome preserves the notion that genomes of closely related organisms have something in common, while the pan genome is in accord with the finding that gene composition differs even among genomes of closely related organisms. In that sense, the core and pan genome concepts begin to address both conceptual problems (e.g., what is a bacterial species?) and applied problems (e.g., how likely is it that an individual of a given bacterial species is a pathogen?). Multiple attempts have been made to estimate the size of pan and core genomes in hopes of quantifying how open or closed a particular set of genomes is to gene exchange [5, 7, 8, 23, 24]. However, estimating the actual list of genes in the pan and core genomes remains intractable.
Thus far, attempts to quantify the size of the core and pan genomes have been based on extrapolations from a limited number of sequenced strains (usually on the order of a dozen or few dozen genomes) to the entire group (generally unknown, but easily upwards of 1012 genomes). Results of such extrapolations have been widely divergent. In the most well-studied case, the pathogen Streptococcus agalactiae, estimates of the pan genome size vary from tens of thousands  to infinite . Extreme variation in estimates of core and pan genome sizes makes it difficult to utilize these measures to quantify or compare the degree of acquisition and loss of gene families within a particular group or to make meaningful biological interpretations of the core and pan genome concepts. One might suspect that robust quantification of core and pan genomes sizes could be achieved with improved statistical estimation methods, combined with increased sequencing coverage. This is not the case. The problem of estimating pan and core genome sizes will not be resolved by gradual improvements in sequencing.
In this paper we demonstrate that current methods to estimate pan and core genome sizes are statistically ill-posed. We do so by demonstrating that sample gene distributions drawn from artificially generated groups of genomes with radically different pan and core genomes sizes are statistically indistinguishable. In contrast, we present an alternative diversity metric, genomic fluidity, whose expected value is equivalent whether estimated from the sample or from the true gene distribution. We then apply a bioinformatics pipeline so as to estimate genomic fluidity within 7 multiply-sequenced bacterial species containing 109 sequenced genomes. We test the robustness of genomic fluidity to changes in the number of sequenced genomes as well as to changes in alignment parameters. In so doing we demonstrate when it is possible to reliably rank order species in terms of genomic fluidity and discuss the implications of our work for inferring information about gene distributions based on subsamples.
where U k and U l are the number of gene families found only in genomes k and l respectively and M k and M l are the total number of gene families found in k and l respectively. Importantly, the same formula for fluidity applies whether N represents the total number of genomes in the population or N represents the total number of genomes in the sample. In other words, genomic fluidity is an estimate of gene-ic dissimilarity, akin to similarity measures used in the study of ecological communities  (see Additional file 1, Figure S2 for a schematic illustration of Eq. (1)). More specifically, genomic fluidity estimates how dissimilar genomes are when evaluated at a gene level. For example, a genomic fluidity of 0.1 represents that a pair of genomes have on average 10% unique genes and share 90% of their genes. As fluidity increases, so too does the probability that gene content differs between genomes in a sample. Genomic fluidity also provides information on novelty in sequencing projects. To see how, note that the best estimate for the probability that a random gene from a newly sequenced genome is not found in a randomly selected prior sequenced genome is simply φ. Importantly, genomic fluidity is robust to small sample size: it can be reliably estimated from a few sampled genomes. For example, in Figure 2 we show how the genomic fluidities for synthetically generated gene distributions are equivalent whether estimated from the true distribution or from a few dozen sampled genomes. In addition, subtle differences in the genomic fluidity between two species can be detected from a small number of sampled genomes. The estimated variance of fluidity was calculated using the jackknife estimate , which is based on leave-one-out statistics (see Methods for more details). In contrast, rarefaction curves used to estimate pan and core genome sizes are statistically indistinguishable for synthetically generated gene distributions, even when the underlying pan and core genome sizes are radically different (see Figure 1C, D).
These results are generally consistent with previous suggestions that B. anthracis has a closed genome, that N. meningitidis may have an open genome due to its natural competence, and that Strep. agalactiae has an open genome . However, now we can describe a group of organisms as being relatively open or closed, instead of being strictly open or strictly closed. In addition, we can utilize variance estimates to suggest when greater sequencing is needed. The comparison of the rank order of φ between species is consistent with recent calls  to utilize the rank, not the absolute magnitude, when comparing the relative diversity of complex ecological communities. This issue is particularly important in the case of gene diversity studies when identification of gene families is strongly depend on thresholds utilized in bioinformatics pipelines. Note that we do not suggest ranking of pan and core genome size estimates, since common genes and genomes do not, in general, inform estimates of rare genes and genomes, respectively.
The proposal that there exists a core and pan genome for bacterial species represents a significant advance in the conceptualization of gene variability within microorganisms . The basic premise of these two concepts have been borne out by the finding that the gene content of bacteria can vary significantly when comparing the sequence of two isolates from a species or genus [5–9, 18, 19, 22]. For example, it is now well established that some genes are found in most, if not all, sequenced genomes of isolates from within a sample. In addition, it is also well established that some genes are found in very few, if only one, sequenced isolate within a sample. However, as we have demonstrated here, efforts to infer the size of the pan and core genomes of an entire species or genus from the frequency distribution of genes within a small sample of sequenced genomes will almost certainly fail. Similarly, efforts to compare the core or pan genomes sizes of bacterial species or genera will be uninformative. The reason is that pan and core genome sizes depend sensitively on the frequency of rare events (such as a rare gene occurring in a genome) whose frequency cannot be accurately estimated from a small sample of sequenced genomes. Instead, we have proposed the use of an alternative diversity metric - genomic fluidity - which is a reliable and robust estimator of the gene dissimilarity amongst a group of sequenced genomes.
This study has a number of key implications for future sequencing efforts. First, it suggests that efforts to understand a single species by sequencing as many isolates as possible may be limited in their ability to comprehensively define the diversity within that species . Clearly, such studies will remain important in their ability to describe expected genomic differences (in contrast to rare genomic differences). Next, our findings also suggest that the expected gene dissimilarity within a given species can be well characterized by sequencing a relatively small number of well-chosen representatives. Sequencing a few dozen genomes is a fairly straightforward task given recent advances in sequencing technology. Finally, perhaps the most far-reaching implication of the work presented here is that we have shown it is possible to compare the relative genomic fluidity of different groups of bacteria (e.g. species, genera, or higher). We have shown that genomic fluidity can reliably distinguish between subtle differences in true gene distributions (in a computational study) as well as determine when it is possible to rank-order a set of 7 species based on the analysis of 109 whole genomes (in a bioinformatics analysis).
Genomic fluidity necessarily varies with the phylogenetic diversity of the group of genomes under consideration. In many cases, this level of diversity is defined through a species or other group definition via observed phenotypic aspects, and not through any account of genome-level divergence. As a result, within-species gene diversity of bacterial species varies greatly. To facilitate fluidity-based comparisons between species, one possibility is to normalize genomic fluidity by the average or median phylogenetic distance between members of the considered group, such as the phylogenetic distance computed using a multiple alignment of housekeeping genes . However, other normalizations are possible and we consider this to be an important target for future research.
Despite its merits, genomic fluidity is not meant to describe all forms of genome variation. Genomic fluidity can provide a reliable estimate for how many new genes additional sequencing is likely to reveal, with respect to a previously sequenced genome. It cannot, however, provide an estimate of the amount of sequencing necessary to cover the gene novelty in the entire group (for reasons similar to why estimates of the pan genome size are impossible). In addition, genomic fluidity restricts itself to one component of genomic difference. There are a variety of forms of genomic differences beyond gene compositional differences or the more classic finding of single-nucleotide polymorphisms. Genomes may differ in terms of gene synteny , copy number variation [31, 32], plasmid and/or prophage presence , codon biases [33, 34], and methylation state . It would be prudent to consider other diversity metrics, in addition to the metric of genomic fluidity studied here, that account for forms of variation in genome state amongst closely related organisms.
Genomic fluidity is an integrated measure of gene diversity within a group of organisms. Genomic fluidity is both estimable given a small number of sequenced genomes and robust to variation in alignment parameters. As such, we recommend that genomic fluidity be used in place of pan and core genome size estimates when assessing gene diversity within a species or a group of closely related organisms. However, the precise relationship between variation in gene composition and genomic fluidity with underlying mechanisms of gene family diversification are yet to be resolved . Recent calls for comparing and contrasting the average overlap of gene content with respect to average nucleotide divergence provide one possible route to disentangling the effects of ecological and genomic structure , but much work remains at the interface of bioinformatics and ecological analysis. For example, the detailed comparison of complete bacterial genomes from closely related biofilm-forming bacteria revealed how and why different organisms have adapted to and shaped their environment . Similarly, genomic analysis of cyanoviruses sampled in the oceans helped uncover photosynthetic pathways which enable the exploitation of a niche distinct from previously cultured E. coli based phages despite sharing many common genes and genome architecture . Genomic fluidity complements the detailed functional comparison of genomes by robustly estimating dissimilarity of genes within groups of genomes and providing insight into their potential evolvability. In so doing, our results highlight the need for continued focus on developing new toolsets for assessing what can be inferred about the genome composition and diversity of prokaryotic species and communities based on analysis of a sub-sample of genomes.
Complete annotated genomes and draft annotated genomes were retrieved from NCBI GenBank in the GenBank format. Genomes were automatically re-annotated without hand-curation using a recently developed infrastructure resulting in new GenBank-formatted files . Automatic re-annotation removes annotation bias arising from variability in annotation methods, depth of curation, and the resulting impact on the list of candidate genes - a similar approach was recently used in the analysis of genomes within a bacterial genus . Following this process, putative protein sequences were extracted from annotated CDS regions and aligned using BLASTP  in all vs. all pairwise amino acid alignment. A pair of genes were considered homologous if the protein alignment covered more than c fraction of each gene's length and identity in the alignment exceeded i. To improve performance, alignments were parallelized between nodes on a compute cluster using the Torque PBS job scheduler.
Next, genes were clustered into gene families using a strict clique requirement, i.e. each new gene considered for inclusion into a family must have an alignment with every member of the family satisfying the minimum criteria described above. In this implementation, we compare all members of a gene family to each other on an equal basis and do not distinguish between orthologs, homologs, or paralogs. This homology-based approach is appropriate, since the fine resolution and gene family structure afforded by true ortholog reconstruction does not affect the inclusion or exclusion of genes with marginal evidence of homology.
Mean and variance of fluidity for the species and conditions examined here are presented in Additional file 2.
Consider two sets of genomes, the first set consisting of n 1 genomes, the second set consisting of n 2 genomes. For each pair of genomes, we determine the fraction of the total number of unique genes and the total number of genes. Averaging over all pairs in the first set gives the fluidity ; in the second set . Suppose and we want to determine whether this inequality is significant.
From the theory of U -statistics it is known that the estimated fluidity has approximately a normal distribution . The mean of this distribution is estimated to be in the first set and in the second set. The variance is estimated (by jackknifing) to be in the first set and in the second set. We use the parameters of the approximate normal distributions to compute the significance of the observed fluidity differences. Formally, this corresponds to a two-sample two-sided z -test with one degree of freedom (the effective number of degrees of freedom are taken into account by the jackknife estimation).
We thank King Jordan, Jessa Lee, Tim Read and Matt Sullivan and two anonymous reviewers for their feedback and suggestions on the manuscript. We thank Anju Varadarajan for her assistance in the implementation of the bioinformatics pipeline. We thank Lee Katz, Scott Sammons, Dhwani Govil, Brian Harcourt, King Jordan and Leonard Mayer for providing access to N. meningitidis genomes sequenced at the Centers for Disease Control and Prevention and for help with their analysis. We acknowledge the support of the Defense Advanced Research Projects Agency under grants HR0011-05-1-0057 and HR0011-09-1-0055. Joshua S. Weitz, Ph.D., holds a Career Award at the Scientific Interface from the Burroughs Wellcome Fund. Andrey Kislyuk was supported, in part, by a Georgia Tech Research and Innovation travel grant.
This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.