If you don't remember your password, you can reset it by entering your email address and clicking the Reset Password button. You will then receive an email that contains a secure link for resetting your password
If the address matches a valid account an email will be sent to __email__ with instructions for resetting your password
Matching metagenomic and/or metatranscriptomic data, currently often under-used, can be useful reference for metaproteomic tandem mass spectra (MS/MS) data analysis. Here we developed a software pipeline for identification of peptides and proteins from metaproteomic MS/MS data using proteins derived from matching metagenomic (and metatranscriptomic) data as the search database, based on two novel approaches Graph2Pro (published) and Var2Pep (new). Graph2Pro retains and uses uncertainties of metagenome assembly for reference-based MS/MS data analysis. Var2Pep considers the variations found in metagenomic/metatranscriptomic sequencing reads that are not retained in the assemblies (contigs). The new software pipeline provides one stop application of both tools, and it supports the use of metagenome assembly from commonly used assemblers including MegaHit and metaSPAdes. When tested on two collections of multi-omic microbiome data sets, our pipeline significantly improved the identification rate of the metaproteomic MS/MS spectra by about two folds, comparing to conventional contig- or read-based approaches (the Var2Pep alone identified 5.6% to 24.1% more unique peptides, depending on the data set). We also showed that identified variant peptides are important for functional profiling of microbiomes. All results suggested that it is important to take into consideration of the assembly uncertainties and genomic variants to facilitate metaproteomic MS/MS data interpretation.
). Researchers started to investigate the therapeutic applications of microbiome, which are very promising as demonstrated in recent case studies: Zhu and colleagues edited gut microbiota on ameliorate colitis (
), single-cell genomics and community metagenomics revealed that Nitrospinae are the most abundant and globally distributed nitrite-oxidizing bacteria in the ocean, whereas metaproteomic and metatranscriptomic analyses suggest that nitrite oxidation is the main pathway of energy production in Nitrospinae. And in another work, Kleiner and colleagues (
) showed that metaproteomics could be used to assess species biomass contributions in microbial communities, which is less prone to the biases in sequencing-based methods. Studies have shown the impact of microbiomes on human health and diseases. Metaproteomics provides a new opportunity of studying the functionality of microbial communities.
Although metaproteomics revealed complementary information to metagenomic and metatranscriptomic studies, it often applies conventional proteomics techniques using mass spectrometry (MS). Metaproteomic data analysis is challenging (
). Successful proteomic database search largely relies on the completeness and specificity of the target protein database. There are two broad categories of approaches for metaproteomic database search: approaches relying on a generic collection of reference proteins (for example, the human/mouse gut bacterial genes (
)), and approaches that use specific proteins/peptides inferred from matching metagenomic and metatranscriptomic data sets of the same microbial community. Different from typical proteomics data analysis, peptide/protein identification from metaproteomics data generally lacks a good reference database of proteins for database searches, except for well-studied microbiomes. On the other hand, although matching metagenomic or metatranscriptomic data sets can provide more adequate reference database, they are often under-used for metaproteomic data analysis. Microbiomes are normally composed of many species, many of those species are closely related, as demonstrated in studies of microbial communities empowered by long sequencing reads (
). Despite of the recent advances in the development of assembly algorithms for metagenomes, it remains a challenge to assemble metagenomes from shotgun metagenomic sequences because of the nature of microbial communities: a large number of species co-exist in a sample; these species are of various abundances; and many species are closely related containing similar genome sequences with variations so that assemblers won't be able to untangle them. As a result, metagenome assemblies often provided limited references for peptide/protein identification from metaproteomic data sets.
Approaches have been proposed to address the above-mentioned limitation of using metagenome assemblies as reference for peptide/protein identification from MS/MS data. One successful effort is from May and colleagues (
), who developed a method to derive metapeptides (short amino acid sequences that may be encoded by multiple organisms) from shotgun metagenomic sequencing reads, and then use the metapeptides as the reference for peptide identification from metaproteomic MS/MS data. They showed that by constructing site-specific metapeptide databases, they were able to detect more than one and a half times as many peptides as by searching against predicted genes from an assembled metagenome and roughly three times as many peptides as by searching against the NCBI environmental proteome database. Also some benchmarking experiments have been optimized to explore better experimental and computational strategies for metaproteomics (
) as well. We have also recognized the limitation of conventional approach of using assembled metagenomes (e.g. using only the contigs) for peptide/protein identification from MS/MS data and developed a graph-centric approach (Graph2Pro) (
) that uses the assembly graph to take into consideration of the assembly uncertainties. Our approach achieved significant improvement of peptide/protein identification from MS/MS data for microbiome research.
Comparing to (meta)genomic sequencing, metaproteomics (and proteomics in general) typically has relatively limited throughput. It is therefore important to develop methods that can extensively use the sequencing reads (in addition to assemblies) for peptide/protein identification from MS/MS data (sequencing reads from low abundant species are unlikely to be assembled into contigs). Further, microbial communities often contain closely related species and strains with genomic variations, many of which are of low abundances (
). Short sequencing reads containing variations are often “ignored” by the assemblers, so as a result, our graph-based approach that uses assembly graphs for MS/MS identification is unable to use these reads for MS/MS identification. The open database search approach (
) in principle can be applied to detect variations, which however will be impractical for metaproteomic identification because of the unknown and extremely large search space. The goal of our new variant-aware approach is to incorporate the potential peptides encoded by the short sequencing reads that are not even assembled into the assembly graph, to improve the identification from MS/MS data. Instead of using variant callers (such as MIDAS (
), our variant-aware approach retains all potential variants containing peptides, and use them as reference for spectral database search. The rationale is that the variants are often of low abundance, so the typical variant detection tools that rely on sequencing coverage to distinguish variants from sequencing errors may not work. If a potential variant is supported by only few reads, but nevertheless is supported by the MS/MS data, it is a confident variant. We also note that, different from the previous method (
), which includes all potential short peptides from short sequencing reads, our variant-aware approach only includes the peptides that are similar to the proteins that have already been identified from MS/MS spectra. We showed that our approach effectively reduced the search space for MS/MS spectral database searching, and identified more peptides from metaproteomic MS/MS data.
We tested our integrated pipeline for peptide/protein identification using the ocean water and wastewater multi-omics data sets. We compared our approaches with contig-based (in which proteins predicted from contigs are used as the reference), read-based approach (using metapeptides derived from reads as the reference), and a simple combination of read- and contig-based approaches (in which peptides were derived using the read- and contig-based approaches separately, and then merged as the final results). Our results showed that by considering variations and uncertainties of assembly, we can significantly improve the MS/MS spectral identification rate, and therefore identify more peptides (and proteins) for downstream analyses including functional profiling. The Var2Pep approach alone identified 5.6% to 24.1% more unique peptides, and our pipeline integrating the Graph2Pro approach and the Var2Pep approach improved the overall identification rate of the metaproteomic MS/MS spectra by about two folds, as compared with conventional contig- or read-based approaches.
Performance of the Meta-proteogenomic Approach on the Seawater Microbiome Data Set
We first tested the pipeline for peptide and protein identification from metaproteomic MS/MS data using the ocean data set. Table I lists the size of the CS and BSt metapeptide databases, as well as the size of Var2Pep databases for comparison. Apparently, the Var2Pep database contains much fewer reads-derived peptides for MS/MS identification. For instance, for the BSt data set, there are about 15.9 million metapeptides, whereas the Var2Pep database only contains about 2.7 million peptides derived from 148 million mismatched or unmapped reads.
Table ISummary of the read-based MS/MS search databases for the ocean data sets
Number of peptides in metapeptides database (obtained from
Table II summarizes the peptide identification for selected data sets (BSt45 and CS51; see results for all data sets in the supplementary Data File S1) using different approaches. Fig. 2 shows the comparison results for all ocean water data sets. Briefly, our approach approximately doubled the number of identified peptides as compared with the contig-based and read-based approaches, when the same FDR (≤ 0.01, at spectrum level) estimated using the target-decoy search approach was applied. When FDR ≤ 0.01 at peptide level was applied, all methods unsurprisingly identified fewer peptide, however, we observed similar trends of improvements by Graph2Pro and Var2Pep: supplemental Fig. S1 shows the comparison of the different approaches on all ocean water data sets, using FDR ≤ 0.01 at peptide level, and supplemental Fig. S2 shows the side-by-side comparison of the results by using FDR at spectrum or peptide level for the CS51 data set, respectively. It was shown in reference (
) that read-based approach using short amino acid sequences predicted from shotgun sequencing reads as the search database identified significantly more peptides as compared with the contig-based approach using genes predicted from contigs assembled by SOAPdenovo version 1.06 (
)). Our results confirmed the observation that when SOAPdenovo2 was used as the assembler, contig-based approach resulted in worse MS/MS identification (see Table II) than the read-based approach. However, when MegaHit was used as the assembler, the trend reversed and the contig-based approach identified more peptides than the read-based approach (see Fig. 2), whereas in both cases adding Graph2Pro and Var2Pep steps further improved the peptide identification.
Table IISummary of peptide identification from MS/MS spectra for selected ocean datasets
BSt45 (90,072 spectra)
CS51 (100,588 spectra)
Contigs based (S*)
Contigs based (M*)
Graph2Pro (M*) + Var2Pep
Note: S* represents SOAPdenovo2, M* represents MegaHit, PSM stands for peptide spectrum match. Bst stands for Bering strait, and CS stands for Chukchi sea. Graph2Pro (S*) and Graph2Pro (M*) represent using assembly graph from SOAPdenovo2 (S*) and MegaHit (M*) as the reference in Graph2Pro, respectively. In all cases, FDR (false discovery rate) was estimated using a target-decoy search approach, and a cutoff of 1% at spectrum level was applied. This table only shows the results for two datasets. See Fig. 2 and supplementary Data File S1 for results of all data sets.
Our pipeline significantly outperformed the simple combination of contig-based and read-based approaches (red bars versus blue bars in Fig. 2). For instance, for the Bst45 data set, reads-based method identified 7795 (8.6% of total) spectra, and the contig-based (using MegaHit) approach identified 8631 (9.5%) of total spectra. Combining contig-based search and reads-based search resulted in the characterization of a total of 16,426 spectra/3813 peptides. By comparison, our pipeline reported the identification of comparable number of spectra (15,913), but many more unique peptides (6240).
Comparison of Our Approach with Grouped and Cascaded Searches
Our pipeline uses separate database searches on the Graph2Pro and Var2Pep databases, and then combines the search results (we called it separate search for comparison purpose). We compared this approach with the other two approaches: a single spectral search on a grouped database combining both databases (called grouped approach), and cascaded search using only rejected spectra from the spectral search against the Graph2Pro database for the subsequent search against Var2Pep database (i.e. the cascaded approach, similar to the approach developed by Kertesz-Farkas and colleagues (
)). Table III summarizes the comparison of separate and cascaded searches at the same FDR level. For example, for CS51, our pipeline based on separate search against the Var2Pep database identified additional 2073 spectra (917 unique peptides); by comparison, the cascaded search against the Var2Pep database identified slightly fewer 1598 spectra (629 unique peptides). Notably, this result appears different from the previous report that cascaded search identified more spectra and unique peptides (
) than the separate search. One possible reason is that removing some high qualify spectra in the cascaded searches may alter the overall score distribution, and thus the estimation of false discovery rate for different searching results. Our approach also outperformed (marginally) the search using combined database. For example, for the CS51 data set, combined database search only identified a total of 20,046 spectra and 7662 unique peptides; by contrast, our approach based on separate searches identified 31,145 spectra and 12,380 peptides, and the cascaded search identified 30,670 spectra and 12,092 peptides from the same spectra data set. Based on our tests, we recommend using either separate or cascaded search, but not using the combined database for spectrum search. Our pipeline supports both separate search (default) and cascaded-search, and the results reported below were based on the separate searches.
Table IIIComparison of additional spectra and unique peptides identified by searching against the Var2Pep database using the separate search approach (our approach) versus cascaded search
Performance of the Meta-proteogenomic Approach on the Wastewater Microbiome Data Set
Next we tested our approach using the wastewater data sets. Fig. 3 illustrates the results of peptide identification using different approaches on this data set using FDR ≤ 0.01 at spectrum level (see supplementary Data File S2 for details; also see supplemental Fig. S3 for the comparison of the different approaches using FDR ≤ 0.01 at peptide level, and supplemental Fig. S4 for side-by-side comparison of the results by using FDR at spectrum versus peptide level for the SD3-MG data set). For comparison purposes, for reads-based approach, we used the Sixgill software (https://github.com/dhmay/sixgill) (
) to predict putative peptides (i.e. metapeptides) from reads for the reads-based approach. Here are the numbers of metapeptides predicted from the different data sets: 266,430 (SD3-MG), 1,640,859 (SD3-MGMT), 2,614,012 (SD6-MG), 4,782,170 (SD6-MGMT), 1,800,892 (SD7-MG), and 3,977,520 (SD7-MGMT) (sequences of the metapeptides can be found at our supplementary website at http://omics.soic.indiana.edu/MP). In all cases, the contig-based approach achieved significantly better performance than the read-based approach, and our approach integrating Graph2Pro and Var2pep significantly outperformed other approaches.
As shown in Fig. 3, significant improvement of peptide identification was achieved by using both metatranscriptomic and metagenomic data sets as the reference as compared with using only metagenomic sequencing data (see supplementary Data File S2 for details). This result is expected, because metatranscriptomic data may provide sequence information for otherwise low-abundant species that may have been missed by metagenomic sequencing. We note that among the wastewater data sets, the Var2Pep approach resulted in the most significant improvements of the peptide identification for the SD3 data set, with an additional 24.1% and 19.5% more identified unique peptides when using matched metagenomic data set alone (SD3-MG), and both the metagenomic and metatranscriptomic data sets (SD3-MGMT).
Considering the increasing popularity of metagenome assembler MetaSPAdes, we also used the wastewater data sets to test if more contigs assembled from metagenome resulted in better identification of metaproteomics data. Table IV shows the comparison of the results for SD3-MG data set using assemblies from MegaHit and MetaSPAdes as the reference in our pipeline. Strikingly, MetaSPAdes produced about three times more contigs (more than doubled the total bases in the contigs) as compared with MegaHit. However, the increasing of contigs did not lead to more peptide identifications.
Table IVComparison of the performance using metagenome assembly from MegaHit and MetaSpades as the reference on SD3-MG dataset
Abundant Known (and Unknown) Proteins Revealed by Metaproteomics
We examined proteins that were identified using our pipeline. For example, for CS51 data set, its first most abundant protein was supported by 270 spectra (0.8% of the total 31,145 identified spectra for this data set). A total of 663 proteins or protein fragments (of minimal length of 60 aa) were each supported by at least 10 spectra (see supplementary Data File S3 for details). Functional annotation of these proteins by myRAST (http://blog.theseed.org/downloads/myRAST-Intel.dmg) (
) revealed 90 ribosomal proteins, 61 heat shock proteins (HSP60), 59 elongation factors (GTP_EFTU), 38 chaperone proteins (DnaK), 30 glutamate aspartate periplasmic binding protein precursors, 26 branched-chain amino acid ABC transporters, 10 TonB-dependent receptors, among other functions (see supplementary Data File S3 for details of the annotations). Our result is consistent with the fact that many of these functions (e.g. ribosomal proteins) are typically the most abundant and highly conserved in prokaryotic genomes.
Among the CS51's 663 proteins each with support of at least 10 spectra, a significant fraction (93, 14%) have no functional annotation based on myRAST annotation, including a protein (ID: NODE_17410_length_722_cov_188.0000_ID_34819_1_456_+) that has the third most spectra support with 257 spectra and many others with strong spectral supports. See supplementary Data File S3 for the complete list of these highly expressed proteins that lack any functional annotations (and their sequences can be found at our supplementary website). We argue that these proteins are most wanted proteins that wait to be further studied experimentally for their functions. Fig. 4 illustrates two examples of proteins supported by a significant number of variant peptides (details of the variant peptides are available at our supplementary website). Protein366 was predicted from the assembly graph by Graph2Pro and was shown to be supported by a total of 27 peptides, including six variant peptides. The myRAST annotation did not reveal any significant hits for this protein. We also did sequence similarity search of this protein using ncbi blast (http://www.ncbi.nlm.nih.gov/blast) against the nr database (as of Nov 18th, 2018), which returned no hits either. Although no significant sequential similarity was detected between Protein366 and proteins in the public databases, a structural comparison (by FATCAT (
)) revealed that it contains a Cystatin-like fold composed of helical segments packed against coiled antiparallel beta-sheet. Protein2811 is another example, which was annotated as a putative peptide/opine/nickel uptake transporter. This protein has a total of 17 supporting peptides, among which nine are variant peptides.
For the wastewater data sets, an interesting example is a protein (Protein264) identified from SD6 data set. This protein has the most spectra support (2266 spectra), and the domain prediction shows that it contains three SLH (S-layer homology) domains and another larger domain that has unknown function (DUF2815). The presence of SLH domains in this protein suggest that it is likely a cell-surface protein, as studies have shown the utilization of SLH domains as vehicles for surface location of function proteins (
); however, the precise function of this protein remains to be determined. Fig. 4 (bottom) shows the locations of the domains and the distribution of identified peptides and variant peptides supporting this protein.
Low Abundant Proteins with Supports Augmented by Variant Peptides
As shown above, metaproteomics is likely limited to revealing highly abundant proteins expressed from highly abundant bacterial species, so metaproteomics-based functional profiling of microbiomes will only provide a shallow glimpse into the functionality of the underlining microbial community. A hope to alleviate this limitation is to augment the identification of proteins using variant peptides. Here we use the wastewater data set (SD3-MG) as an example to illustrate this approach. We used BlastKOALA (https://www.kegg.jp/blastkoala/) (
) to annotate the functions of identified proteins from this data set. A total of 1401 KEGG families (i.e. the KO families) were identified using proteins with at least one peptide support; however, this number increased to 1938 if variant peptides were also considered (a 38% improvement). We acknowledge that the number of families is still small, because of the shallow functional profiling of microbiomes based on metaproteomics.
We note that the variant peptides identified by our pipeline have both nucleotide level (from metagenomic and/or metatranscriptomic sequences) and mass spectral level supports. Here, we examined the types of substitutions found in the variant peptides. Not surprisingly, we found that most variations are the substitutions likely to be found in homologous proteins, with positive BLOSUM scores (BLOSUM matrices contain log-odd scores of all possible substitution pairs of amino acids, where positive scores represent favorable substitutions observed in homologous proteins such as Asp to Asn substitution and negative scores represent unfavorable substitutions such as Arg to Asp) (
). As shown in Fig. 5, the average BLOSUM (BLOSUM62) score between two amino acids is −2 (see “all” box), whereas the average BLOSUM score of substitutions in variant peptides is positive for Protein366. This result indicates that Var2Pep likely detected biologically meaningful variations.
DISCUSSION AND CONCLUSIONS
We developed the Var2Pep algorithm, which combined with Graph2pro achieved drastic improvement of peptide identification from metaproteomics data. Var2Pep makes use of shotgun sequencing reads for MS/MS identification, and it reduces the search space imposed by the large number of sequencing reads by only keeping peptides that share similarity with proteins already identified by Graph2Pro. As shown in Table I, there were 15.9 million putative peptides in BSt data set from the mismatched/unmatched reads, and Var2Pep only contained 2.7 million peptides for MS/MS database search. This process significantly reduced the number of possible decoys in the database to permit the database searches, which significantly contribute to increasing the identification rate of metaproteomic MS/MS data. Similar ideas have been applied to improve metaproteomics data identification: for example MetaLab (
) uses a reduced gut reference gene catalogue to improve MS/MS identification). We note Var2Pep may still be used when assembly graph is not available (i.e. only the contigs are available). In such a case, peptides encoded in short sequencing reads that share similarity with MS/MS supported peptides identified from the contigs can be included in the database searching for MS/MS analysis.
The main goal of our approaches is to optimize the use of matching metagenomic and/or metatranscriptomic data for metaproteomics data analysis. Our integrated pipeline results in two databases for MS/MS identification, one from Graph2Pro (Graph2Pro DB), and the other from Var2Pep (Var2Pep DB). Our pipeline provides one way of using these databases as shown in Fig. 1. Users can make use of these databases in different ways. In order to get high confident peptide identification, as well as variant peptides, it is highly recommended that users apply strict false discovery rate control, either in spectral level or peptide level. A reasonable concern on our approach of combining the results from two separated searches each at 1% FDR is that the actually FDR would be higher. We re-calculated the FDRs for combined results, based on the unique sets of identified spectra for all the data sets we tested. The results showed that the FDRs estimated using this approach were slightly higher, but still reasonably low (at about 1.5%).
It is possible to combine our approaches with those that are recently developed, including the utilization of spectral clustering to speed up the search as shown in (
) can also be incorporated into the pipeline. Our significantly improved results indicate there is still a huge space for algorithmic improvement in the field of metaproteomics and proteomics in general.
The goal of metaproteomics is not only to identify proteins expressed in the microbial community, but also to estimate their abundances (i.e. their expression levels) under different conditions. Nevertheless, a protein can be quantified only if it can be identified by using the metaproteomic data. Therefore, the methods presented here that increase the coverage of protein identification will also help the subsequent steps for protein quantification. Furthermore, identification of peptide variants may contribute to accurate quantification of the proteins from which these peptides are derived. Tools such as MetaGOmics (
) can be applied to infer functional and taxonomic contents of the microbial communities based on identified proteins/peptides. Better identification of peptides and/or proteins from metaproteomic MS/MS data can improve the functional and pathway annotation of the corresponding microbial communities, as we have shown in (
We showed that metaproteomics (and matching metagenomics/metatranscriptomics) data can be used to identify hypothetical proteins that sometimes are very abundant in the microbial communities. These MS/MS supported hypothetical proteins provide interesting candidates for further studies of their functions using computational and experimental approaches. However, our analyses of the identified proteins also suggested that unless we significantly increase the depth of the metaproteomic experiments, what we identify are limited to those highly abundant (often well studied) proteins.
In conclusion, we developed a pipeline for metaproteomic MS/MS data analysis using matching metagenomic and/or metatranscriptomic sequencing data as the reference. Tests of our pipeline using publicly available meta-omics data sets showed that it is important to consider the assembly uncertainties (captured in assembly graph) and genomic variants to maximize the utilization of metagenomes and/or metatranscriptomes as the reference for metaproteomic data interpretation.