Evaluation and comparison of bioinformatic tools for the enrichment analysis of metabolomics data

1Enrichment techniques for ‘omics’ data are key tools for understanding complex biological systems. These tools reduce the complexity of the data, improve interpretation and understanding of biological systems, and help to generate hypotheses. Although the number of tools for ‘omics’ is rapidly growing, suitable tools for metabolomics are still scarce. Most of the available tools for metabolomics data have been previously developed for other ‘omics’ technologies. These tools have been described in detail elsewhere.

2Enrichment tools denote any analytic technique that benefits from molecular pathway or network information to gain insight into a biological system. The most widely used methodology for performing such analysis is termed functional enrichment or over-representation analysis (ORA). This analysis looks for keywords or descriptors of the set of molecules of interest (e.g. those over-expressed) with respect to a background reference set (e.g. the whole genome/transcriptome/proteome/metabolome or the set of molecules detected by the technology employed). Classical enrichment analyses employ Fisher’s exact test, but many other enrichment methods have derived from it, e.g. hypergeometric, Kolmogorov–Smirnov or Wilcoxon statistical tests.

3To the best of our knowledge studies evaluating the performance of enrichment tools for metabolite sets do not exist yet. The aim of the present work will be to dissect, for the first time, these techniques. First of all, we have carried out an exploratory multivariate analysis of the state-of-the-art of bioinformatic tools for metabolomics sets to visualize their diversity. Then, we have examined the completeness of metabolite databases, the performance of ORA methods and accuracy of disease-based analyses. For these purposes, we have used datasets from metabolomic repositories, whose results have been already published in peer-reviewed journals. In addition, we have enriched selected metabolic pathways and then compared the outputs of these tools when using real datasets or enriched data. Thus the present study provides a global insight of the current status of bioinformatic tools for the analysis and interpretation of metabolite sets from metabolomic studies.

4The list of metabolites used in this work refers to five datasets from metabolomics studies in humans, already published in peer-reviewed journals, whose raw data, study information and the list of identified metabolites are available in MetabolomeXchange, an online portal of metabolomics repositories including MetaboLights, Metabolomics Workbench and Metabolomic Repository Bordeaux. A brief summary of the datasets is shown in Table 1. These datasets correspond to the following publications: 1) Lanza et al.; 2) Fiehn et al.; 3) Kaluarachchi et al.; 4) Hart et al.; and 5) Zhu et al.. P- and adjusted p-values were obtained from the original papers and only metabolites with an adjusted p-value < 0.05 were used for tools comparison. The list of metabolites is shown as (Additional file 1: Table S1).

5Bioinformatic tools for enrichment analysis require the metabolite name or code (identifier) from a metabolite database. Although Kyoto Encyclopaedia of Genes and Genomes (KEGG Compound) identifiers are the most commonly used in metabolomics, some tools prefer other database identifiers such as PubChem, BioCyc/HumanCyc (hereinafter only referred as HumanCyc) or Chemical Entities of Biological Interest (ChEBI).

6We analysed the current completeness of the following metabolite libraries: Human Metabolome Database (HMDB), KEGG, PubChem, HumanCyc, ChEBI, ChemSpider, the metabolic reconstruction Recon2, METLIN and Lipid Metabolites and Pathways Strategy (LipidMAPS).

7The list of significant metabolites from was used to assess the completeness of these nine databases. The identification of metabolites had been carried out by original authors in all the datasets, and in some cases KEGG and HMDB identifiers were already provided by authors. Since the HMDB website provides links to other metabolite databases, we started the search of the identifiers on this site and then we extended this search to the LipidMAPS website and MetExplore (for Recon2 codes). All the identifiers were then double-checked in the corresponding metabolite databases. If more than one metabolite identifier was found (e.g. in PubChem or ChEBI databases), we took all those identifiers and checked in ConsensusPathDB which ones were recognized by the tool (not shown). When the stereochemistry of the metabolite was not specified, the most common chemical configuration was assumed. The complete list of metabolite identifiers is shown in Additional file 1: Table S1.

8Most of the bioinformatic tools for the enrichment of metabolomics datasets accept a list of identifiers as output, while a lower number require quantitative data, e.g. concentration, fold change or peak intensity. Therefore we decided to work with a list of metabolites (name or identifier) from real datasets and enriched data to compare these tools. Although this approach do not allow us to assess some of the available tools, i.e. 3omics or PAPi, the use of simulated o synthetic data would have allowed us to examine a lower number of tools.

