MOTIVATION: Several applications in constraint-based modelling can be mathematically formulated as cardinality optimization problems involving the minimization or maximization of the number of nonzeros in a vector. These problems include testing for stoichiometric consistency, testing for flux consistency, testing for thermodynamic flux consistency, computing sparse solutions to flux balance analysis problems and computing the minimum number of constraints to relax to render an infeasible flux balance analysis problem feasible. Such cardinality optimization problems are computationally complex, with no known polynomial time algorithms capable of returning an exact and globally optimal solution. RESULTS: By approximating the zero-norm with nonconvex continuous functions, we reformulate a set of cardinality optimization problems in constraint-based modelling into a difference of convex functions. We implemented and numerically tested novel algorithms that approximately solve the reformulated problems using a sequence of convex programs. We applied these algorithms to various biochemical networks and demonstrate that our algorithms match or outperform existing related approaches. In particular, we illustrate the efficiency and practical utility of our algorithms for cardinality optimization problems that arise when extracting a model ready for thermodynamic flux balance analysis given a human metabolic reconstruction. AVAILABILITY AND IMPLEMENTATION: Open source scripts to reproduce the results are here https://github.com/opencobra/COBRA.papers/2023_cardOpt with general purpose functions integrated within the COnstraint-Based Reconstruction and Analysis toolbox: https://github.com/opencobra/cobratoolbox.

MOTIVATION: Constraint-Based Reconstruction and Analysis (COBRA) is a widely used approach for the interrogation and stratification of microbiome samples, yet applications to large-scale cohorts are hampered by limited scalability and efficiency of simulations. RESULTS: We substantially improved the computation speed and scalability of a previous implementation for the construction and interrogation of personalized constraint-based microbiome models as well as implemented additional functionalities for analysis and visualization. AVAILABILITY AND IMPLEMENTATION: Microbiome Modelling Toolbox and tutorials are freely available as part of the COBRA Toolbox at https://git.io/microbiomeModelingToolbox.

MOTIVATION: Genome-scale metabolic reconstructions have been assembled for thousands of organisms using a wide range of tools. However, metabolite annotations, required to compare and link metabolites between reconstructions, remain incomplete. Here, we aim to further extend metabolite annotation coverage using various databases and chemoinformatic approaches. RESULTS: We developed a COBRA toolbox extension, deemed MetaboAnnotator, which facilitates the comprehensive annotation of metabolites with database independent and dependent identifiers, obtains molecular structure files, and calculates metabolite formula and charge at pH 7.2. The resulting metabolite annotations allow for subsequent cross-mapping between reconstructions and mapping of, e.g., metabolomic data. AVAILABILITY AND IMPLEMENTATION: MetaboAnnotator and tutorials are freely available at https://github.com/opencobra. SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online.

Abstract Motivation Manual curation of genome-scale reconstructions is laborious, yet existing automated curation tools do not typically take species-specific experimental and curated genomic data into account. Results We developed Data-drivEn METabolic nEtwork Refinement (DEMETER), a Constraint-Based Reconstruction and Analysis (COBRA) Toolbox extension, which enables the efficient, simultaneous refinement of thousands of draft genome-scale reconstructions, while ensuring adherence to the quality standards in the field, agreement with available experimental data and refinement of pathways based on manually refined genome annotations. Availability and implementation DEMETER and tutorials are freely available at https://github.com/opencobra. Supplementary information Supplementary data are available at Bioinformatics online.

