The recent progress in genome sequencing techniques has led to the development of genome-level models of metabolism that have been analyzed using constraint-based approaches, such as flux-balance analysis (FBA) [1, 2]. The success of FBA stems from the fact that, unlike kinetic models, FBA aims to identify optimal metabolic steady-state activity patterns that satisfy constraints imposed by mass balance, the metabolic network structure, and the availability of nutrients. The most common cellular task to be optimized (the system’s objective function) is that of growth, although other choices are possible depending on the selective environment of the cell [3, 4]. The FBA framework has been applied to many genome-level models (see e.g., [5–11]) with great success, as well as the systematic prediction of genetic knockout phenotypes [12, 13], the global organization of metabolic fluxes , and the discovery of novel regulatory interactions . However, fulfillment of systems biology’s goal to generate models that integrate data from all cellular levels (genomic, transcriptomic, proteomic, metabolomic, etc.), and can accurately predict metabolic phenomena under different environmental conditions has hitherto been hampered by minimal application of regulatory constraints.
According to the central dogma of biology, information flows from DNA to mRNA and ultimately to enzymes which catalyze and regulate various cellular functions. Hence, one might envision a fully “hierarchical” regulation of metabolism where expression levels of mRNA correlate directly with the amount of enzymes and thus with the flux through associated reactions. For some conditions, this simplified assumption can be used for the purpose of modeling metabolic activity . However, this type of hierarchical control does not take place in general since there are several levels of flux regulation which operate separately from the purely genetic. These mechanisms include variations in protein translation, protein activation/inactivation and metabolite regulation of enzymatic activity. Studies have shown that even within one pathway there may exist a variety of flux regulatory mechanisms for each reaction that range from purely hierarchical to fully metabolic control [17–20].
The varying role of hierarchical regulation for network reactions has limited the utilization of gene-expression data to improve predictions of genome-scale metabolic models. The earliest attempt at imposing transcriptional regulation on constraint-based models was conducted by Palsson and coworkers who developed regulatory flux-balance analysis (rFBA) [21–25] where, using Boolean logic, a transcriptional regulatory network was superimposed on an FBA model. rFBA can be used to predict a form of quasi-dynamic flux profile (i.e., series of steady-state flux profiles) in a changing environment. The time course of an experiment is divided into a number of successive short intervals and at each time step, new regulations based on metabolic steady state of the previous time is formulated. Next, FBA is used to predict a steady state flux that is consistent with the set regulatory rules at that moment.
Later, Nielsen and coworkers  further developed the idea of combined regulatory metabolic control by implementing gene-expression data as a Boolean switch to block the activity of any reaction for which the responsible mRNA was not expressed. Further progress on this methodology was made when Becker and Palsson  introduced the Gene Inactivity Moderated by Metabolism and Expression (GIMME) algorithm which uses a set of pre-determined thresholds for transition of each gene from “on” to “off”. The user selects a priori a minimally acceptable outcome for the FBA models and GIMME iteratively activates genes that were initially turned “off” in order to ensure that the FBA model achieves its required metabolic functionalities.
Another method dubbed E-Flux  uses gene-expression values to relatively regulate the flux that reactions in a model can carry. In a process akin to “setting the width of pipes” in a network, E-flux uses gene-expression data for different conditions to set normalized relative upper flux limits on affected reactions and then optimize a previously chosen objective function. Although the method is innovative in that it utilizes the actual gene-expression data, it is still limited in that a) it requires a pre-determined objective function for the condition associated with the gene-expression data, and b) the flux limit for each reaction is purely determined by the value of gene-expression values, and hence is unlikely to account for metabolic regulation. All subsequent advances involve utilizing mixed-integer linear programming (MILP) to identify cellular states that optimally adhere to both regulatory and metabolic regulations.
The introduction of steady-state regulatory flux balance analysis (SR-FBA)  which utilizes MILP to maximize biomass growth while concurrently trying to adhere to the maximum number of regulatory constraints, allowed a detailed quantification of the extent to which metabolic and transcriptional regulation control the metabolic behavior of a cell. Jensen and Papin further improved this mode of analysis by developing the Metabolic Adjustment by Differential Expression (MADE) methodology . This method, unlike GIMME, does not require a prior selection of expression thresholds and instead uses MILP and the statistical significance of changes in gene-expression to develop a metabolic model that recreates the measured expression dynamics while ensuring that the FBA model maintains previously determined threshold functionality. Although these methods have been useful in qualitatively predicting gene-expression patterns and metabolic adjustments between different conditions, they are limited by the fact that they require an a priori user-defined objective function and also do not fully make use of the predictions of FBA models; thus, a significant portion of the available data is not fully utilized.
Further work by Shlomi et al.  that has been incorporated in the iMAT algorithm  uses gene-expression data and a Boolean gene-to-reaction mapping to impose hierarchical regulation on a metabolic model. Here, affected reactions are classified based on associated gene-expression data as either highly expressed (RH), moderately expressed or lowly expressed (RL). iMAT utilizes MILP to identify a possible steady-state flux distribution among those that maximize the number of reactions with predicted flux consistent with the gene-expression data as well as the model’s stoichiometric and thermodynamic constraints. Thus, the goal of iMAT is to maximize the sum of the number of reactions in RL that carry a flux of zero, and the number of reactions in RH that carry a flux greater than an arbitrarily chosen threshold . Consequently, iMAT maximizes only the pattern of hierarchical regulation. Although the method has been successfully applied to model different human tissues (e.g., [33, 34]) and other multi-cellular organisms , the utility of the method is limited since ensuring that active reactions carry a minimum flux does not necessarily ensure that the model can predict correct cellular objective flux(es). Despite these deficiencies, iMAT has a strong advantage over other methods [26, 27, 29, 30], in that it does not need a predefined set of required metabolic functionalities and an FBA objective function.
Here, we present a new approach that uses gene-expression data to optimize not only the pattern of hierarchical regulation, but also the level of differential gene-expression within the rigid framework of metabolic constraints placed on a system by the connectivity of the reaction network. Although our steady-state based method does not account for capacity limitations in various enzymes, and thus the beneficial, deleterious or regulatory role of metabolite concentrations, the model's adherence to conservation of mass balance and network connectivity imposes pseudo-metabolic regulation. The coupled interaction of this absolute form of metabolic control with optimal hierarchical control in gene-expression FBA (GX-FBA) improves our theoretical capabilities for analyses of a wide range of phenomena, such as cellular responses to environmental perturbations which traditionally have been considered outside the realm of FBA.
A recently published approach by Lee et al.  is also focused on using actual gene-expression levels to guide metabolic flux prediction. However, this method differs significantly from GX-FBA in that it minimizes the absolute difference between metabolic fluxes and gene-expression data from RNA-seq experiment.
To illustrate the utility of GX-FBA, we have analyzed the genome-scale metabolic model for the etiological agent of bubonic plague, the gram-negative bacterium Yersinia pestis (YP) . We have studied YP’s genome-scale metabolic response in physiologically important conditions: temperature shifts known to induce virulence in low calcium media [38, 39], as well as its response to stress induced separately by the antibiotics streptomycin and chloramphenicol. Our analyses open windows into the metabolic workings of this bacterium while it survives within macrophages following initial introduction into a mammalian host, proliferates in the blood, and attempts to resist therapeutic efforts. Our analyses indicate that majority of cellular metabolic changes associated with response to stress is unique to the type of perturbation. The only common adaptive response to all four types of stress was for YP to initiate a series of energy saving measures.