https://orcid.org/0000-0002-7811-6711
Abstract
Bacteria and Archaea display a variety of phenotypic traits and can adapt to diverse ecological niches. However, systematic annotation of prokaryotic phenotypes is lacking. We have therefore developed ProTraits, a resource containing ∼545 000 novel phenotype inferences, spanning 424 traits assigned to 3046 bacterial and archaeal species. These annotations were assigned by a computational pipeline that associates microbes with phenotypes by text-mining the scientific literature and the broader World Wide Web, while also being able to define novel concepts from unstructured text. Moreover, the ProTraits pipeline assigns phenotypes by drawing extensively on comparative genomics, capturing patterns in gene repertoires, codon usage biases, proteome composition and co-occurrence in metagenomes. Notably, we find that gene synteny is highly predictive of many phenotypes, and highlight examples of gene neighborhoods associated with spore-forming ability. A global analysis of trait interrelatedness outlined clusters in the microbial phenotype network, suggesting common genetic underpinnings. Our extended set of phenotype annotations allows detection of 57 088 high confidence gene-trait links, which recover many known associations involving sporulation, flagella, catalase activity, aerobicity, photosynthesis and other traits. Over 99% of the commonly occurring gene families are involved in genetic interactions conditional on at least one phenotype, suggesting that epistasis has a major role in shaping microbial gene content.
INTRODUCTION
Bacteria and Archaea inhabit a wide spectrum of ecological niches, including growth in extreme environments and association to plant or animal hosts, whether mutualistic, commensal or parasitic. This is made possible by a plethora of physiological adaptations observed in prokaryotes, such as the use of different carbon sources and electron acceptors, resistance to stressors and molecular interactions with host cells. Broadly construed, the notion of a microbial phenotypic trait encompasses all of the above facilities – the ability to colonize ecological niches and the underlying physiological features. A deeper characterization of such traits can be obtained by linking them to the genetic makeup of the microbes (1,2). Statistical associations of genes to phenotypes can implicate certain proteins or pathways, providing insight into the mechanistic basis of phenotypic traits (3,4), as demonstrated for adaptation to stress (5,6), host-association (7,8), pathogenesis (9,10), drug resistance (11,12) and relevance to biotechnological applications (13,14).
The amount of prokaryotic genomes is increasing rapidly (15), aided by single-cell sequencing (16) and metagenomics (17). However, efforts to obtain systematic, high-quality phenotype annotation of microbes are not keeping pace, meaning that the potential for comprehensive gene-trait association studies cannot be realized (1,18). A thorough characterization of phenotypes of phylogenetically diverse microbial taxa could implicate genes in phenotypic traits by capturing the evolutionary signal across the prokaryotic tree of life. Moreover, multiple traits could be considered jointly, as they have been for mammalian phenotypes (19,20), thus boosting statistical power and elucidating relationships between the traits (21). Comprehensive resources that link phenotypes to genes and pathways exist for several eukaryotic model organisms (22,23), but not so for prokaryotes.
Current databases of microbial genome sequencing projects (24,25) do supply a certain amount of annotated traits, however, these important resources do not focus on phenotype data and thus the coverage is not extensive. In contrast, the scientific literature abounds with trait descriptions stored as unstructured text, which are therefore not directly accessible to automated analyses. Relying on manual curation to organize these data does not scale with the increasing volume of scientific publications. Motivated by the above, we have developed the ProTraits resource, a comprehensive atlas of 424 microbial phenotypes, covering 3046 Bacterial and Archaeal species, each receiving on average 165 novel high-confidence annotations in addition to 23 previously known labels; available online at http://protraits.irb.hr/.
The pipeline behind ProTraits relies on automated text mining of diverse corpora of biological literature, annotating microbes with existing phenotypic traits, while additionally being able to define novel phenotypic concepts from free text in an unsupervised manner. Furthermore, our approach draws extensively on genomic data, including novel methods of automated phenotype inference from conserved gene neighborhoods and from evolution of synonymous codon biases. The broad coverage with phenotype annotations allows us to systematically discover thousands of high-confidence statistical associations that link genes to traits, thereby informing about their genetic basis. Finally, our analyses suggest that epistatic interactions are ubiquitous within bacterial and archaeal gene repertoires, affecting almost all commonly-occurring gene families.
MATERIALS AND METHODS
Definitions
Here, we adopt a broad definition of a phenotypic trait (or, simply, trait) of a microorganism. Firstly, this encompasses the common meanings of this term in the field of microbiology (metabolic capabilities, morphology, growth conditions). Secondly, this also encompasses the microbe's ability to colonize certain ecological niches, which includes the association to a particular host and an organ/tissue thereof, and moreover the type of this association (pathogenic, symbiotic). While a trait is any attribute of an organism which matches the above definition, the term ‘phenotype’ also implies a value of this attribute. Upon initial data collection (described below), all traits were binarized, meaning that multi-valued categorical traits were split into three or more binary traits. In other words, as defined herein, there are two phenotypes each trait can exhibit: presence or absence of the trait.
Phenotype data collection
The initial set of phenotype assignments was combined from several sources. First, we merged the phenotype data in the NCBI microbial genome projects list (‘lproks0’ table, now retired) with the largely overlapping set of traits from the BacMap database (24). We collapsed together synonymous traits from the two databases and furthermore the data on causal roles in diseases for various hosts from BacMap was manually categorized (Supplementary Table S1). In rare cases of traits exhibiting discordant phenotype labels between databases, we recorded that trait/species combination as a missing value (n = 147, list in Supplementary Table S1).
Second, we included descriptions of the ecosystems where microbes were isolated from, as provided in the GOLD database of genome sequencing projects (25). Since GOLD provides only positive assertions, we have provisionally annotated as negative examples for a certain ecosystem all those organisms which were not explicitly assigned to that ecosystem type. For instance, organisms annotated as ‘marine’ were used as provisional negatives for ‘soil’ or for ‘thermal springs’. In GOLD, each organism can have more than one assigned value, e.g. being annotated as both ‘marine’ and ‘freshwater’ and thus receiving positive labels for these two ecosystems, and negative labels for all remaining ones. Such provisional negative annotations were used to train the classification models (see below).
Third, we further considered a set of biochemical phenotypes, which here implies the experimental data we manually collected from 265 articles describing new microbial strains, published in the IJSEM journal (26) during the years 2013 and 2014. These publications describe series of laboratory tests meant to discriminate novel species, such as the ability to grow on a particular substrate. Next, we collapsed together synonymous biochemical traits, further retaining the trait only if at least 30 species were covered with annotations (full list of synonymous names in Supplementary Table S1). We merged this data with a smaller set of experimental measurements of growth on various substrates for 40 species (27), measured using Biolog phenotype arrays, a technique with broad application to metabolic modeling (28).
Text sources and preprocessing
We analyzed text documents describing 1640 bacterial and archaeal species from (i) Wikipedia, (ii) MicrobeWiki, (iii) HAMAP proteomes, (iv) PubMed abstracts retrieved upon searching for the name of each species (setting ‘sort by relevance’; first 100 items kept); (v) a set of PubMed Central publications that were looked up via the ‘Reference’ field in the KEGG Organisms repository; and (vi) a mixed collection of smaller resources that were pooled together (Bacmap, Genoscope, Joint Genome Institute, KEGG, NCBI and Karyn's Genomes). From the HAMAP source we downloaded all available organisms, for Wikipedia and MicrobeWiki we searched for bacterial and archaeal species names as defined in NCBI Taxonomy, and for the ‘Pubmed abstracts’ corpus (which was retrieved last) we searched for all species for which we already had at least one text document. All documents describing bacterial/archaeal strains were mapped to the corresponding species-level taxon, as defined by the NCBI Taxonomy. The documents were pre-processed by removing reference parts, English language stop-words (29) and all words occurring less than four times. Porter stemming (30) was used to reduce words to their root form.
Unsupervised discovery of phenotypic concepts
We applied non-negative matrix factorization (NMF) (31), a statistical method commonly used for text analysis (32,33) in order to discover novel phenotypic concepts which may not be represented in existing databases. To this end, we constructed a standard ‘bag-of-words’ representation: a matrix where rows correspond to words and columns to organisms for each of the five corpora separately (excluding the mixed collection). In order to enforce consistency between corpora, for the NMF analysis we used only the words that appeared in all corpora (see Supplementary Methods for an exception). Matrix element values were the term frequency-inverse document frequency (tf-idf) weights, often used as a measure of how important a word is to a document, given a collection of documents (34). Briefly, tf-idf normalizes data so as to give greater weights to words which occur rarely, and thus are presumably more informative, while common words are down-weighted. In the NMF analysis, we filtered out names of bacterial and archaeal taxa, as well as a custom list of very frequent words (Supplementary Table S1).
We performed NMF in Matlab 2011b using alternating non-negativity constrained least squares (35) on each data matrix separately, with the number of hidden factors set to 50 and to 100 (in two separate NMF runs). We grouped together similar NMF factors across corpora, while requiring that a similar phenotypic concept had to be consistently discoverable in at least three corpora. The similarity between factors was calculated as the Pearson correlation coefficient between the NMF weights of top 20 ranked words. We further summarized the group by finding the median of word weights across factors in the group, so that each group (phenotypic concept) was represented with the top ranked words of the summary. These groups were then manually examined and those that were of high quality, in terms of consistency and content of words relevant for phenotypes, were further retained and treated in the same way as the other phenotypic traits collected from databases (see above). The groups and the constituent NMF factors are listed in Supplementary Table S2. Several runs of the NMF algorithm were performed to maximize coverage of discovered concepts (Supplementary Methods).
Predicting phenotypes by text mining
The training and the unlabeled data sets were generated for each trait and for each of the six text corpora separately, wherein learning examples were species, each described with the document(s) assigned to that species in a given text corpus. Again, we a used a standard bag-of-words representation with tf-idf weighting of word frequencies across texts; unlike the NMF analysis, here all words were used. All species with known positive or negative trait labels were part of the training set, while all species without known annotations for that trait were in the unlabeled data set, which later received predictions. A support vector machine (SVM) classifier (36) with a linear kernel was trained for each trait and for each text corpus separately using LibSVM (37). The regularization parameter C was optimized in five runs of 4-fold cross-validation (C = 2−15,2−14,…,25), retaining the value yielding the best average area-under-the-curve (AUC) score. We then used a single SVM run of 10-fold cross-validation with the optimal C to estimate confidence scores of predictions, which were converted to precision (or positive predictive value) scores; equivalent to 1-false discovery rate (FDR). The precision scores were obtained separately for the positive and the negative class by thresholding the precision-recall curves obtained via cross-validation; see Supplementary Methods. If a training set corresponding to a trait-corpus pair did not have AUC > 0.6, or did not have at least one example that could be assigned at FDR < 0.25, we considered it as ‘unlearnable’ and did not use those predictions. We retained 424 of 522 traits that were learnable from at least one text corpus. In case of the GOLD ecosystems, provisionally annotated negatives (see above) were used for SVM training, but they were also included in the unlabeled data set to which the trained models were applied, meaning these provisional annotations could potentially be reassigned.
Inferring phenotypes from genomes and metagenomes
Prokaryotic genome sequences and gene annotations were from NCBI Entrez Genomes, while COG/NOG gene families were from eggNOG 3 (38). We only considered species with available genomes having a quality score ≥ 0.9 in (39). We further complemented inferences from text mining using genomic data, by predicting traits from: (i) the proteome composition, encoded as relative frequencies of amino acids (40,41); (ii) the gene repertoire, encoded as presence/absence indicators of COG gene families in a genome (42,43); (iii) co-occurrence of species across environmental sequencing data sets (44); (iv) gene neighborhoods (45,46), encoded as pairwise chromosomal distances between commonly occurring COGs; and (v) genomic signatures of translation efficiency in gene families (6,47), encoded as the MILC codon bias measure (48). For all five genomic representations, learning examples were species and class labels were the phenotypic traits. For the details of how the genomic features were calculated, see Supplementary Methods. For each representation and each trait, the Random Forest (RF) classification model (49) was trained to learn to distinguish between positive and negative trait annotations, using the FastRandomForest (50) implementation with 500 trees. A single run of 10-fold cross-validation was used to determine the performance and to estimate the precision/FDR scores of individual predictions in the same manner as for the text mining (Supplementary Methods). Again, we required cross-validation AUC > 0.6 and at least one organism labelled at FDR < 0.25 to consider a trait-genomic representation pair as learnable. For the purposes of determining individual genomic features informative of traits, we used the RF Gini importance measure (as implemented in scikit-learn; 200 trees; missing values imputed by the median).
Gene-trait associations and epistatic interactions
The Supplementary Methods describe the logistic regression-based statistical methodology we used to detect associations between the occurrence of 80 576 prokaryotic COG/NOG gene families and each phenotype, as well as the genetic interactions between pairs of COGs conditional on a phenotype, while controlling for phylogenetic relatedness of organisms.
Gene function analyses
The lists of genes known to be involved in endospore formation were from the B. subtilis genetic screen in (51) and from the collection of differentially expressed genes during sporulation (52). For the analysis of conserved gene neighborhoods associated with sporulation, we considered only organisms belonging to the Firmicutes and Actinobacteria, as almost all known spore-forming bacteria are contained within. There is a known confounding effect of genome size with the other genetic determinants of spore-forming ability (52). We have therefore used the propensity score matching methodology (53) to control for this effect, as implemented in the R package Matching, with parameters replace = false and caliper = 0.5. In brief, this algorithm resamples genomes to create two sets, one sporulating and other non-sporulating, which are matched by their genome size distribution. For every pair of examined COGs, we repeated this matching only on the set of genomes containing both COGs. The pairwise distances of the COGs in the sporulating versus the genome size-matched set of non-sporulating bacteria were compared using Mann–Whitney tests on all 44 850 pairs of 300 ubiquitous COGs (see above).
RESULTS
Integrating the known phenotype annotation databases
The ProTraits methodology infers microbial phenotypes from free-text and from genomic data by using supervised machine learning algorithms, in particular the Support Vector Machine (SVM) and Random Forest classifiers (Materials and Methods). In order to operate, such classifiers require an initial set of labelled examples – here, organisms with known phenotypes. These were obtained from four sources (Figure 1A): (i) 97 microbial traits collected from NCBI Genomes and from the BacMap genome atlas (24); (ii) 105 ecosystems where microbes were isolated from, provided by the GOLD database of genome sequencing projects (25); (iii) 113 phenotypic concepts that we inferred de novo by a topic modeling of biological literature, see below; and (iv) 109 biochemical phenotypic traits that we manually curated from experimental measurements in scientific publications (Materials and Methods), including phenotype array data from reference (27). Considered together, these sources span a wide range of microbial traits, however, they are often assigned only to a limited number of microbes (Figure 1B). In particular, out of 3046 prokaryotic species with sufficient coverage with text or high-quality genome sequences (Materials and Methods), a trait is currently annotated to a median of 112 species (Q1–Q3: 47–327), counting both the positive and the negative assignments. The modest coverage presents an opportunity for automated inference methods that draw on unstructured text sources, such as the ProTraits pipeline.
Inferring microbial phenotypic traits from free text and genomic data. (A) A schematic overview of the ProTraits pipeline. (B) The coverage with currently available trait labels in various databases. Inlay shows number of bacterial and archaeal species with available free-text and high-quality genome data we used for phenotype inference.
Defining novel phenotypic concepts from free text
We collected text documents describing microbes from six sources: (i) Wikipedia, (ii) MicrobeWiki, (ii) HAMAP proteomes (54), (iv) PubMed abstracts, (v) PubMedCentral publications and (vi) a collection of smaller text sources. All six corpora were encoded using a standard ‘bag-of-words’ approach, wherein each microbe was represented by the normalized word frequencies in the texts describing the organism (Materials and Methods). In order to discover novel phenotypic traits from this extensive collection, we applied non-negative matrix factorization algorithm (NMF (31)) to model phenotypic concepts across the texts (55). Each individual corpus was processed with repeated NMF runs and the results were clustered to uncover trends that are consistent across at least three text corpora (Figure 2A; Materials and Methods). Such novel phenotypic trait candidates (example in Figure 2B; list in Supplementary Table S2) further passed manual curation to ensure that the clustered words conform to a common theme, thus resulting in 113 non-redundant traits discovered from biological literature de novo. Expectedly, this method has recovered some of the phenotype annotations in existing databases: based on the overlap in microbes assigned to the known versus NMF-inferred traits, we find 9 traits highly similar to known ones (Spearman correlation of species’ NMF weights >0.7), serving to validate the methodology. We have compared this NMF-based approach against a previous methodology that used principal components analysis (PCA) followed by hierarchical clustering (HC) to discover ‘word sets’ informative of various traits (56). We performed a blinded evaluation, where two human curators consistently judged the NMF concepts to be preferable over the PCA+HC word sets, according to subjective criteria (sign-test P-value ranging from 10-8 to 10-5 across evaluators and text corpora; outputs of the methods and the evaluation results are in Supplementary Table S2).
Unsupervised discovery of phenotypic concepts from text by non-negative matrix factorization (NMF). (A) An example phenotype ‘beer/spoilage/spoil/lactic/beverage’ captures lactic acid bacteria that cause beer spoilage. NMF weights for selected species (listed left) are shown separately for various keywords from three text sources: Wikipedia, MicrobeWiki and PubMed abstracts. Gray squares are missing data. (B and C) The final NMF-inferred phenotype consists of (B) keywords and of (C) species that were consistently highly weighted by NMF across the text corpora; shown weights are from medoids of clusters of NMF factors.
Automated annotation of microbes with phenotypes
In order to annotate bacterial and archaeal taxa with new phenotypes, we have input the bag-of-words data sets and the known phenotype labels to a SVM classifier. This was done separately for each trait and for each text corpus, and cross-validation was used to estimate the accuracy for the 2272 trained SVM models. Of note, the ecosystem classification of GOLD does not provide explicit negative annotations (Figure 1B) and we thus modeled using provisional negative labels (Materials and Methods). We also considered four other algorithms for this text classification task, but the SVM consistently outperformed them in accuracy on held-out data, across many different traits (Supplementary Figure S1; Supplementary Table S3).
In total, 424 (of 522) phenotypic traits could produce one or more highly-confident novel predictions (Methods), and this ‘learnable’ set of traits was considered further. In many cases, the text-based SVM models were highly predictive: median AUC scores for the traits in NCBI+Bacmap, GOLD and NMF phenotypic concepts of 0.87, 0.86 and 0.93, respectively (Figure 3C). Expectedly, the novel concepts appear easiest to predict, since they were discovered from the same text corpora using a different method (NMF). Biochemical phenotypes were challenging to infer from text in an automated manner (Figure 3C), but there were individual examples with good accuracy.
Annotating microbes with phenotypic traits from free text. (A) An example precision-recall curve for predicting halophilicity from text of MicrobeWiki pages using a Support Vector Machine (SVM) classifier (left). Arrows highlight microbes recovered at certain threshold values of precision (equivalent to 1-false discovery rate, FDR). The proportion of recovered halophilic microbes (recall) increases with more permissive precision thresholds (shown right). (B) The keywords important for predicting halophilicity in the six text corpora, as estimated by the SVM. (C) Accuracy of 2272 SVM models trained on all combinations of text corpora (columns) and phenotypic trait groups (rows), measured as the area-under-curve (AUC) score, where 1.0 is perfect performance and 0.5 indicates random guessing. Plots are histograms, showing the total number of classification models (traits) in a certain accuracy range. Data shown in all panels are from cross-validation.
Next, we applied the trained SVM models to predict phenotypes for microbes that had no previous label for the given trait, while obtaining FDR estimates for each prediction using crossvalidation precision-recall curves (57). An example for the halophilic phenotype is shown in Figure 3A, where progressively more relaxed FDR thresholds allow better coverage with phenotypic annotations, while trading off precision. Importantly, even the less-accurate SVM models can still annotate some microbes confidently (at low FDR), although the coverage is limited in such cases (Supplementary Figure S2A and B). The AUC scores for the individual classification models and the estimated coverage with annotations at several FDR thresholds is provided in Supplementary Table S3.
Orthogonal predictions from comparative genomics
We next complement these text-mining inferences by using an orthogonal approach, which relies on genome sequences and encompasses five independent data sources. These include three data types known to be informative of microbial phenotypes, in particular (i) the amino acid content of the proteome, e.g. (41), (ii) the gene repertoire of the genome, e.g. (42), here encoded as the presence/absence profile of COG/NOG gene families; and (iii) the patterns of co-occurrence of prokaryotic taxa across environmental sequencing data sets (44). Such classification models have previously been applied to predict small sets of up to 10 phenotypic traits, e.g. (40,43). Here, we extend this approach to a comprehensive set of 424 traits, thereby creating many broadly accurate models (median cross-validation AUC 0.83, 0.85 and 0.69, for the proteome composition, gene repertoires and metagenome co-occurrence, respectively; Figure 4A). While the biochemical phenotypes were, overall, challenging to infer by text mining (see above; Figure 3C), the genome data is more informative in this respect, thus complementing the text-based annotations (median AUC from text is 0.54 versus 0.72 from genomes; Figure 4A). For this task of phenotype inference from comparative genomics, we used a Random Forest classifier, which we found to outperform the SVM and three other algorithms, providing more accurate predictions on a held-out data set (Supplementary Figure S1; Supplementary Table S3).
Inferring microbial phenotypes from comparative genomics. (A) Accuracy of RF classification models in predicting groups of phenotypic traits (rows) from five independent representations of genomic data (columns). This includes known methods: proteome composition (blue), gene repertoires (green) and metagenomic co-occurrence (red), as well as two novel methods: gene neighborhoods (yellow) and translation efficiency (estimated via codon usage biases; orange). Histograms show frequency of RF models of certain accuracy, estimated as the cross-validation AUC score, ranging from 0.5 (random guessing) to 1.0 (perfect performance). (B) Examples of gene synteny patterns associated with the spore-forming phenotype, which involve ribosomal proteins or translation factors (columns) and known sporulation factors (rows). Overlaid percentages are FDRs for significant shift in distribution of gene distances in sporulating microbes by Mann–Whitney test. (C) Gene neighborhoods involving a ribosomal protein and a sporulation gene are more often associated with spore-forming ability than other pairs (P < 0.001 for the first distribution of P-values versus each other distribution, Mann–Whitney test). Dashed line is the FDR = 1% threshold. (D) Precision-recall curves, in cross-validation, for two example phenotypes where the gene neighborhoods (top) and codon usage biases (bottom) have prediction accuracy that compares favorably to the gene repertoire method.
Gene synteny patterns accurately predict phenotypes
Evolutionarily conserved gene neighborhoods in prokaryotes (45) are known to reflect gene functional interactions (58,59). Given that a phenotype emerges from a network of many such interactions, we hypothesized that presence of certain gene neighborhoods in a genome may be used to infer the phenotype. One known example are the syntenic regions associated with cold tolerance in the genus Shewanella (46). Here, we systematically look for links between gene neighborhoods and various traits. We have thus described each genome using the pairwise chromosomal distances of 300 COGs that commonly occur across organisms (Materials and Methods), and indeed found that this can predict phenotypes accurately (median AUC = 0.80, Q1–Q3: 0.72–0.87), comparable to the well-established gene repertoire method. In some cases, the conserved gene neighborhoods were particularly accurate, e.g. for predicting chemolithotrophic organisms (Figure 4D). Thus, gene synteny patterns appear broadly associated with many microbial phenotypic traits.
We examined the spore-forming phenotype in more depth, since the involved genes were well-characterized experimentally (51,60,61) and via analyses of genomic occurrence (62–64). A core set of genes with homologs in nearly all sporulating Firmicutes bacteria was outlined previously (52); curiously, many of these genes also have homologs in their non-sporulating relatives, suggesting they may be retained because of roles in other biological processes (52,62). We hypothesized that their involvement in sporulation, or lack thereof, may be reflected in the genomic neighborhoods of such genes. In our data, the ability to form spores was highly predictable from genomic clustering of commonly occurring COGs (AUC = 0.98), prompting us to examine the association of the genomic proximity between the individual COG pairs with the sporulation phenotype. This yielded 3010 significant COG pairs (Mann–Whitney test; Bonferroni-adjusted P < 0.01; Figure 4C). The best-supported results were enriched in COG pairs involving a known sporulation gene (51) and a ribosomal protein (RP) gene or a translation factor (P < 2 × 10-16, Mann–Whitney test on distribution of P-values; Figure 4C). RP gene operons are well-known to cluster in prokaryotes (65,66) with other genes of apparently unrelated function ‘hitchhiking’ in the vicinity (45). We find many sporulation gene COGs that are proximal to RPs specifically in sporulating, but not in non-sporulating bacteria. For instance, COGs with spoVFB, spoIVB and spoIIIE genes form statistically significant clusters with eight RPs and translation factors (Figure 4B, Supplementary Figure S2C; median FDR across all pairs <1e-6). The known sporulation gene spoIIIJ clusters with a tRNA modification enzyme trmE—which exhibits a strong sporulation knockout phenotype (51)—and with RP genes S18, L9 and L31 (Supplementary Figure S2D). These examples illustrate how gene neighborhoods can be conserved preferentially in organisms exhibiting a particular phenotype, conceivably due to co-regulation (45) or differential gene essentiality (67).
Systematic inference of phenotypes from codon adaptation
Next, we considered the possibility that evolution of codon usage biases within gene families may be predictive of microbial phenotypes. Highly expressed genes are known to be enriched with optimal codons across many bacterial and archaeal genomes (68,69), facilitating efficient and accurate protein translation. This ‘codon adaptation’ of orthologous genes may however change in evolution (70), driving phenotypic divergence (47,71). Such genomic signatures of translation efficiency within certain gene families have been used to infer the adaptive value of individual genes to various environmental niches (6) and we thus hypothesized that the overall pattern of codon adaptation across many genes of an organism can predict its phenotype (Materials and Methods). Indeed, such translation efficiency profiles predicted phenotypes with median AUC = 0.75 (Q1–Q3: 0.67–0.81), when used as features in Random Forest models. Again, a number of models based on codon biases were similarly informative as gene repertoires (Figure 4B; example in Figure 4D), consistent with past work linking evolution of translation efficiency to a broad range of gene functions and phenotypic traits; reviewed in (18).
For each of the five comparative genomics approaches described above, we supply the sets of most important features (Materials and Methods) as Supplementary Table S4, thereby providing information about the functioning of our genome-based models. Importantly, these models draw on statistical associations that are robust and predictive on held-out data, but may not necessarily be representative of the biological mechanisms underlying particular phenotypes.
Validating the accuracy of novel annotations
In summary, we provide predictions for a comprehensive set of 424 microbial traits using six independent text corpora, and five independent sources of genomic data; available at http://protraits.irb.hr/. Of these ∼545 000 novel annotations (at FDR < 10%, including both positive and negative annotations; Supplementary Figure S3), ∼308 000 are supported in two or more independent data sources (Figure 5A). Next, we evaluated a random sample of 2489 phenotype predictions by a literature search performed by two curators. Overall, our FDR estimates (Materials and Methods) appear trustworthy, reflecting researcher judgement well, particularly when requiring agreement of any two independent predictions ('two-votes') for each phenotype (Figure 5C; Supplementary Figure S4; detailed statistics in Supplementary Table S5). Of note, combining the 11 individual predictors using the two-votes scheme consistently provided increased coverage over each individual predictor (Supplementary Figure S5), without trading off accuracy.
The number and accuracy of inferred annotations at different confidence levels. (A) Histograms show the number of organisms covered with traits, tallying both positive and negative labels. Shown for (from left to right): known phenotype labels, knows labels including NMF-inferred traits; using predicted annotations at different stringency levels, requiring agreement in three, two or only a single data source (including both text and genomic data sources). All predictions are at FDR < 10%. (B) The relative proportions of predictions accessible only with certain integration schemes, from the most stringent ("Known") to the least stringent ("1 vote"). (C) The precision (equivalent to 1-FDR) estimated by manual curation of 2489 predictions, obtained at different nominal precision thresholds. Bins by nominal precision on x-axis; observed precision calculated using Wilson point estimate on y-axis. Error bars are 95% C.I. (adjusted Wald). Top panel compares schemes for integrating predictions across the individual data sources; middle panel shows accuracy of the six text sources; bottom panel shows accuracy of the five genomic methods (PC stands for proteome composition, GR for gene repertoires, GN for gene neighborhoods, TE for translation efficiency and MC for metagenome co-occurrence). In all cases, predictions for microbes without previously assigned phenotypes are evaluated (see Supplementary Figure S4 for data on previously known labels, indicating putative mis-annotations); detailed statistics are in Supplementary Table S5.
When examining sets of traits individually, in case of the phenotype labels obtained via the BacMap/NCBI databases, the two-votes annotations with nominal FDR < 10% will have an actual FDR of 5.8%, and similarly so for the GOLD ecosystem annotations (FDR = 7.5%; Supplementary Figure S4A). The phenotypic concepts we inferred from free text are also supported in validation, with FDR = 12.9%. These estimates are highly consistent between the two curators (Supplementary Figure S4B). Furthermore, upon breaking down the annotations by the prediction methodology, we find that the five comparative genomics methods have the actual FDR ranging from 6.6% for the proteome composition to 8.9% for the metagenome co-occurrence, and similarly so for the six text mining corpora (6.8–8.2%; Figure 5C); all values given for the nominal FDR < 10%. At more permissive thresholds that afford a broader coverage with annotations, the nominal FDRs again match the observed FDR found via curation (Figure 5C; Supplementary Figure S4; Supplementary Table S5).
Importantly, we have also applied the models to organisms that already had known labels for a particular phenotypic trait (in crossvalidation; Materials and Methods), following the rationale that some of these initial labels may be incorrect. Indeed, we found that the apparent false positives for models that predict existing NCBI/BacMap phenotypes with high confidence (<10% nominal FDR) are, in fact, in many cases genuine positives (55.1% accuracy from curation; Supplementary Figure S4C). We suggest that automated phenotype prediction methods that use free text and genomic data can be useful in highlighting misannotated phenotypes in existing databases.
Additionally, we examined the infrequent cases of annotations that differed between the source databases (listed in Supplementary Table S1), and found the predictive accuracy to be roughly comparable to the general set of annotations (Supplementary Figure S4C).
The network of co-occurring microbial phenotypic traits
The ProTraits atlas of 424 phenotypic traits annotated to thousands of taxa presents an opportunity for a global, unbiased analysis of the interrelatedness among microbial phenotypic traits. Here, we estimate the distance between any two traits as the overlap between the sets of species that display the one or the other trait, accounting for both the positive and the high-confidence negative phenotype annotations (Supplementary Methods). For instance, free-living bacteria are more often flagellated (odds ratio, OR = 4.2, 95% CI: [3.4, 5.2], P = 10-44) while host-associated ones are associated with ability to grow on sucrose (OR = 4.8, 95% CI: [1.0, 23.5], P = 0.045). We further represented all such pairwise similarities as a phenotype network that describes the convergence between groups of microbial traits (Figure 6, Supplementary Figure S6; Supplementary Methods).
A network of microbial phenotypic traits. Similarity of trait pairs was estimated using the F1-measure (shown as line thickness) that accounts for overlap in both the positive and the negative labels. A clustering was used to partition the network (colors), and selected sections are highlighted. The nodes labeled with comma-separated keywords are the NMF-inferred phenotypic concepts. In Supplementary Figure S6, the full network is shown with labels for all nodes.
Clustering the network highlights well-known associations between phenotypes, such as a densely interconnected cluster of human and animal pathogenicity-related traits (Figure 6, top), or grouping of traits related to food fermentation, lactic acid production and vaginal microflora (Figure 6, right). This serves as a validation of our approach, and suggests that some less-obvious associations between phenotypes might be informative. For instance, secondary metabolite production appears strongly associated with the aerobic lifestyle (Figure 6, right; OR = 124.5, 95% CI: [17.3, 896.6], P = 10-30), as was noted in the literature, even though the basis of this association is elusive (72,73). Moreover, the microaerophilic and the facultative aerobic phenotypes are associated with mentions of bacterial infectivity and pathogenicity (from NMF concept discovery; Figure 6, top): OR = 8.2, P = 0.01; and OR = 8.9, P = 10-58, respectively, consistent with prior evidence linking the transition to microoxic conditions with expression of virulence genes (74). Of note, this phenotypic trait network pertains to the ∼3000 genomes covered by ProTraits. An important caveat is that this network may not necessarily reflect universal trends of microbial trait co-occurrence, since the representation of prokaryotic phyla among the currently sequenced set of genomes is not unbiased (15).
Increased power to detect genes associated with traits
An important application of a comprehensive atlas of phenotypic traits would be to search for genetic underpinnings of the traits of interest. Indeed, past analyses have sought associations between patterns of presence/absence of homologous genes across microbial genomes, and certain phenotypes (3,4,75–79), while focusing on selected sets of traits (n ≤ 10) and genomes (n ≤ 300). Here, we extend each of these dimensions by approximately an order of magnitude, greatly increasing statistical power to detect associations, which is particularly important for rarely occurring phenotypes or gene families.
We tested for statistically significant associations of ∼80 000 COG or NOG gene families to 332 traits, while controlling for confounding effects of the phylogeny (principal components of the 16S rRNA evolutionary tree were included as covariates in logistic regression; see Supplementary Methods) and of genomic %G+C and genome size (44,80). Previously known phenotypic labels (Figure 1) can retrieve 20 348 associations across all COGs, while requiring FDR < 10% (t-test for significance of regression coefficient; Figure 7A) and OR ≥ 4 or ≤ 0.25. However, when combining the previously known with the newly annotated phenotypes, this increases approximately 6-fold, to 116 639 gene-trait associations. About half of those (n = 57 088) are also supported at a stringent FDR < 1%. Of note, phenotype labels used in this test were inferred only from text sources and not from genomics (Materials and Methods). This new, broader set of associations includes, among others, 648 and 1187 very high-confidence associations of genes to pathogenicity to plants or to mammals, respectively, while original phenotypic labels would recover less than one-tenth of these links (61 and 57, respectively; at FDR ≤ 1%). Other traits of high general interest are also represented in the list of significant associations (Supplementary Table S6).
The extended collection of phenotypes reveals many novel associated genes. (A) The number of significant gene-trait associations for various traits (rows), using the previously known phenotypes (dark green) and the newly predicted phenotypes at different stringency levels (lighter shades; broken down by FDR cutoffs). (B) The two gene families encoding catalase genes katG (COG376) and katE (COG753) are associated with the biochemical phenotype describing catalase activity. Size of tiles in the mosaic plot is proportional to the number of organisms with (1) or without (0) a catalase homolog, and those exhibiting (1) or not exhibiting (0) catalase activity. (C) Validation of the predicted gene-phenotype interactions using two known sets of sporulation genes (in B. subtilis) and a known set of flagellar genes. Mosaic plots as above. (D) Same as in (A), but counting the number of epistatic interactions involving pairs of genes and a phenotype. (E) A network of epistatic interactions involving the sporulating phenotype, known B. subtilis sporulation genes (green) and novel sporulation genes (red) that were recently discovered in a B. subtilis genetic screen (51). Edges represent synergistic epistatic interactions, in all cases FDR < 1% and ORinter > 9.2. Edge thickness denotes log OR. (F) Examples of an antagonistic epistatic interaction between the pair of catalase genes katG (COG376) and katE (COG753) and the ‘catalase activity’ phenotype (left), and a synergistic epistatic interaction between the pair of sporulation genes spoIIE (COG3854) and spoIIIAA (COG2208) and the sporulation phenotype (right). Mosaic plots as above. (G) Same as in (A), but counting the number of significantly enriched Gene Ontology (GO) terms for various traits.
Overall, the ProTraits resource yields a 6.6-fold increase in the coverage of microbial gene-trait associations (at FDR < 1%). This suggests that further efforts toward systematically annotating phenotypes are yet to uncover many novel genes underlying the corresponding traits. To investigate this, we simulated sets of phenotype annotations with varying coverage (Materials and Methods). The tally of highly confident gene-trait associations can be approximated by a linear fit to number of available phenotype labels (Supplementary Figure S7) and does not generally appear to saturate when considering all currently available annotations. The slope of the linear fits suggests how many genes might be linked to traits with future increases in phenotype coverage: for every additional microbial species labeled, we estimate that a median of 0.67 novel gene families will be associated with the trait (Q1–Q3: 0.18–1.06; at <1% FDR and odds ratio ≥ 4), thus extending our knowledge of the genetic basis of microbial phenotypes.
Validation of novel gene-trait associations
As a validation, we considered a set of genes with sporulation phenotypes in knockout strains of B. subtilis (51) and also a broader set of genes regulated during sporulation (52). We additionally considered a set of genes previously found to co-occur with the sporulation phenotype via phylogenetic profiling (64). Since sporulation is well-investigated experimentally, these known genes can serve to benchmark our methodology for predicting whether an organism sporulates (S) or does not sporulate (NS). In particular, if the ProTraits pipeline infers this phenotype correctly, its predictions should enhance the statistical power to recover known genes. Gene-phenotype associations using original phenotypic labels (104 S and 524 NS species) can capture 18, 24 and 10 genes from the three validation sets. This increases to 25, 34 and 13 genes (Figure 7C), respectively, upon including additional phenotype annotations predicted from text (FDR < 10%; 164 S and 1156 NS). Importantly, the enrichment with known genes was not significantly different between original and extended sets of phenotype labels (P = 0.7–0.8, Z-test for difference of log OR, given for the three validation sets; Supplementary Figure S8A). This implies that the newly obtained gene-trait associations are similarly accurate as the original ones, while affording broader coverage. Next, we evaluated associations of the gene families known to encode the catalase enzymes, with the ‘catalase activity’ biochemical phenotype. In COG0376 (katG) and COG0753 (katE), the statistical support for the associations increased upon annotating a broader set of microbes as being catalase positive or negative (Figure 7B). Finally, we also examined the COG gene families that appear significantly associated with the ‘flagellated’ phenotype only upon introducing the additional annotations obtained herein, by validating against known B. subtilis flagellar genes (81) (Figure 7C). Again, the enrichment was similar for the original and the extended set of labels: OR = 16.1 and OR = 14.6, respectively (Supplementary Figure S8A). The novel labels allowed us to recover an additional flagellar gene (the chaperone motE) and moreover two flagella-related genes not given in the validation set: the regulator of flagellin synthesis flbT and a chemotaxis signal transduction gene cheF (82). Using only the original set of phenotypic labels, these genes would not have been retrieved.
Next, the genes we linked to individual phenotypic traits were examined for consistency of biological function. For instance, the trait ‘strictly aerobic’ has phenotypic annotations that we found to be associated with 1187 gene families (t-test on logistic regression coefficient, FDR < 10%), which are enriched in the Gene Ontology (GO) terms aerobic respiration (FDR = 6 × 10-5, Fisher's exact test) and the more specific functions tricarboxylic acid cycle (0.3%) and electron transport chain (3%). Next, the gene families associated with the label ‘energy source = photosynthetic’ in our extended set of annotations are enriched in the GO term carbon fixation (FDR = 6%). These expected functional enrichments validate our general approach. Importantly, using previously known labels, the linked genes were enriched in average of 14.7 GO categories per trait (FDR < 10% by Fisher's exact test), increasing to 22.7 GO categories with the extended set of phenotypic labels (Figure 7G). In other words, many previously inaccessible trait-GO associations were revealed. We highlight some examples in Table 1, and moreover provide a comprehensive set of 8180 statistically supported links between microbial phenotypic traits and gene functions or pathways (Supplementary Table S6), which may boost further efforts to elucidate the mechanisms underlying particular microbial phenotypes.
Example GO terms significantly enriched with phenotype-associated gene families
| Phenotype | GO term | FDR original | FDR inferred | OR original | OR inferred |
|---|---|---|---|---|---|
| gelatin hydrolysis | peptidase activity | n.s. | 5% | n.s. | 86.4 |
| iron-reducer | iron ion binding | 22% | 1% | 12.5 | 11.2 |
| halophilic | sodium ion transport | n.s. | 7×10-5 | n.s. | 38.9 |
| energy source = photosynthetic | pigment biosynthesis | 20% | 1% | 5.0 | 9.3 |
| pathogenic in mammals | pathogenesis | 29% | 3% | 4.7 | 8.6 |
| lactic/cheese/food/ferment/milk | galactose metabolic process | NA | 1% | NA | 52.6 |
| legum/symbiosis/rhizobia/nod/stem | nitrogen fixation | NA | 3×10-6 | NA | 1362.6 |
| Phenotype | GO term | FDR original | FDR inferred | OR original | OR inferred |
|---|---|---|---|---|---|
| gelatin hydrolysis | peptidase activity | n.s. | 5% | n.s. | 86.4 |
| iron-reducer | iron ion binding | 22% | 1% | 12.5 | 11.2 |
| halophilic | sodium ion transport | n.s. | 7×10-5 | n.s. | 38.9 |
| energy source = photosynthetic | pigment biosynthesis | 20% | 1% | 5.0 | 9.3 |
| pathogenic in mammals | pathogenesis | 29% | 3% | 4.7 | 8.6 |
| lactic/cheese/food/ferment/milk | galactose metabolic process | NA | 1% | NA | 52.6 |
| legum/symbiosis/rhizobia/nod/stem | nitrogen fixation | NA | 3×10-6 | NA | 1362.6 |
Association of genes to phenotypes is determined by a t-test on the logistic regression coefficient, requiring FDR < 10%. Table shows the examples in which the GO-phenotype association was not significant in the sets of gene families retrieved using the original set of phenotypic labels, but became significant as the inferred phenotypic labels were also considered. The FDR column is by Fisher's exact test for COG-GO association by Fisher's exact test, one-tailed. ‘OR’, odds ratio. ‘n.s.’ denotes associations that did not have statistically significant gene-trait associations; Supplementary Methods). ‘NA’ denotes that NMF phenotypic concepts were not available in the original data.
| Phenotype | GO term | FDR original | FDR inferred | OR original | OR inferred |
|---|---|---|---|---|---|
| gelatin hydrolysis | peptidase activity | n.s. | 5% | n.s. | 86.4 |
| iron-reducer | iron ion binding | 22% | 1% | 12.5 | 11.2 |
| halophilic | sodium ion transport | n.s. | 7×10-5 | n.s. | 38.9 |
| energy source = photosynthetic | pigment biosynthesis | 20% | 1% | 5.0 | 9.3 |
| pathogenic in mammals | pathogenesis | 29% | 3% | 4.7 | 8.6 |
| lactic/cheese/food/ferment/milk | galactose metabolic process | NA | 1% | NA | 52.6 |
| legum/symbiosis/rhizobia/nod/stem | nitrogen fixation | NA | 3×10-6 | NA | 1362.6 |
| Phenotype | GO term | FDR original | FDR inferred | OR original | OR inferred |
|---|---|---|---|---|---|
| gelatin hydrolysis | peptidase activity | n.s. | 5% | n.s. | 86.4 |
| iron-reducer | iron ion binding | 22% | 1% | 12.5 | 11.2 |
| halophilic | sodium ion transport | n.s. | 7×10-5 | n.s. | 38.9 |
| energy source = photosynthetic | pigment biosynthesis | 20% | 1% | 5.0 | 9.3 |
| pathogenic in mammals | pathogenesis | 29% | 3% | 4.7 | 8.6 |
| lactic/cheese/food/ferment/milk | galactose metabolic process | NA | 1% | NA | 52.6 |
| legum/symbiosis/rhizobia/nod/stem | nitrogen fixation | NA | 3×10-6 | NA | 1362.6 |
Association of genes to phenotypes is determined by a t-test on the logistic regression coefficient, requiring FDR < 10%. Table shows the examples in which the GO-phenotype association was not significant in the sets of gene families retrieved using the original set of phenotypic labels, but became significant as the inferred phenotypic labels were also considered. The FDR column is by Fisher's exact test for COG-GO association by Fisher's exact test, one-tailed. ‘OR’, odds ratio. ‘n.s.’ denotes associations that did not have statistically significant gene-trait associations; Supplementary Methods). ‘NA’ denotes that NMF phenotypic concepts were not available in the original data.
Widespread epistasis in microbial gene repertoires
Evolutionary innovation in prokaryotes often occurs by gene losses and by gains via horizontal transfer, which are subject to constraints wherein the ability to gain a gene depends on the prior gene content (83). We thus examine genetic interactions involving presence/absence patterns of gene families across genomes, focusing on those which become apparent only in the context of a particular phenotypic trait. This bears a certain analogy to large-scale yeast experiments that systematically examined phenotypes of double gene knockouts (84), but instead relies on comparative genomics to gauge the relevance of epistasis for complex traits (85). Crucially, our extended set of annotations for many microorganisms affords more statistical power to search for gene interactions associated with diverse traits. We focus on pairs of 2663 commonly occurring COG/NOG gene families, in relation to traits defined above, while controlling for confounding effects of phylogeny (covariates in logistic regression; Supplementary Methods). Overall, we find that epistatic interactions are pervasive in prokaryotic gene repertoires and appear broadly associated with microbial traits. In particular, we find 3.9 × 106 significant three-way (gene-gene-trait) interactions (at FDR < 1%) after having included our phenotypic annotations found by text mining; this is a 6.8-fold increase over the coverage that could be obtained with original annotations only (Figure 7D). Many of the tested traits (76%, 254 out of 332) have at least one significant epistatic interaction associated with them, and 230 have five or more interactions (Supplementary Figure S8B). Conversely, almost all (99.4%) of the tested COGs are involved in an epistatic interaction with respect to at least one trait. Importantly, the phenotypic concepts we inferred from free text de novo also provided a context for the genetic interactions: 53 of 73 tested NMF concepts (Materials and Methods) have ≥5 interactions; see above. This suggests that the bottleneck in determining the genetic underpinnings of microbial traits is not only the number of microbes annotated with common phenotypes, but also the incompleteness of our current dictionary of phenotypic traits.
Epistatic interactions associate genes with diverse phenotypes
Of the detected genetic interactions, 56% are instances of antagonistic epistasis, where the phenotype is associated with gene mutual exclusivity, while the remaining 44% are synergistic epistasis, here meaning that the phenotype is associated with gene co-occurrence. An example of the former are the gene families of the catalases katG and katE, where their joint occurrence is less strongly associated with the catalase phenotype than expected from individual effects (Figure 7F; P = 2 × 10−5, t-test for significance of regression coefficient of gene interaction term). This can be interpreted as either one or the other COG being sufficient to attain a measurable catalase activity, with little additional benefit in having both. Next, many prominent examples of synergistic epistasis were observed among sporulation genes, where, for instance, the spoIIE and spoIIIAA gene families are jointly associated with the spore-forming phenotype much more strongly than expected from their individual effects (Figure 7F, P = 3 × 10−8, t-test on interaction term; logistic regression OR = 6.8 and 3.7 for individual COGs and ORinter = 45.2 for the interaction term).
Such synergistic relationships suggest that the two proteins may function as subsequent steps in the same molecular pathway. We found further instances of epistatic interactions that involve one gene known to be involved in sporulation, and other genes where such a role was not previously described. In particular, the gene minC has a known role in septum formation during B. subtilis sporulation, and it is strongly linked with, among others, the gene families containing the transcription factor ydiH (or rex), the pyridoxine biosynthesis enzyme yaaD (or pdxS) and the oligoribonuclease nrnA; in all three cases, P < 5 × 10−5 and ORinter > 28. A recent transposon mutagenesis screen in B. subtilis has indeed found ydiH, yaaD and nrnA to exhibit very strong sporulation phenotypes (51). We also detected epistatic relationships between well-known sporulation genes and three other novel genes (Figure 7E) that were validated in experimental data (51). These examples suggest that signatures of epistasis in gene repertoires may be broadly useful in assigning novel roles to poorly characterized prokaryotic genes, similarly as large-scale experimental genetic interaction screens have proven to be for model eukaryotes (86,87). Moreover, the comparative genomics approach can in principle screen many phenotypes at once, thus highlighting functional links between genes and additionally associating the linked genes to a particular phenotype.
DISCUSSION
There is an unmet need for automated methods that systematically collect and compute over phenotypic data (1,88). Microbes represent a particular challenge in that respect, given the staggering diversity of traits observed across the bacterial and archaeal domains of life. A wealth of phenotype information is contained within the text of scientific articles and other online sources, which is not directly amenable to statistical analysis – implying, most prominently, elucidating the genetic determinants of various traits. Indeed, mining of free-text shows great promise for establishing gene-trait associations, as previously demonstrated by linking individual gene families to PubMed keywords (56). We aimed to further exploit these free-text resources by employing machine learning to accurately assign known phenotypes to novel microbes, as well as to infer phenotypic concepts de novo from free text, while focusing on robust trends that replicate across individual text corpora.
In addition to the abundance of text, a further opportunity for systematic inference of microbial phenotypes is presented by the availability of many complete genomes, which can provide independent support for the text-based predictions. An example is the Genome Properties system (89,90), which applies curated rules to detect activity of biochemical pathways or other molecular subsystems, based on the occurrence of critical gene families. The ability of such rule-based systems to serve as accurate predictors for certain (mechanistically well-understood) traits provides the rationale for more generally attempting to predict a broader spectrum of phenotypes from gene content. To this end, we have employed a general-purpose statistical learning method, the Random Forest. The algorithm chooses the combination of gene families that best explains a particular phenotype, in terms of high statistical support, and uses that combination to annotate the phenotype—or absence thereof—in novel genomes. Since this does not require human input, it scales well with the increase in number of phenotypic traits and sequenced genomes. The ProTraits pipeline can thus generate new predictive models in an automated manner, not requiring expert knowledge on the genetic basis of the phenotypes, which is still lacking for many complex traits.
Further generalizing this well-known gene content-based methodology, our work demonstrates how other comparative genomics approaches normally used to predict gene function—here, synteny patterns (58,59) and codon adaptation (6,70)—can be efficiently repurposed into phenotype predictors. Here, given that the resulting models may be complex and thus challenging to interpret, it is critical to ensure that the predictions they provide are trustworthy and that they robustly classify novel data. We systematically gauged the accuracy of the FDR estimates provided for the inferences in ProTraits by manually evaluating a large sample of the predictions by literature review (Figure 5C; Supplementary Figure S4, Supplementary Table S5). Since the supplied confidence estimates have a probabilistic interpretation which we have validated, they enable the users of the ProTraits resource to make informed decisions about how best to use this massive data set in their work.
We provide a summary of the ∼545 000 phenotypic labels assigned to 3046 microbial species via the ProTraits pipeline by outlining the structure of the microbial phenotype co-occurrence network. Furthermore, we perform a systematic search for the genetic underpinnings of various prokaryotic traits. The ∼57 000 significant gene family-trait links represent a 6.6-fold increase over the genes that could be implicated by using only the previously available databases. Our data suggest that the availability of phenotypic labels is, in many cases, still limiting for elucidating the genetic basis of the traits, which commonly involves epistatic interactions between genes. Thus, it is imperative to direct further effort toward systematically annotating microbial phenotypes, thereby matching the strides made toward systematic sequencing of prokaryotic genomes (15). We anticipate that future developments in the natural language processing methods that can annotate semi-structured or free-text (91–94) will make important contributions towards organizing the legacy phenotype data scattered throughout the scientific literature into structured, computable formats.
SUPPLEMENTARY DATA
Supplementary Data are available at NAR Online.
ACKNOWLEDGEMENTS
Authors contributions: M.B. collected and processed the data sets and performed all statistical analyses. M.P. developed the ProTraits Web application. V.V. contributed genomic data. A.K. and F.S. curated phenotypic predictions for validation. F.S. conceived and supervised the project. T.S. supervised the project. F.S. and M.B. interpreted the data and wrote the manuscript.
FUNDING
This work was funded by the Croatian Science Foundation grants HRZZ-9623 (DescriptiveInduction) and HRZZ-5660 (MultiCaST) and by the European Union FP7 grants ICT-2013-612944 (MAESTRA) and REGPOT-2012-2013-1-316289 (InnoMol). FS acknowledges the support of the Spanish Ministry of Economy and Competitiveness, ‘Centro de Excelencia Severo Ochoa 2013-2017’ (SEV-2012-0208) and the FP7 project 4DCellFate (277899). Funding for open access charge: EU FP7 FET ICT-2013-612944 (MAESTRA).
Conflict of interest statement. None declared.
Comments