Ethanol is one of the most widely used recreational substances in the world and due to its ubiquitous use, ethanol abuse has been the cause of over 3.3 million deaths each year. In addition to its effects, ethanol's primary metabolite, acetaldehyde, is a carcinogen that can cause symptoms of facial flushing, headaches, and nausea. How strongly ethanol or acetaldehyde affects an individual depends highly on the genetic polymorphisms of certain genes. In particular, the genetic polymorphisms of mitochondrial aldehyde dehydrogenase, ALDH2, play a large role in the metabolism of acetaldehyde. Thus, it is important to characterize how genetic variations can lead to different exposures and responses to ethanol and acetaldehyde. While the pharmacokinetics of ethanol metabolism through alcohol dehydrogenase have been thoroughly explored in previous studies, in this paper, we combined a base physiologically-based pharmacokinetic (PBPK) model with a whole-body genome-scale model (WBM) to gain further insight into the effect of other less explored processes and genetic variations on ethanol metabolism. This combined model was fit to clinical data and used to show the effect of alcohol concentrations, organ damage, ALDH2 enzyme polymorphisms, and ALDH2-inhibiting drug disulfiram on ethanol and acetaldehyde exposure. Through estimating the reaction rates of auxiliary processes with dynamic Flux Balance Analysis, The PBPK-WBM was able to navigate around a lack of kinetic constants traditionally associated with PK modelling and demonstrate the compensatory effects of the body in response to decreased liver enzyme expression. Additionally, the model demonstrated that acetaldehyde exposure increased with higher dosages of disulfiram and decreased ALDH2 efficiency, and that moderate consumption rates of ethanol could lead to unexpected accumulations in acetaldehyde. This modelling framework combines the comprehensive steady-state analyses from genome-scale models with the dynamics of traditional PK models to create a highly personalized form of PBPK modelling that can push the boundaries of precision medicine.

MOTIVATION: The application of constraint-based modeling to functionally analyze metagenomic data has been limited so far, partially due to the absence of suitable toolboxes. RESULTS: To address this gap, we created a comprehensive toolbox to model (i) microbe-microbe and host-microbe metabolic interactions, and (ii) microbial communities using microbial genome-scale metabolic reconstructions and metagenomic data. The Microbiome Modeling Toolbox extends the functionality of the constraint-based reconstruction and analysis toolbox. AVAILABILITY AND IMPLEMENTATION: The Microbiome Modeling Toolbox and the tutorials at https://git.io/microbiomeModelingToolbox.

BACKGROUND: The maturation of the female germ cell, the oocyte, requires the synthesis and storing of all the necessary metabolites to support multiple divisions after fertilization. Oocyte maturation is only possible in the presence of surrounding, diverse, and changing layers of somatic cells. Our understanding of metabolic interactions between the oocyte and somatic cells has been limited due to dynamic nature of ovarian follicle development, thus warranting a systems approach. RESULTS: Here, we developed a genome-scale metabolic model of the mouse ovarian follicle. This model was constructed using an updated mouse general metabolic model (Mouse Recon 2) and contains several key ovarian follicle development metabolic pathways. We used this model to characterize the changes in the metabolism of each follicular cell type (i.e., oocyte, granulosa cells, including cumulus and mural cells), during ovarian follicle development in vivo. Using this model, we predicted major metabolic pathways that are differentially active across multiple follicle stages. We identified a set of possible secreted and consumed metabolites that could potentially serve as biomarkers for monitoring follicle development, as well as metabolites for addition to in vitro culture media that support the growth and maturation of primordial follicles. CONCLUSIONS: Our systems approach to model follicle metabolism can guide future experimental studies to validate the model results and improve oocyte maturation approaches and support growth of primordial follicles in vitro.

Gastrointestinal side effects are among the most common classes of adverse reactions associated with orally absorbed drugs. These effects decrease patient compliance with the treatment and induce undesirable physiological effects. The prediction of drug action on the gut wall based on in vitro data solely can improve the safety of marketed drugs and first-in-human trials of new chemical entities. We used publicly available data of drug-induced gene expression changes to build drug-specific small intestine epithelial cell metabolic models. The combination of measured in vitro gene expression and in silico predicted metabolic rates in the gut wall was used as features for a multilabel support vector machine to predict the occurrence of side effects. We showed that combining local gut wall-specific metabolism with gene expression performs better than gene expression alone, which indicates the role of small intestine metabolism in the development of adverse reactions. Furthermore, we reclassified FDA-labeled drugs with respect to their genetic and metabolic profiles to show hidden similarities between seemingly different drugs. The linkage of xenobiotics to their transcriptomic and metabolic profiles could take pharmacology far beyond the usual indication-based classifications.

