StressME: Unified computing framework of Escherichia coli metabolism, gene expression, and stress responses.

Generalist microbes have adapted to a multitude of environmental stresses through their integrated stress response system. Individual stress responses have been quantified by E. coli metabolism and expression (ME) models under thermal, oxidative and acid stress, respectively. However, the systematic quantification of cross-stress & cross-talk among these stress responses remains lacking. Here, we present StressME: the unified stress response model of E. coli combining thermal (FoldME), oxidative (OxidizeME) and acid (AcidifyME) stress responses. StressME is the most up to date ME model for E. coli and it reproduces all published single-stress ME models. Additionally, it includes refined rate constants to improve prediction accuracy for wild-type and stress-evolved strains. StressME revealed certain optimal proteome allocation strategies associated with cross-stress and cross-talk responses. These stress-optimal proteomes were shaped by trade-offs between protective vs. metabolic enzymes; cytoplasmic vs. periplasmic chaperones; and expression of stress-specific proteins. As StressME is tuned to compute metabolic and gene expression responses under mild acid, oxidative, and thermal stresses, it is useful for engineering and health applications. The modular design of our open-source package also facilitates model expansion (e.g., to new stress mechanisms) by the computational biology community.

Novel context-specific genome-scale modelling explores the potential of triacylglycerol production by Chlamydomonas reinhardtii.

Gene expression data of cell cultures is commonly measured in biological and medical studies to understand cellular decision-making in various conditions. Metabolism, affected but not solely determined by the expression, is much more difficult to measure experimentally. Finding a reliable method to predict cell metabolism for expression data will greatly benefit metabolic engineering. We have developed a novel pipeline, OVERLAY, that can explore cellular fluxomics from expression data using only a high-quality genome-scale metabolic model. This is done through two main steps: first, construct a protein-constrained metabolic model (PC-model) by integrating protein and enzyme information into the metabolic model (M-model). Secondly, overlay the expression data onto the PC-model using a novel two-step nonconvex and convex optimization formulation, resulting in a context-specific PC-model with optionally calibrated rate constants. The resulting model computes proteomes and intracellular flux states that are consistent with the measured transcriptomes. Therefore, it provides detailed cellular insights that are difficult to glean individually from the omic data or M-model alone. We apply the OVERLAY to interpret triacylglycerol (TAG) overproduction by Chlamydomonas reinhardtii, using time-course RNA-Seq data. We show that OVERLAY can compute C. reinhardtii metabolism under nitrogen deprivation and metabolic shifts after an acetate boost. OVERLAY can also suggest possible 'bottleneck' proteins that need to be overexpressed to increase the TAG accumulation rate, as well as discuss other TAG-overproduction strategies.

Metabolic and genetic basis for auxotrophies in Gram-negative species.

Auxotrophies constrain the interactions of bacteria with their environment, but are often difficult to identify. Here, we develop an algorithm (AuxoFind) using genome-scale metabolic reconstruction to predict auxotrophies and apply it to a series of available genome sequences of over 1,300 Gram-negative strains. We identify 54 auxotrophs, along with the corresponding metabolic and genetic basis, using a pangenome approach, and highlight auxotrophies conferring a fitness advantage in vivo. We show that the metabolic basis of auxotrophy is species-dependent and varies with 1) pathway structure, 2) enzyme promiscuity, and 3) network redundancy. Various levels of complexity constitute the genetic basis, including 1) deleterious single-nucleotide polymorphisms (SNPs), in-frame indels, and deletions; 2) single/multigene deletion; and 3) movement of mobile genetic elements (including prophages) combined with genomic rearrangements. Fourteen out of 19 predictions agree with experimental evidence, with the remaining cases highlighting shortcomings of sequencing, assembly, annotation, and reconstruction that prevent predictions of auxotrophies. We thus develop a framework to identify the metabolic and genetic basis for auxotrophies in Gram-negatives.

Computation of condition-dependent proteome allocation reveals variability in the macro and micro nutrient requirements for growth.