9For data enrichment, the dataset with the most metabolites was selected (colorectal cancer, ST000284). The list of significant metabolites of this dataset (n = 42, obtained from) was analysed with MetaboAnalyst, using the option ‘pathway analysis’. MetaboAnalyst’s output was examined and the three KEGG pathways that presented the lowest false discovery rate (FDR), based on the Benjamini-Hochberg procedure, were chosen for pathway enrichment: 1) Alanine, aspartate and glutamate metabolism; 2) Aminoacyl-tRNA biosynthesis, and 3) Arginine and proline metabolism.

10The R package KEGGREST (v.1.17.0) was employed to build an adjacency matrix which linked the metabolites of the dataset (n = 113) with their corresponding KEGG pathways. One was assigned if the metabolite was part of that particular pathway, or 0 if not. Then five metabolites of each pathway were randomly sampled. Enriched data are shown in Additional file 2: Table S2.

11The most commonly used tools for metabolomics data enrichment were chosen for similarity analysis. This selection was formed by 3omics, BioCyc/HumanCyc, ConsensusPathDB, IMPaLA, Ingenuity® Pathway Analysis (IPA®, QIAGEN, Redwood City, CA), KEGG, MassTRIX, MBRole, MetaboAnalyst, Metabox, MetaCore™ (Thomson Reuters Inc., Carlsbad, CA), MetaMapp, MetExplore, MetScape, MPEA, PaintOmics, PAPi, PathVisio, Reactome, Small Molecule Pathway Database (SMPDB), WikiPathways and XCMS.

12The main features of these tools were summarized on a binary matrix (Yes/No responses) including whether they 1) perform ORA, integration with other ‘omics’ or other enrichment analyses; 2) visualization of pathways, networks or other types of visualization; 3) use KEGG, BioCyc, Reactome, Wikipathways, SMPDB or other pathway databases; 4) are databases, programmable, open-source or online tools (Additional file 3: Table S3).

13The similarity analysis was performed with the R package vegan (v.2.4–4). First, the Yes/No responses were transformed to 1 and 0, respectively. Then Jaccard’s coefficients were calculated and a non-metric multidimensional scaling (NMDS) was performed. This method plots dissimilar objects far apart in the ordination space and similar objects close to one another preserving ordering relationships among them.

14The performance of the tools that perform ORA in metabolomics datasets was assessed with the list of significant metabolites of the colorectal cancer dataset (ST000284) (n = 42) and of enriched data. The comparative analysis of ORA tools was performed on tools employing KEGG, Reactome and HumanCyc as pathway database. The selected tools were ConsensusPathDB, HumanCyc, IMPaLA, MBRole, MetaboAnalyst, Metabox, MetExplore, MPEA, PathVisio and Reactome and the pathway mapping tool KEGGREST. Table 2 summarizes the main features of these tools and type of identifiers used. Analyses were performed following the guidelines of each tool.

15The output of these tools was examined for the following metabolic pathways: 1) KEGG: the three aforementioned pathways; 2) Reactome: Metabolism of amino acids and derivates pathway; and 3) HumanCyc: tRNA charging pathway. Ranking (position in the list of pathways sorted by significance), total number of metabolites/pathway, number of hits/pathway, p- and adjusted p-value (generally FDR, calculated by the tools) were recorded from each output.

16Disease-based enrichment analyses were performed by using the list of significant metabolites of the five datasets on: 1) MetaboAnalyst (SMPDB disease pathway database); 2) MBRole (HMDB disease database); 3) IPA® (Ingenuity® disease database); and 4) MetaCore™ (MeSH and OMIM disease databases). Disease, ranking (position according to their p-value), total number of metabolites/disease, number of hits/disease, p- and adjusted p-values were recorded from each output.

17Figure 1 displays a similarity plot of the most commonly used bioinformatic tools. Tools were distributed all along the two dimensions revealing their diversity. The first dimension mainly separated tools that: 1) perform ORA and are non-open source, 2) perform ORA and are open-source, and 3) are a metabolite database. On the other hand, the second dimension mainly separated tools that: 1) perform metabolite identification, 2) perform ORA and are not programmable, and 3) perform ORA and are programmable. MetScape and MetaMapp, which only carry out data visualization, were distant in the plot.