SUMMARY: In constraint-based metabolic modelling, physical and biochemical constraints define a polyhedral convex set of feasible flux vectors. Uniform sampling of this set provides an unbiased characterization of the metabolic capabilities of a biochemical network. However, reliable uniform sampling of genome-scale biochemical networks is challenging due to their high dimensionality and inherent anisotropy. Here, we present an implementation of a new sampling algorithm, coordinate hit-and-run with rounding (CHRR). This algorithm is based on the provably efficient hit-and-run random walk and crucially uses a preprocessing step to round the anisotropic flux set. CHRR provably converges to a uniform stationary sampling distribution. We apply it to metabolic networks of increasing dimensionality. We show that it converges several times faster than a popular artificial centering hit-and-run algorithm, enabling reliable and tractable sampling of genome-scale biochemical networks. AVAILABILITY AND IMPLEMENTATION: https://github.com/opencobra/cobratoolbox . CONTACT: [email protected] or [email protected]. SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online.

Motivation: Flux balance analysis and its variants are widely used methods for predicting steady-state reaction rates in biochemical reaction networks. The exploration of high dimensional networks with such methods is currently hampered by software performance limitations. Results: DistributedFBA.jl is a high-level, high-performance, open-source implementation of flux balance analysis in Julia. It is tailored to solve multiple flux balance analyses on a subset or all the reactions of large and huge-scale networks, on any number of threads or nodes. Availability and Implementation: The code is freely available on github.com/opencobra/COBRA.jl. The documentation can be found at opencobra.github.io/COBRA.jl. Contact: [email protected]. Supplementary information: Supplementary data are available at Bioinformatics online.

Motivation: A genome-scale reconstruction of human metabolism, Recon 2, is available but no interface exists to interactively visualize its content integrated with omics data and simulation results. Results: We manually drew a comprehensive map, ReconMap 2.0, that is consistent with the content of Recon 2. We present it within a web interface that allows content query, visualization of custom datasets and submission of feedback to manual curators. Availability and Implementation: ReconMap can be accessed via http://vmh.uni.lu , with network export in a Systems Biology Graphical Notation compliant format released under a Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 International License. A Constraint-Based Reconstruction and Analysis (COBRA) Toolbox extension to interact with ReconMap is available via https://github.com/opencobra/cobratoolbox . Contact: [email protected].

The metabolic phenotype of cancer cells is reflected by the metabolites they consume and by the byproducts they release. Here, we use quantitative, extracellular metabolomic data of the NCI-60 panel and a novel computational method to generate 120 condition-specific cancer cell line metabolic models. These condition-specific cancer models used distinct metabolic strategies to generate energy and cofactors. The analysis of the models' capability to deal with environmental perturbations revealed three oxotypes, differing in the range of allowable oxygen uptake rates. Interestingly, models based on metabolomic profiles of melanoma cells were distinguished from other models through their low oxygen uptake rates, which were associated with a glycolytic phenotype. A subset of the melanoma cell models required reductive carboxylation. The analysis of protein and RNA expression levels from the Human Protein Atlas showed that IDH2, which was an essential gene in the melanoma models, but not IDH1 protein, was detected in normal skin cell types and melanoma. Moreover, the von Hippel-Lindau tumor suppressor (VHL) protein, whose loss is associated with non-hypoxic HIF-stabilization, reductive carboxylation, and promotion of glycolysis, was uniformly absent in melanoma. Thus, the experimental data supported the predicted role of IDH2 and the absence of VHL protein supported the glycolytic and low oxygen phenotype predicted for melanoma. Taken together, our approach of integrating extracellular metabolomic data with metabolic modeling and the combination of different network interrogation methods allowed insights into the metabolism of cells.