Sustaining a robust metabolic network requires a balanced and fully functioning proteome. In addition to amino acids, many enzymes require cofactors (coenzymes and engrafted prosthetic groups) to function properly. Extensively validated resource allocation models, such as genome-scale models of metabolism and gene expression (ME-models), have the ability to compute an optimal proteome composition underlying a metabolic phenotype, including the provision of all required cofactors. Here we apply the ME-model for Escherichia coli K-12 MG1655 to computationally examine how environmental conditions change the proteome and its accompanying cofactor usage. We found that: (1) The cofactor requirements computed by the ME-model mostly agree with the standard biomass objective function used in models of metabolism alone (M-models); (2) ME-model computations reveal non-intuitive variability in cofactor use under different growth conditions; (3) An analysis of ME-model predicted protein use in aerobic and anaerobic conditions suggests an enrichment in the use of peroxyl scavenging acids in the proteins used to sustain aerobic growth; (4) The ME-model could describe how limitation in key protein components affect the metabolic state of E. coli. Genome-scale models have thus reached a level of sophistication where they reveal intricate properties of functional proteomes and how they support different E. coli lifestyles.

Adaptive evolution reveals a tradeoff between growth rate and oxidative stress during naphthoquinone-based aerobic respiration.

Evolution fine-tunes biological pathways to achieve a robust cellular physiology. Two and a half billion years ago, rapidly rising levels of oxygen as a byproduct of blooming cyanobacterial photosynthesis resulted in a redox upshift in microbial energetics. The appearance of higher-redox-potential respiratory quinone, ubiquinone (UQ), is believed to be an adaptive response to this environmental transition. However, the majority of bacterial species are still dependent on the ancient respiratory quinone, naphthoquinone (NQ). Gammaproteobacteria can biosynthesize both of these respiratory quinones, where UQ has been associated with aerobic lifestyle and NQ with anaerobic lifestyle. We engineered an obligate NQ-dependent γ-proteobacterium, Escherichia coli ΔubiC, and performed adaptive laboratory evolution to understand the selection against the use of NQ in an oxic environment and also the adaptation required to support the NQ-driven aerobic electron transport chain. A comparative systems-level analysis of pre- and postevolved NQ-dependent strains revealed a clear shift from fermentative to oxidative metabolism enabled by higher periplasmic superoxide defense. This metabolic shift was driven by the concerted activity of 3 transcriptional regulators (PdhR, RpoS, and Fur). Analysis of these findings using a genome-scale model suggested that resource allocation to reactive oxygen species (ROS) mitigation results in lower growth rates. These results provide a direct elucidation of a resource allocation tradeoff between growth rate and ROS mitigation costs associated with NQ usage under oxygen-replete condition.

Cellular responses to reactive oxygen species are predicted from molecular mechanisms.

Catalysis using iron-sulfur clusters and transition metals can be traced back to the last universal common ancestor. The damage to metalloproteins caused by reactive oxygen species (ROS) can prevent cell growth and survival when unmanaged, thus eliciting an essential stress response that is universal and fundamental in biology. Here we develop a computable multiscale description of the ROS stress response in Escherichia coli, called OxidizeME. We use OxidizeME to explain four key responses to oxidative stress: 1) ROS-induced auxotrophy for branched-chain, aromatic, and sulfurous amino acids; 2) nutrient-dependent sensitivity of growth rate to ROS; 3) ROS-specific differential gene expression separate from global growth-associated differential expression; and 4) coordinated expression of iron-sulfur cluster (ISC) and sulfur assimilation (SUF) systems for iron-sulfur cluster biosynthesis. These results show that we can now develop fundamental and quantitative genotype-phenotype relationships for stress responses on a genome-wide basis.

OxyR Is a Convergent Target for Mutations Acquired during Adaptation to Oxidative Stress-Prone Metabolic States.

Oxidative stress is concomitant with aerobic metabolism. Thus, bacterial genomes encode elaborate mechanisms to achieve redox homeostasis. Here we report that the peroxide-sensing transcription factor, oxyR, is a common mutational target using bacterial species belonging to two genera, Escherichia coli and Vibrio natriegens, in separate growth conditions implemented during laboratory evolution. The mutations clustered in the redox active site, dimer interface, and flexible redox loop of the protein. These mutations favor the oxidized conformation of OxyR that results in constitutive expression of the genes it regulates. Independent component analysis of the transcriptome revealed that the constitutive activity of OxyR reduces DNA damage from reactive oxygen species, as inferred from the activity of the SOS response regulator LexA. This adaptation to peroxide stress came at a cost of lower growth, as revealed by calculations of proteome allocation using genome-scale models of metabolism and macromolecular expression. Further, identification of similar sequence changes in natural isolates of E. coli indicates that adaptation to oxidative stress through genetic changes in oxyR can be a common occurrence.

Synthesizing Systems Biology Knowledge from Omics Using Genome-Scale Models.