18Metabolites of the five datasets were used to assess the completeness of the metabolite and pathway databases. Almost all the metabolites presented PubChem (97%), ChEBI (91%), METLIN (91%), KEGG (88%), ChemSpider (87%) and HMDB (86%) identifiers, and the 97% of the lipid subset had LipidMAPS identifiers. In some cases, KEGG, HumanCyc and Recon2 provided chemical class identifiers instead of a single identifier to certain metabolites, especially to the lipid subset (Additional file 1: Table S1 and Additional file 4: Table S4).

19In general, ORA methods yielded consistent results using both real and enriched data in all the range of tools tested (Tables 3 and 4). Also similar results were obtained in paired analyses/tools such as MetaboAnalyst hypergeometric test - Fisher’s exact test, MBRole full - Homo sapiens database, MPEA top down - bottom up analyses, or ConsensusPathDB - IMPaLA tools, as expected (Tables 3 and 4).

20Minor differences in the total number of metabolites/pathway and number hits/pathway were found. For instance, MPEA (Table 3) and MBRole (Table 4) presented a higher number of metabolites/pathway than the other tools. Other divergences were also observed, e.g. MPEA provided higher adjusted p-values values (nearly 1 in all the cases) than other tools, or not all the tools mapped the same metabolites of the dataset onto the queried pathways (not shown).

21The significant metabolites of the five datasets were used to analyse the accuracy of the SMPDB, HMDB, IPA®, MeSH and OMIM disease-based libraries. Outputs revealed that the diseases queried (diabetes type 1 and 2, obesity, respiratory alterations and breast and colorectal cancer) were not successfully identified by these tools, as they appeared in a low position in the list of potential diseases and most of the times they presented a p > 0.05 (Table 5).

22Interpretation of metabolomic data is much less straightforward than that with genomic and proteomic datasets. In the present work we have described the diversity of bioinformatic tools for metabolite sets and have evaluated their performance by exploring three features: the completeness of metabolite databases, ORA approaches and disease-based analyses. To that end, we have used five metabolite sets of blood biomarkers of different diseases obtained from LC-MS, GC-MS and NMR metabolomics approaches. This approach allowed minimizing the possible bias introduced by a given metabolomic platform and thus working with a wide range of metabolites.

23Metabolomics is a developing field, thus bioinformatic tools designed to perform enrichment of metabolomics datasets are being developed and released by various groups using diverse statistical tests. Our exploratory multivariate analysis mirrors the high diversity of the currently available tools for the analysis of metabolite sets.

24To date about 30,000 endogenous metabolites have already been identified, but this number is rapidly increasing due to advances in high-throughput technologies. Current metabolite databases do not have the full potential to quickly absorb these advances in the description of the endogenous metabolome yet, as not a single metabolite database used in this work covered the full list of significant metabolites of the five datasets. Among all the metabolites databases, PubChem was the one that covered more metabolites from the datasets. However, PubChem is a crowded compound database and presents duplicated metabolite entries, which might produce a larger number of false positives than searching against the KEGG database. To address the low metabolite coverage of metabolite databases, some of them such as KEGG and HumanCyc assign chemical class identifiers to certain types of compounds, especially lipids such as phosphatidylcholines, sphingomyelins or triglycerides. For instance, KEGG coded phosphatidylcholines and sphingomyelins as ‘C00157’ and ‘C00550’, respectively, and HumanCyc as ‘PHOSPHATIDYLCHOLINE’ and ‘Sphyngomyelin (class)’.

25Missing, ambiguous or redundant entries have been commonly found in public repositories. Indeed metabolites with more than one PubChem, HMDB or ChEBI identifiers were found in this work, which reduce enrichment analyses’ accuracy. Several on-going initiatives on identifiers standardization such as BridgeDB and the Chemical Translation Service are trying to overcome redundancy. Some tools such as MetaboAnalyst, ConsensusPathDB or PathVisio embrace these initiatives and accept different types of identifiers, which are then transformed into an internal identifier prior to the enrichment analysis. However, this approach also presents different pitfalls. For instance, these tools usually transform the input code into KEGG identifiers, and thus certain types of metabolites such as lipids lose their uniqueness and become a chemical class KEGG identifier. Consequently bioinformatic tools analyse these lipids as a single entity, thereby losing the diversity of these metabolites.