Recent advances focusing on the metabolic interactions within and between cellular populations have emphasized the importance of microbial communities for human health. Constraint-based modeling, with flux balance analysis in particular, has been established as a key approach for studying microbial metabolism, whereas individual-based modeling has been commonly used to study complex dynamics between interacting organisms. In this study, we combine both techniques into the R package BacArena (https://cran.r-project.org/package=BacArena) to generate novel biological insights into Pseudomonas aeruginosa biofilm formation as well as a seven species model community of the human gut. For our P. aeruginosa model, we found that cross-feeding of fermentation products cause a spatial differentiation of emerging metabolic phenotypes in the biofilm over time. In the human gut model community, we found that spatial gradients of mucus glycans are important for niche formations which shape the overall community structure. Additionally, we could provide novel hypothesis concerning the metabolic interactions between the microbes. These results demonstrate the importance of spatial and temporal multi-scale modeling approaches such as BacArena.

MOTIVATION: Genome-scale metabolic reconstructions summarize current knowledge about a target organism in a structured manner and as such highlight missing information. Such gaps can be filled algorithmically. Scalability limitations of available algorithms for gap filling hinder their application to compartmentalized reconstructions. RESULTS: We present fastGapFill, a computationally efficient tractable extension to the COBRA toolbox that permits the identification of candidate missing knowledge from a universal biochemical reaction database (e.g. Kyoto Encyclopedia of Genes and Genomes) for a given (compartmentalized) metabolic reconstruction. The stoichiometric consistency of the universal reaction database and of the metabolic reconstruction can be tested for permitting the computation of biologically more relevant solutions. We demonstrate the efficiency and scalability of fastGapFill on a range of metabolic reconstructions. AVAILABILITY AND IMPLEMENTATION: fastGapFill is freely available from http://thielelab.eu. SUPPLEMENTARY INFORMATION: Supplementary data are available at Bioinformatics online.

BACKGROUND: Biological processes such as metabolism, signaling, and macromolecular synthesis can be modeled as large networks of biochemical reactions. Large and comprehensive networks, like integrated networks that represent metabolism and macromolecular synthesis, are inherently multiscale because reaction rates can vary over many orders of magnitude. They require special methods for accurate analysis because naive use of standard optimization systems can produce inaccurate or erroneously infeasible results. RESULTS: We describe techniques enabling off-the-shelf optimization software to compute accurate solutions to the poorly scaled optimization problems arising from flux balance analysis of multiscale biochemical reaction networks. We implement lifting techniques for flux balance analysis within the openCOBRA toolbox and demonstrate our techniques using the first integrated reconstruction of metabolism and macromolecular synthesis for E. coli. CONCLUSION: Our techniques enable accurate flux balance analysis of multiscale networks using off-the-shelf optimization software. Although we describe lifting techniques in the context of flux balance analysis, our methods can be used to handle a variety of optimization problems arising from analysis of multiscale network reconstructions.

MOTIVATION: In flux balance analysis of genome scale stoichiometric models of metabolism, the principal constraints are uptake or secretion rates, the steady state mass conservation assumption and reaction directionality. Here, we introduce an algorithmic pipeline for quantitative assignment of reaction directionality in multi-compartmental genome scale models based on an application of the second law of thermodynamics to each reaction. Given experimental or computationally estimated standard metabolite species Gibbs energy and metabolite concentrations, the algorithms bounds reaction Gibbs energy, which is transformed to in vivo pH, temperature, ionic strength and electrical potential. RESULTS: This cross-platform MATLAB extension to the COnstraint-Based Reconstruction and Analysis (COBRA) toolbox is computationally efficient, extensively documented and open source. AVAILABILITY: http://opencobra.sourceforge.net.

MOTIVATION: Genome-scale metabolic networks are widely used in systems biology. However, to date no freely available tool exists that ensures quality control during the reconstruction process. RESULTS: Here, we present a COBRA toolbox extension, rBioNet, enabling the construction of publication-level biochemical networks while enforcing necessary quality control measures. rBioNet has an intuitive user interface facilitating the reconstruction process for novices and experts. AVAILABILITY: The rBioNet is freely available from http://opencobra.sourceforge.net.

BACKGROUND: Flux variability analysis is often used to determine robustness of metabolic models in various simulation conditions. However, its use has been somehow limited by the long computation time compared to other constraint-based modeling methods. RESULTS: We present an open source implementation of flux variability analysis called fastFVA. This efficient implementation makes large-scale flux variability analysis feasible and tractable allowing more complex biological questions regarding network flexibility and robustness to be addressed. CONCLUSIONS: Networks involving thousands of biochemical reactions can be analyzed within seconds, greatly expanding the utility of flux variability analysis in systems biology.

Recent advances in reconstruction and analytical methods for signaling networks have spurred the development of large-scale models that incorporate fully functional and biologically relevant features. An extended reconstruction of the human Toll-like receptor signaling network is presented herein. This reconstruction contains an extensive complement of kinases, phosphatases, and other associated proteins that mediate the signaling cascade along with a delineation of their associated chemical reactions. A computational framework based on the methods of large-scale convex analysis was developed and applied to this network to characterize input-output relationships. The input-output relationships enabled significant modularization of the network into ten pathways. The analysis identified potential candidates for inhibitory mediation of TLR signaling with respect to their specificity and potency. Subsequently, we were able to identify eight novel inhibition targets through constraint-based modeling methods. The results of this study are expected to yield meaningful avenues for further research in the task of mediating the Toll-like receptor signaling network and its effects.

Metabolic network reconstructions represent valuable scaffolds for '-omics' data integration and are used to computationally interrogate network properties. However, they do not explicitly account for the synthesis of macromolecules (i.e., proteins and RNA). Here, we present the first genome-scale, fine-grained reconstruction of Escherichia coli's transcriptional and translational machinery, which produces 423 functional gene products in a sequence-specific manner and accounts for all necessary chemical transformations. Legacy data from over 500 publications and three databases were reviewed, and many pathways were considered, including stable RNA maturation and modification, protein complex formation, and iron-sulfur cluster biogenesis. This reconstruction represents the most comprehensive knowledge base for these important cellular functions in E. coli and is unique in its scope. Furthermore, it was converted into a mathematical model and used to: (1) quantitatively integrate gene expression data as reaction constraints and (2) compute functional network states, which were compared to reported experimental data. For example, the model predicted accurately the ribosome production, without any parameterization. Also, in silico rRNA operon deletion suggested that a high RNA polymerase density on the remaining rRNA operons is needed to reproduce the reported experimental ribosome numbers. Moreover, functional protein modules were determined, and many were found to contain gene products from multiple subsystems, highlighting the functional interaction of these proteins. This genome-scale reconstruction of E. coli's transcriptional and translational machinery presents a milestone in systems biology because it will enable quantitative integration of '-omics' datasets and thus the study of the mechanistic principles underlying the genotype-phenotype relationship.

ABSTRACT: BACKGROUND: The set of extreme pathways (ExPa), {pi}, defines the convex basis vectors used for the mathematical characterization of the null space of the stoichiometric matrix for biochemical reaction networks. ExPa analysis has been used for a number of studies to determine properties of metabolic networks as well as to obtain insight into their physiological and functional states in silico. However, the number of ExPas, p = |{pi}|, grows with the size and complexity of the network being studied, and this poses a computational challenge. For this study, we investigated the relationship between the number of extreme pathways and simple network properties. RESULTS: We established an estimating function for the number of ExPas using these easily obtainable network measurements. In particular, it was found that log [p] had an exponential relationship with log[ summation operatori=1Rd-id+ici] MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacyGGSbaBcqGGVbWBcqGGNbWzdaWadaqaamaaqadabaGaemizaq2aaSbaaSqaaiabgkHiTmaaBaaameaacqWGPbqAaeqaaaWcbeaakiabdsgaKnaaBaaaleaacqGHRaWkdaWgaaadbaGaemyAaKgabeaaaSqabaGccqWGJbWydaWgaaWcbaGaemyAaKgabeaaaeaacqWGPbqAcqGH9aqpcqaIXaqmaeaacqWGsbGua0GaeyyeIuoaaOGaay5waiaaw2faaaaa@4414@, where R = |Reff| is the number of active reactions in a network, d-i MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqWGKbazdaWgaaWcbaGaeyOeI0YaaSbaaWqaaiabdMgaPbqabaaaleqaaaaa@30A9@ and d+i MathType@MTEF@5@5@+=feaafiart1ev1aaatCvAUfKttLearuWrP9MDH5MBPbIqV92AaeXatLxBI9gBaebbnrfifHhDYfgasaacH8akY=wiFfYdH8Gipec8Eeeu0xXdbba9frFj0=OqFfea0dXdd9vqai=hGuQ8kuc9pgc9s8qqaq=dirpe0xb9q8qiLsFr0=vr0=vr0dc8meaabaqaciaacaGaaeqabaqabeGadaaakeaacqWGKbazdaWgaaWcbaGaey4kaSYaaSbaaWqaaiabdMgaPbqabaaaleqaaaaa@309E@ the incoming and outgoing degrees of the reactions ri in Reff, and ci the clustering coefficient for each active reaction. CONCLUSION: This relationship typically gave an estimate of the number of extreme pathways to within a factor of 10 of the true number. Such a function providing an estimate for the total number of ExPas for a given system will enable researchers to decide whether ExPas analysis is an appropriate investigative tool.