Omic technologies have enabled the complete readout of the molecular state of a cell at different biological scales. In principle, the combination of multiple omic data types can provide an integrated view of the entire biological system. This integration requires appropriate models in a systems biology approach. Here, genome-scale models (GEMs) are focused upon as one computational systems biology approach for interpreting and integrating multi-omic data. GEMs convert the reactions (related to metabolism, transcription, and translation) that occur in an organism to a mathematical formulation that can be modeled using optimization principles. A variety of genome-scale modeling methods used to interpret multiple omic data types, including genomics, transcriptomics, proteomics, metabolomics, and meta-omics are reviewed. The ability to interpret omics in the context of biological systems has yielded important findings for human health, environmental biotechnology, bioenergy, and metabolic engineering. The authors find that concurrent with advancements in omic technologies, genome-scale modeling methods are also expanding to enable better interpretation of omic data. Therefore, continued synthesis of valuable knowledge, through the integration of omic data with GEMs, are expected.

Adaptations of Escherichia coli strains to oxidative stress are reflected in properties of their structural proteomes.

BACKGROUND: The reconstruction of metabolic networks and the three-dimensional coverage of protein structures have reached the genome-scale in the widely studied Escherichia coli K-12 MG1655 strain. The combination of the two leads to the formation of a structural systems biology framework, which we have used to analyze differences between the reactive oxygen species (ROS) sensitivity of the proteomes of sequenced strains of E. coli. As proteins are one of the main targets of oxidative damage, understanding how the genetic changes of different strains of a species relates to its oxidative environment can reveal hypotheses as to why these variations arise and suggest directions of future experimental work. RESULTS: Creating a reference structural proteome for E. coli allows us to comprehensively map genetic changes in 1764 different strains to their locations on 4118 3D protein structures. We use metabolic modeling to predict basal ROS production levels (ROStype) for 695 of these strains, finding that strains with both higher and lower basal levels tend to enrich their proteomes with antioxidative properties, and speculate as to why that is. We computationally assess a strain's sensitivity to an oxidative environment, based on known chemical mechanisms of oxidative damage to protein groups, defined by their localization and functionality. Two general groups - metalloproteins and periplasmic proteins - show enrichment of their antioxidative properties between the 695 strains with a predicted ROStype as well as 116 strains with an assigned pathotype. Specifically, proteins that a) utilize a molybdenum ion as a cofactor and b) are involved in the biogenesis of fimbriae show intriguing protective properties to resist oxidative damage. Overall, these findings indicate that a strain's sensitivity to oxidative damage can be elucidated from the structural proteome, though future experimental work is needed to validate our model assumptions and findings. CONCLUSION: We thus demonstrate that structural systems biology enables a proteome-wide, computational assessment of changes to atomic-level physicochemical properties and of oxidative damage mechanisms for multiple strains in a species. This integrative approach opens new avenues to study adaptation to a particular environment based on physiological properties predicted from sequence alone.

Visualizing metabolic network dynamics through time-series metabolomic data.

BACKGROUND: New technologies have given rise to an abundance of -omics data, particularly metabolomic data. The scale of these data introduces new challenges for the interpretation and extraction of knowledge, requiring the development of innovative computational visualization methodologies. Here, we present GEM-Vis, an original method for the visualization of time-course metabolomic data within the context of metabolic network maps. We demonstrate the utility of the GEM-Vis method by examining previously published data for two cellular systems-the human platelet and erythrocyte under cold storage for use in transfusion medicine. RESULTS: The results comprise two animated videos that allow for new insights into the metabolic state of both cell types. In the case study of the platelet metabolome during storage, the new visualization technique elucidates a nicotinamide accumulation that mirrors that of hypoxanthine and might, therefore, reflect similar pathway usage. This visual analysis provides a possible explanation for why the salvage reactions in purine metabolism exhibit lower activity during the first few days of the storage period. The second case study displays drastic changes in specific erythrocyte metabolite pools at different times during storage at different temperatures. CONCLUSIONS: The new visualization technique GEM-Vis introduced in this article constitutes a well-suitable approach for large-scale network exploration and advances hypothesis generation. This method can be applied to any system with data and a metabolic map to promote visualization and understand physiology at the network level. More broadly, we hope that our approach will provide the blueprints for new visualizations of other longitudinal -omics data types. The supplement includes a comprehensive user's guide and links to a series of tutorial videos that explain how to prepare model and data files, and how to use the software SBMLsimulator in combination with further tools to create similar animations as highlighted in the case studies.