26KEGG and HumanCyc are the most used pathway libraries in metabolomics and Reactome is widely used in other ‘omics’ studies. Thus we have evaluated and compared outputs of ORA methods that employ these pathway libraries. Some limitations prior to ORA analysis were found. For instance, despite the fact that almost all the metabolites of ST000284 dataset had a KEGG code, not all of them were mapped in a KEGG pathway. However, these compounds (e.g. 5-hydroxytryptophan and salicylurate) were mapped in other pathway databases such as Reactome, Wikipathways and SMPDB (not shown). In addition, the KEGG code for glutamic acid (C00025) was not recognized by MetaboAnalyst and the alternative suggested by the tool corresponded to the compound amphetamine (C07514).

27The number of total metabolites and hits per pathway varied according to the tool used and those tools that employ the newer database versions (Table 2) presented the higher number of metabolites, as expected. Surprisingly, KEGGREST, a R package that provides an updated client interface to the KEGGREST server, did not provide the highest number of total metabolites among the tested KEGG pathways. Despite regular updates to some pathway databases, such as KEGG or Reactome, being carried out, most of the tools evaluated do not use up-to-date database versions (Table 2). Wadi et al. performed an elegant review on the impact of outdated annotations on pathway enrichment analysis, which revealed that many software tools use functional information not updated for years, thereby strongly affecting the quality of the analyses.

28We can conclude that current ORA methods, despite their differences, provide consistent, robust and reproducible results regardless of their analytic approach (statistical test, p-value adjustment or pathway database used), despite the limitations and small differences found between outputs. The most discordant result was obtained with MPEA, probably due to the fact that it employs a different method to handle many-to-many relationships that may occur between the query compounds and metabolite annotations.

29Although we cannot recommend one tool over the others, we suggest choosing those tools that employ updated metabolite/pathway databases in order to obtain more complete results. Nevertheless, we also consider that the enrichment analysis must not be restricted to a single database or tool. The combined use of libraries such as KEGG, Reactome, HumanCyc or WikiPathways will increase the metabolome coverage and the statistical power of the enrichment analysis.

30Disease-based enrichment analysis did not yield accurate results. Although we only used serum/plasma biomarkers, results with other types of biological samples would have been similar. On one hand, metabolite disease sets are not up-to-date. For instance, MetaboAnalyst and MBRole (SMPDB and HMDB disease databases, respectively) base their searches of literature dated between 1975 and 2008, as stated in the outputs of these tools. Since 2008, advances in high-throughput techniques have remarkably improved metabolomics analyses and, consequently, more knowledge about these diseases is available. As previously discussed, the use of not updated annotation sets strongly affect the quality of the analyses. On the other hand, metabolites can overlap between unrelated physiopathological events since similar metabolic processes are altered. This fact could complicate the development and accuracy of background sets for disease-based enrichment analysis.

31Although extensive work in developing bioinformatic tools for metabolite sets has been carried out in recent years, more effort in improving metabolite/pathway databases and tools is still needed. On one hand, metabolite databases have to rapidly absorb new information from unstoppable advances in high-throughput technologies. On the other hand, enrichment methods should include a wider range of metabolite identifiers (e.g. LipidMAPS, ChemSpider or METLIN) and metabolite pathway databases in order to increase the metabolome coverage. For instance, the LipidMAPS Structure Database contains about 30,000 human endogenous lipids and 12,000 plant lipids, but also databases based on lipid metabolism and signalling pathways, MS/MS spectra and protein-related data. ChemSpider is a general chemical database and offers access to information for almost 25 million experimentally determined structures of natural and synthetic compounds. However, similarly to PubChem, ChemSpider may lead to a high number of false positives. The METLIN database includes nearly 1,000,000 molecules, ranging from lipids, steroids, plant & bacteria metabolites, small peptides, carbohydrates, exogenous drugs/metabolites and central carbon metabolites, and more than 200,000 MS/MS spectra. Including these information sources in current bioinformatic tools would also involve more effort in the improvement of metabolite identifiers converters. Therefore, there is still a long way ahead to achieve complete metabolite and pathway databases and thus accurate enrichment analyses of metabolite sets.

32We have extensively reviewed, for the first time, the state-of-the-art of bioinformatic tools for the enrichment of metabolite sets from metabolomics studies, visualized their diversity, and examined their performance. The redundancy of identifiers, the use of chemical class identifiers and the incompleteness of metabolite databases and disease metabolite sets limit the extent of the analyses and reduce their accuracy. In general, ORA tools provided consistent results among tools revealing that these analyses are robust and reproducible regardless of their analytic approach. However, more work in the completeness of metabolite/pathway databases is required to get more accurate and global insights of the metabolome.