Genome-scale model of metabolism and gene expression provides a multi-scale description of acid stress responses in Escherichia coli.

Response to acid stress is critical for Escherichia coli to successfully complete its life-cycle by passing through the stomach to colonize the digestive tract. To develop a fundamental understanding of this response, we established a molecular mechanistic description of acid stress mitigation responses in E. coli and integrated them with a genome-scale model of its metabolism and macromolecular expression (ME-model). We considered three known mechanisms of acid stress mitigation: 1) change in membrane lipid fatty acid composition, 2) change in periplasmic protein stability over external pH and periplasmic chaperone protection mechanisms, and 3) change in the activities of membrane proteins. After integrating these mechanisms into an established ME-model, we could simulate their responses in the context of other cellular processes. We validated these simulations using RNA sequencing data obtained from five E. coli strains grown under external pH ranging from 5.5 to 7.0. We found: i) that for the differentially expressed genes accounted for in the ME-model, 80% of the upregulated genes were correctly predicted by the ME-model, and ii) that these genes are mainly involved in translation processes (45% of genes), membrane proteins and related processes (18% of genes), amino acid metabolism (12% of genes), and cofactor and prosthetic group biosynthesis (8% of genes). We also demonstrated several intervention strategies on acid tolerance that can be simulated by the ME-model. We thus established a quantitative framework that describes, on a genome-scale, the acid stress mitigation response of E. coli that has both scientific and practical uses.

BOFdat: Generating biomass objective functions for genome-scale metabolic models from experimental data.

Genome-scale metabolic models (GEMs) are mathematically structured knowledge bases of metabolism that provide phenotypic predictions from genomic information. GEM-guided predictions of growth phenotypes rely on the accurate definition of a biomass objective function (BOF) that is designed to include key cellular biomass components such as the major macromolecules (DNA, RNA, proteins), lipids, coenzymes, inorganic ions and species-specific components. Despite its importance, no standardized computational platform is currently available to generate species-specific biomass objective functions in a data-driven, unbiased fashion. To fill this gap in the metabolic modeling software ecosystem, we implemented BOFdat, a Python package for the definition of a Biomass Objective Function from experimental data. BOFdat has a modular implementation that divides the BOF definition process into three independent modules defined here as steps: 1) the coefficients for major macromolecules are calculated, 2) coenzymes and inorganic ions are identified and their stoichiometric coefficients estimated, 3) the remaining species-specific metabolic biomass precursors are algorithmically extracted in an unbiased way from experimental data. We used BOFdat to reconstruct the BOF of the Escherichia coli model iML1515, a gold standard in the field. The BOF generated by BOFdat resulted in the most concordant biomass composition, growth rate, and gene essentiality prediction accuracy when compared to other methods. Installation instructions for BOFdat are available in the documentation and the source code is available on GitHub (https://github.com/jclachance/BOFdat).

Systematic discovery of uncharacterized transcription factors in Escherichia coli K-12 MG1655.

Transcriptional regulation enables cells to respond to environmental changes. Of the estimated 304 candidate transcription factors (TFs) in Escherichia coli K-12 MG1655, 185 have been experimentally identified, but ChIP methods have been used to fully characterize only a few dozen. Identifying these remaining TFs is key to improving our knowledge of the E. coli transcriptional regulatory network (TRN). Here, we developed an integrated workflow for the computational prediction and comprehensive experimental validation of TFs using a suite of genome-wide experiments. We applied this workflow to (i) identify 16 candidate TFs from over a hundred uncharacterized genes; (ii) capture a total of 255 DNA binding peaks for ten candidate TFs resulting in six high-confidence binding motifs; (iii) reconstruct the regulons of these ten TFs by determining gene expression changes upon deletion of each TF and (iv) identify the regulatory roles of three TFs (YiaJ, YdcI, and YeiE) as regulators of l-ascorbate utilization, proton transfer and acetate metabolism, and iron homeostasis under iron-limited conditions, respectively. Together, these results demonstrate how this workflow can be used to discover, characterize, and elucidate regulatory functions of uncharacterized TFs in parallel.

Thermosensitivity of growth is determined by chaperone-mediated proteome reallocation.

Maintenance of a properly folded proteome is critical for bacterial survival at notably different growth temperatures. Understanding the molecular basis of thermoadaptation has progressed in two main directions, the sequence and structural basis of protein thermostability and the mechanistic principles of protein quality control assisted by chaperones. Yet we do not fully understand how structural integrity of the entire proteome is maintained under stress and how it affects cellular fitness. To address this challenge, we reconstruct a genome-scale protein-folding network for Escherichia coli and formulate a computational model, FoldME, that provides statistical descriptions of multiscale cellular response consistent with many datasets. FoldME simulations show (i) that the chaperones act as a system when they respond to unfolding stress rather than achieving efficient folding of any single component of the proteome, (ii) how the proteome is globally balanced between chaperones for folding and the complex machinery synthesizing the proteins in response to perturbation, (iii) how this balancing determines growth rate dependence on temperature and is achieved through nonspecific regulation, and (iv) how thermal instability of the individual protein affects the overall functional state of the proteome. Overall, these results expand our view of cellular regulation, from targeted specific control mechanisms to global regulation through a web of nonspecific competing interactions that modulate the optimal reallocation of cellular resources. The methodology developed in this study enables genome-scale integration of environment-dependent protein properties and a proteome-wide study of cellular stress responses.

Global transcriptional regulatory network for Escherichia coli robustly connects gene expression to transcription factor activities.

Transcriptional regulatory networks (TRNs) have been studied intensely for >25 y. Yet, even for the Escherichia coli TRN-probably the best characterized TRN-several questions remain. Here, we address three questions: (i) How complete is our knowledge of the E. coli TRN; (ii) how well can we predict gene expression using this TRN; and (iii) how robust is our understanding of the TRN? First, we reconstructed a high-confidence TRN (hiTRN) consisting of 147 transcription factors (TFs) regulating 1,538 transcription units (TUs) encoding 1,764 genes. The 3,797 high-confidence regulatory interactions were collected from published, validated chromatin immunoprecipitation (ChIP) data and RegulonDB. For 21 different TF knockouts, up to 63% of the differentially expressed genes in the hiTRN were traced to the knocked-out TF through regulatory cascades. Second, we trained supervised machine learning algorithms to predict the expression of 1,364 TUs given TF activities using 441 samples. The algorithms accurately predicted condition-specific expression for 86% (1,174 of 1,364) of the TUs, while 193 TUs (14%) were predicted better than random TRNs. Third, we identified 10 regulatory modules whose definitions were robust against changes to the TRN or expression compendium. Using surrogate variable analysis, we also identified three unmodeled factors that systematically influenced gene expression. Our computational workflow comprehensively characterizes the predictive capabilities and systems-level functions of an organism's TRN from disparate data types.

COBRAme: A computational framework for genome-scale models of metabolism and gene expression.

Genome-scale models of metabolism and macromolecular expression (ME-models) explicitly compute the optimal proteome composition of a growing cell. ME-models expand upon the well-established genome-scale models of metabolism (M-models), and they enable a new fundamental understanding of cellular growth. ME-models have increased predictive capabilities and accuracy due to their inclusion of the biosynthetic costs for the machinery of life, but they come with a significant increase in model size and complexity. This challenge results in models which are both difficult to compute and challenging to understand conceptually. As a result, ME-models exist for only two organisms (Escherichia coli and Thermotoga maritima) and are still used by relatively few researchers. To address these challenges, we have developed a new software framework called COBRAme for building and simulating ME-models. It is coded in Python and built on COBRApy, a popular platform for using M-models. COBRAme streamlines computation and analysis of ME-models. It provides tools to simplify constructing and editing ME-models to enable ME-model reconstructions for new organisms. We used COBRAme to reconstruct a condensed E. coli ME-model called iJL1678b-ME. This reformulated model gives functionally identical solutions to previous E. coli ME-models while using 1/6 the number of free variables and solving in less than 10 minutes, a marked improvement over the 6 hour solve time of previous ME-model formulations. Errors in previous ME-models were also corrected leading to 52 additional genes that must be expressed in iJL1678b-ME to grow aerobically in glucose minimal in silico media. This manuscript outlines the architecture of COBRAme and demonstrates how ME-models can be created, modified, and shared most efficiently using the new software framework.

Quantitative time-course metabolomics in human red blood cells reveal the temperature dependence of human metabolic networks.

The temperature dependence of biological processes has been studied at the levels of individual biochemical reactions and organism physiology (e.g. basal metabolic rates) but has not been examined at the metabolic network level. Here, we used a systems biology approach to characterize the temperature dependence of the human red blood cell (RBC) metabolic network between 4 and 37 °C through absolutely quantified exo- and endometabolomics data. We used an Arrhenius-type model (Q10) to describe how the rate of a biochemical process changes with every 10 °C change in temperature. Multivariate statistical analysis of the metabolomics data revealed that the same metabolic network-level trends previously reported for RBCs at 4 °C were conserved but accelerated with increasing temperature. We calculated a median Q10 coefficient of 2.89 ± 1.03, within the expected range of 2-3 for biological processes, for 48 individual metabolite concentrations. We then integrated these metabolomics measurements into a cell-scale metabolic model to study pathway usage, calculating a median Q10 coefficient of 2.73 ± 0.75 for 35 reaction fluxes. The relative fluxes through glycolysis and nucleotide metabolism pathways were consistent across the studied temperature range despite the non-uniform distributions of Q10 coefficients of individual metabolites and reaction fluxes. Together, these results indicate that the rate of change of network-level responses to temperature differences in RBC metabolism is consistent between 4 and 37 °C. More broadly, we provide a baseline characterization of a biochemical network given no transcriptional or translational regulation that can be used to explore the temperature dependence of metabolism.

Biomarkers are used to predict quantitative metabolite concentration profiles in human red blood cells.

Deep-coverage metabolomic profiling has revealed a well-defined development of metabolic decay in human red blood cells (RBCs) under cold storage conditions. A set of extracellular biomarkers has been recently identified that reliably defines the qualitative state of the metabolic network throughout this metabolic decay process. Here, we extend the utility of these biomarkers by using them to quantitatively predict the concentrations of other metabolites in the red blood cell. We are able to accurately predict the concentration profile of 84 of the 91 (92%) measured metabolites (p < 0.05) in RBC metabolism using only measurements of these five biomarkers. The median of prediction errors (symmetric mean absolute percent error) across all metabolites was 13%. The ability to predict numerous metabolite concentrations from a simple set of biomarkers offers the potential for the development of a powerful workflow that could be used to evaluate the metabolic state of a biological system using a minimal set of measurements.

Systems biology definition of the core proteome of metabolism and expression is consistent with high-throughput data.

Finding the minimal set of gene functions needed to sustain life is of both fundamental and practical importance. Minimal gene lists have been proposed by using comparative genomics-based core proteome definitions. A definition of a core proteome that is supported by empirical data, is understood at the systems-level, and provides a basis for computing essential cell functions is lacking. Here, we use a systems biology-based genome-scale model of metabolism and expression to define a functional core proteome consisting of 356 gene products, accounting for 44% of the Escherichia coli proteome by mass based on proteomics data. This systems biology core proteome includes 212 genes not found in previous comparative genomics-based core proteome definitions, accounts for 65% of known essential genes in E. coli, and has 78% gene function overlap with minimal genomes (Buchnera aphidicola and Mycoplasma genitalium). Based on transcriptomics data across environmental and genetic backgrounds, the systems biology core proteome is significantly enriched in nondifferentially expressed genes and depleted in differentially expressed genes. Compared with the noncore, core gene expression levels are also similar across genetic backgrounds (two times higher Spearman rank correlation) and exhibit significantly more complex transcriptional and posttranscriptional regulatory features (40% more transcription start sites per gene, 22% longer 5'UTR). Thus, genome-scale systems biology approaches rigorously identify a functional core proteome needed to support growth. This framework, validated by using high-throughput datasets, facilitates a mechanistic understanding of systems-level core proteome function through in silico models; it de facto defines a paleome.

solveME: fast and reliable solution of nonlinear ME models.

BACKGROUND: Genome-scale models of metabolism and macromolecular expression (ME) significantly expand the scope and predictive capabilities of constraint-based modeling. ME models present considerable computational challenges: they are much (>30 times) larger than corresponding metabolic reconstructions (M models), are multiscale, and growth maximization is a nonlinear programming (NLP) problem, mainly due to macromolecule dilution constraints. RESULTS: Here, we address these computational challenges. We develop a fast and numerically reliable solution method for growth maximization in ME models using a quad-precision NLP solver (Quad MINOS). Our method was up to 45 % faster than binary search for six significant digits in growth rate. We also develop a fast, quad-precision flux variability analysis that is accelerated (up to 60× speedup) via solver warm-starts. Finally, we employ the tools developed to investigate growth-coupled succinate overproduction, accounting for proteome constraints. CONCLUSIONS: Just as genome-scale metabolic reconstructions have become an invaluable tool for computational and systems biologists, we anticipate that these fast and numerically reliable ME solution methods will accelerate the wide-spread adoption of ME models for researchers in these fields.