Abstract
Abnormal metabolism is a hallmark of cancer, yet its regulation remains poorly understood. Cancer cells were considered to utilize primarily glycolysis for ATP production, referred to as the Warburg effect. However, recent evidence suggests that oxidative phosphorylation (OXPHOS) plays a crucial role during cancer progression. Here we utilized a systems biology approach to decipher the regulatory principle of glycolysis and OXPHOS. Integrating information from literature, we constructed a regulatory network of genes and metabolites, from which we extracted a core circuit containing HIF-1, AMPK, and ROS. Our circuit analysis showed that while normal cells have an oxidative state and a glycolytic state, cancer cells can access a hybrid state with both metabolic modes coexisting. This was due to higher ROS production and/or oncogene activation, such as RAS, MYC, and c-SRC. Guided by the model, we developed two signatures consisting of AMPK and HIF-1 downstream genes, respectively, to quantify the activity of glycolysis and OXPHOS. By applying the AMPK and HIF-1 signatures to The Cancer Genome Atlas patient transcriptomics data of multiple cancer types and single-cell RNA-seq data of lung adenocarcinoma, we confirmed an anticorrelation between AMPK and HIF-1 activities and the association of metabolic states with oncogenes. We propose that the hybrid phenotype contributes to metabolic plasticity, allowing cancer cells to adapt to various microenvironments. Using model simulations, our theoretical framework of metabolism can serve as a platform to decode cancer metabolic plasticity and design cancer therapies targeting metabolism. Cancer Res; 77(7); 1564–74. ©2017 AACR.
Major Findings
We developed a theoretical framework for modeling gene regulation of cancer metabolism. We found that, in addition to glycolytic and oxidative metabolism in normal cells, cancer cells have a new hybrid phenotype in which both metabolic modes coexist. Cells in the hybrid state have enhanced metabolic plasticity to facilitate tumorigenesis and metastasis. Guided by the modeling, we developed the AMPK and HIF-1 signatures to quantify the activity of OXPHOS and glycolysis, respectively. We propose a new cancer therapeutic strategy by targeting the hybrid state.
Introduction
Cells can utilize multiple metabolic pathways for energy production and biosynthesis depending on the requirements for cellular function and the availability of metabolites. In the presence of glucose, cells typically uptake glucose and convert it to pyruvate inside the cytosol by glycolysis. Under normoxic conditions, pyruvate is further transported into mitochondria where it undergoes oxidative phosphorylation (OXPHOS) to produce ATP through the tricarboxylic acid (TCA) cycle and the electron transport chain (ETC). On the other hand, under hypoxia, cells utilize anaerobic glycolysis instead, which converts pyruvate into lactate, and produces ATP in a much faster but less efficient way. If fatty acid or glutamine is available, cells under normoxia can also undergo fatty acid oxidation (also called β-oxidation) or glutamine oxidation.
Quick Guide to Equations and Assumptions
In this study, we constructed a regulatory network of cell metabolism that consists of regulatory proteins and metabolites. To understand the behavior of the network, we coarse-grained the whole network into a minimalist circuit of three components —AMPK, HIF-1, and ROS. Here we assumed that the essential feature of the network could be captured by this circuit, while the effects of many other relevant genes can be regarded as inputs to this circuit. The dynamic behaviors of the regulatory circuit AMPK:HIF-1:ROS (Fig. 2) were analyzed by solving nonlinear differential rate equations. Typically, the deterministic rate equation for a protein or metabolite has a generic form:
Here, 
represents the level of the protein or the metabolite. 
and 
represent the basal production and degradation rates of 
, respectively. 
and 
represent the regulation on the production and degradation of 
, respectively.
Overall, we assumed that the regulation of 
by a component 
can be described by a nonlinear function, namely the shifted Hill function (1, 2),
where 
is the threshold level, 
is the Hill coefficient, and 
is the fold change.
is the inhibitory Hill function, defined as
is the excitatory Hill function, defined as
represents an inhibitory regulation and 
represents an excitatory regulation.
When the production of 
is regulated by 
. When the degradation of 
is regulated by 
.
When the production (or degradation) of 
is regulated by two components 
and 
simultaneously, 
can be expressed as:
Details of the functional form of 
can be found in Supplementary Information Section S3.
In our modeling, AMPK and HIF-1 competitively regulate the production of ROS (both mitochondria ROS and cytosol ROS); ROS and oxygen competitively regulate the degradation of HIF-1. All the other combined regulatory processes are assumed to be independent. The parameters were estimated per experimental evidence, as shown in Supplementary Tables S3–S6. Furthermore, the model allowed evaluation of the effects of oncogenes and metabolic therapies. Details can be found in Supplementary Information Section S3–S6, S9.
Abnormal metabolism is a hallmark of cancer (3–5). Unlike normal cells, cancer cells largely depend on glycolysis to produce energy even in the presence of oxygen, which is referred to as the Warburg effect (6, 7) or aerobic glycolysis. Aerobic glycolysis is an aggressive metabolic phenotype in that it has the advantage to produce ATP at a high rate and prepare biomass for amino acids and fatty acid synthesis that are required for rapid cell proliferation (7). Although aerobic glycolysis is generally regarded as a dominant metabolism in cancer cells (8), recent evidence suggests that mitochondrial OXPHOS is also utilized by cancer cells (9–12). Indeed, the influence of mitochondria on cancer cells is well documented (13–19), suggesting that mitochondria are actively functioning in cancer cells. More importantly, it has been revealed that OXPHOS contributes to cancer metastasis (11–15, 19–21). Yet, it is still unclear how these metabolic modes are regulated in cancer cells, particularly the advantages of having each one of them.
To shed light on these puzzles, we utilize a computational systems biology approach to study the genetic regulation of glycolysis and OXPHOS (including glucose oxidation, fatty acid oxidation, and glutamine oxidation). Certainly, cancer metabolism has been studied by several computational modeling approaches. However, each of these earlier works mostly focused on only a part of the interplay between glycolysis and OXPHOS. For example, there are studies of the Warburg effect, including modeling of the regulation by reactive oxygen species (ROS) of hypoxia-inducible factor 1 (HIF-1) in response to hypoxia (22), and modeling of the genome-scale metabolic network based on metabolic flux balance analysis (23). Some other studies focused on mitochondrial functions, for example, in the study of the cooperative effects of TrxR and Prx in defending oxidative stresses (24). The evolutionary advantage of metabolic heterogeneity of cancer has also been investigated (25) with emphasis on glycolysis and acidity. However, we still lack a comprehensive picture of how different metabolic pathways contribute to the survival and oncogenic potential of cancer cells. Therefore, there is an urgent need to explore in detail how genes and metabolites regulate both aerobic glycolysis and OXPHOS.
In this study, we establish a theoretical framework for modeling genetic regulation of the interplay between glycolysis and OXPHOS in cancer metabolism. On the basis of an extensive literature survey, we construct a metabolic network model for both glycolysis and OXPHOS. The network includes both regulatory genes and metabolites that are involved in these metabolic pathways. To capture the two regimes of cancer metabolism, we further coarse-grained the network into a core regulatory circuit, composed of HIF-1, 5′ AMP-activated protein kinase (AMPK) and ROS. The circuit contains HIF-1 and AMPK, as they are the master regulators of glycolysis and OXPHOS, respectively (7); and ROS is included because both cytosolic and mitochondrial ROS (mtROS) mediate the interplay between the regulatory genes and metabolites (20, 26). This reduced representation provides a clear picture of the interplay between glycolysis and OXPHOS. Similar HIF-1:AMPK:ROS representation has also been proposed in the study of longevity in Caenorhabditis elegans (C. elegans; ref. 26), suggesting that the circuit is evolutionarily robust. The core regulatory circuit is modulated by several oncogenic signaling pathways, including MYC, RAS, and c-SRC (27–29). Many metabolic drugs (30) directly target the genes that are associated with specific metabolic reactions, such as GLUT1, HK and FASN, or regulatory genes in the current study, that is, AMPK and HIF-1. Therefore, the regulatory circuit allows us to model the effects of the oncogenes and metabolic drugs on cellular metabolism. As we show later, this model can explain important experimental observations, and is sufficient to capture the main differences between cancer and normal cells.
Computational modeling of our reduced circuit shows that in general cancer cells can have three stable states—a Warburg state (W: high HIF-1, low pAMPK), an oxidative state (O: low HIF-1, high pAMPK), and a hybrid state (W/O: high HIF-1, high pAMPK). Here, the AMPK activity is represented by the level of phosphorylated AMPK (pAMPK) at threonine-172 of the α subunit. Compared with cancer cells, normal cells usually lack the hybrid phenotype due to their lower mtROS production and higher rate of HIF-1 degradation. The activation of oncogenic pathways is an alternative way to induce cells to lie in the hybrid state. The discovery of the hybrid metabolic state for cancer cells from the network analysis has been supported by recent experimental evidence (summarized in Supplementary Table S1; refs. 10, 11, 16, 19, 20). The association of cancer metabolism and oncogenic pathways are also supported by extensive The Cancer Genome Atlas (TCGA) and single-cell transcriptomics data analysis of multiple cancer types under this theoretical framework. Particularly, we develop two signatures to quantify AMPK and HIF-1 activities and find that AMPK activity anticorrelates with HIF-1 activity in several tumors but not in the others. We propose that the W/O hybrid state enhances metabolic plasticity, therefore allowing the cancer cells to adapt to the tumor microenvironment and promote cancer progression and metastasis. Guided by the simulations of targeting AMPK or HIF-1, we make predictions on the effectiveness of various therapeutic strategies, such as hyperbaric oxygen, metformin, and 3-bromopyrvate (3BP), in reducing the metabolic plasticity of cancer cells. All told, our model provides not only a theoretical framework for charactering the regulatory mechanism of glycolysis and OXPHOS and its connection to tumorigenesis, but also new insights into the optimal therapeutic strategies that specifically target cancer metabolism.
Results
The regulatory network of glycolysis and OXPHOS
After an extensive literature analysis, we constructed a comprehensive network (Fig. 1) featuring the regulation of oxidative respiration and glycolysis by both genes and metabolites. The network contains the following four types of regulatory interactions. First, the three major metabolic pathways, aerobic glycolysis, glucose oxidation and fatty acid oxidation (green rectangles), directly inhibit each other because they compete for shared metabolites (see Supplementary Information Section S1 for details). Second, the activity of these metabolic pathways is directly regulated by specific genes (red ovals). Most prominently, the fatty acid oxidation and glucose oxidation are regulated by an energy sensor gene, AMPK and the glycolytic pathway is regulated by the hypoxia-inducible factor, HIF-1. Third, some metabolites produced by the metabolic pathways (yellow rectangles) in turn regulate the activities of some regulatory genes. HIF-1 is stabilized by both NOX-derived ROS (noxROS) from cytosol and OXPHOS-derived ROS from mitochondria (mtROS). Meanwhile, ROS induces phosphorylation and activation of AMPK, yet excessive production of ATP by either metabolic pathway could block the AMPK activity (see Supplementary Table S1 for details). Fourth, the regulatory genes are also coupled to several oncogenic pathways (blue ovals), including RAS, MYC, and c-SRC (Supplementary Table S2).
Schematic illustration of the regulatory network of metabolism. The network includes both regulatory genes (ellipses) and metabolites (rectangles). The arrows represent positive regulation and the bars represent negative regulation. AMPK and HIF-1 play a central role in regulating different metabolic pathways, while the metabolic pathways regulate activities of AMPK and HIF-1 partially through ROS. The oncogenic pathways directly regulate the activity of AMPK, HIF-1, and noxROS.
The core AMPK:HIF-1:ROS regulatory circuit
To capture the principles of how genes and metabolites modulate metabolism, we coarse-grained the extensive metabolic network (Fig. 1) to a minimalist regulatory circuit consisting of AMPK, HIF-1, and ROS (including noxROS and mtROS; Fig. 2). These components were chosen, because they play critical roles in regulating the decision-making of both glycolysis and OXPHOS. Moreover, the core regulatory circuit captures the main features of the more complete network. As we show later, the reduced circuit is sufficient to explain important experimental observations of glycolysis and OXPHOS. It is worth noting that the core circuit is directly coupled with both the oncogenic pathways and metabolic pathways. Its dynamic behavior can shed insight into the interplay of the various metabolic modes in cancer.
The AMPK:HIF-1:ROS regulatory circuit. AMPK and HIF-1 are the master regulators of OXPHOS and glycolysis, respectively. ROS represents both mtROS and noxROS. RAS, MYC, and c-SRC modulates the balance of glycolysis and OXPHOS (details in Fig. 4).
As shown in Fig. 2, the regulatory links among these four components, AMPK, HIF-1, noxROS, and mtROS, are derived from either direct or indirect regulatory interactions (see Supplementary Information Section S2 for details). Various oncogenic pathways are modeled as input signals to the circuit, and the activities of OXPHOS and glycolysis are modeled as the readout of the circuit. The detailed experimental evidence for the regulatory interactions of the full network and the associated reduced circuit are listed in the Supplementary Table S2. The corresponding chemical rate equations are provided in Supplementary Information Section S3. The parameters and their experimental justifications are provided in Supplementary Tables S3–S6.
The metabolic circuit allows multiple cell phenotypes
To identify the possible cell states that are allowed by the core circuit, we performed an analysis of the circuit equations with two sets of parameter that correspond to normal and cancer cells, respectively. Here we illustrate the results by considering the differences between cancer cells and normal cells in two aspects. First, cancer cells have higher mtROS production due to reprogrammed mitochondria (20). Second, cancer cells have more stable HIF-1 because less oxygen is available to each cell due to abnormally rapid proliferation. The effects of the oncogenic pathways will be evaluated in the next section. Therefore, as compared with normal cells, cancer cells produce more mtROS in response to the AMPK activation (represented by 
, maximum fold change of mtROS by the AMPK activation), and have a lower HIF-1 degradation rate (denoted by 
). As mentioned above, the activity of AMPK was quantified by the concentration of its phosphorylated form, that is, pAMPK.
Our simulations show that normal cells have two stable steady states, which correspond to the high HIF-1, low pAMPK state and the low HIF-1, high pAMPK state (Fig. 3A), which are associated with the Warburg effect (W) and oxidative respiration (O) respectively (the “W” state and “O” state in Fig. 3A). This result is consistent with the fact that while cells usually use glucose oxidation to produce energy, they switch to glycolysis during anaerobic exercise. Next, we performed an analogous analysis for the cancer cells, as reflected by larger 
and lower 
. Interestingly, we found that cancer cells have a new hybrid state – (high pAMPK, high HIF-1; “W/O” state in Fig. 3B) in addition to the “W” and “O” states. The new hybrid state is found to be a robust feature of the cancer cells by parameter sensitivity analysis (Supplementary Fig. S1). In the hybrid state, the levels of pAMPK and HIF-1 are both high because the cancer cells maintain high levels of ROS that both stabilizes HIF-1 and activates AMPK. The existence of the hybrid metabolic phenotype means that cancer cells have the flexibility to simultaneously utilize both glycolysis and mitochondrial OXPHOS.
The nullclines and steady states in the phase space of AMPK and HIF-1. The red line represents the nullcline of 
, and the blue line represents the nullcline of 
(see Supplementary Information Eq. S1). The green solid dots denote stable steady states and the green hallow dots denote unstable steady state. Each stable state is associated with a metabolic phenotype. For normal cells (A), the circuit allows a Warburg state, denoted by “W,” and an OXPHOS state, denoted by “O.” For cancer cells (B), the circuit allows an additional hybrid state, denoted as “W/O.”
The role of mtROS production and ontogenetic pathways in cancer metabolism
We further evaluated the effects of mtROS production and several oncogenes on the dynamic behavior of the AMPK:HIF-1:ROS circuit. In Fig. 4A, we show that, when the rest of the parameters remain unchanged, the increase of 
shifts the cells from bistability (with two stable steady states) to tristability (with three stable steady states) and gives rise to the hybrid metabolic state – “W/O,” in addition to “W” and “O” phenotypes. This is because, once again, higher production of mtROS by mitochondria stimulates both AMPK activation and HIF-1 stabilization, which is consistent with the experimental results that partial inhibition of ETC leads to the coexistence of glycolysis and OXPHOS by producing excessive mtROS (20). In addition to ROS production, stabilization of HIF-1, as reflected by lower kh, can also give rise to the hybrid W/O state (Supplementary Fig. S2).
Bifurcation diagrams of the levels of phosphorylated AMPK (pAMPK) and HIF-1 driven by mtROS production (A), MYC (B), c-SRC (C), and RAS (D). 
represents the fold change of mtROS by AMPK activation. The blue, green, and magenta solid lines represent the phenotypes O, W/O, and W, respectively. The red dashed line represents the unstable states. Different background colors represent different phases.
Similarly, we simulated the response of the core circuit to the changes in the activity of the following oncogenic pathways, that is, MYC (Fig. 4B), c-SRC (Fig. 4C), and RAS (Fig. 4D; See Supplementary Information Section S4–S6 for the modeling details). Interestingly, in all these cases, we found that moderate activation of RAS, MYC, or c-SRC induces the circuit from a two-state system to a three-state system containing the hybrid state “W/O.” Strong activation of RAS can drive the circuit to either the “W” or the “W/O” phenotype, while strong activation of c-SRC can drive the circuit to either the “O” or the “W/O” phenotype. However, strong activation of MYC can drive the circuit to the “W” state. The results are consistent with the experimental evidences on the role of some oncogenic pathways in modulating both glycolysis and OXPHOS. Previously, we show that cells with high c-SRC activity have substantial levels of glycolysis and fatty acid oxidation (31). In addition, MYC overexpression has been shown to promote fatty acid oxidation in triple-negative (TN) breast cancer (32). Moreover, RAS can simultaneously activate glycolysis and OXPHOS (33). It is likely that cancer cells can be in the hybrid phenotype when the effects of the excessive production of ROS or the activation of the oncogenes outweigh the mutual inhibitions between AMPK and HIF-1.
Quantification of metabolic state by the activities of AMPK and HIF-1
To validate the metabolic circuit model, we performed an extensive data analysis on the TCGA data from eight cancer types. Here we quantified the activities of AMPK and HIF-1 by evaluating the expression of downstream targets of both AMPK and HIF-1 (a total of 33 genes for AMPK, and 23 genes for HIF-1, see Supplementary Information Section S7 for details). For each cancer type, we performed principal component analysis (PCA) on the RNA-seq data independently for either the AMPK downstream genes or the HIF-1 downstream genes, from which we assigned the first principal component as the axis to quantify the activities of AMPK or HIF-1. The AMPK or the HIF-1 axes (signatures) derived from different cancer types are surprisingly similar (Supplementary Table S7 and S8; Supplementary Fig. S3), indicating the consistent regulatory functions of AMPK and HIF-1 across these tissues and cancer types. Most of the genes have positive contributions to the respective principal component (Supplementary Fig. S4), indicating that they are positively regulated by AMPK or HIF-1. For a few genes that have substantial negative contributions, we identified experimental evidence for the negative regulation (Supplementary Table S9).
Strikingly, even though the AMPK and HIF-1 gene sets are independent and the corresponding principal axes were obtained independently, we observed strong anticorrelations between the HIF-1 and AMPK activities for several major cancer types—liver hepatocellular carcinoma (HCC), lung adenocarcinoma (LUAD), breast invasive carcinoma (Fig. 5A–C). We obtained similar results for three other types of cancer—stomach adenocarcinoma, acute myeloid leukemia, and pancreatic adenocarcinoma (PAC; Supplementary Fig. S5A–S5C), but not in clear cell renal cell carcinoma (CCRCC), prostate adenocarcinoma, and colorectal adenocarcinoma (Supplementary Fig. S5D–S5F). From the data analysis, we observed a spectrum of cases with “W,” “W/O” and “O” phenotypes. Notably, we did not observe two clearly separated clusters representing glycolysis and OXPHOS respectively; therefore, the AMPK and HIF-1 activities are strongly anticorrelated but not mutually exclusive in cancer. However, the TCGA data analysis was not sufficient to explain the existence of the hybrid “W/O” state due to the mixture of normal and cancer cells in the tissues and possible heterogeneity of metabolic states in the cell populations. Thus, we looked into single-cell RNA-seq data of LUAD (34), which shows a similar anticorrelation between the AMPK and HIF-1 activities and mainly a bimodal distribution with substantial amount of hybrid cells (Fig. 5D; see Supplementary Figs. S8–S11 for more details on the single-cell RNA-seq data analysis). This finding suggests that the hybrid metabolic phenotype can also be observed for individual cancer cell.
Evaluation of the AMPK and HIF-1 activities using TCGA patient data and single-cell RNA-seq data of lung adenocarcinoma in HCC (n = 373; A), LUAD (n = 517; B), and breast-invasive carcinoma (n = 971; C). Left, each point represents the AMPK and HIF-1 activities of one sample. For each dataset, the standard k-mean analysis was applied to group the cases into the W, W/O, and O states. Middle, the fraction of the top 20% samples with high oncogene activities in each metabolic group (W, W/O, and O). χ2 test was used to test the significance (See Supplementary Information Section S10). *, P < 0.05; **, P < 0.01; ***, P < 0.001. Right, Kaplan–Meier (KM) overall survival curves of patients that were stratified by their metabolic states. D, Evaluating the AMPK and HIF-1 activities of single LUAD cells (n = 77). Left, each point represents the AMPK and HIF-1 activities of a single cell. Histogram represents the distribution of single cells projected to the first principal component of their AMPK and HIF-1 signatures. Three Gaussians are required to fit to the histogram (Supplementary Fig. S11). Right, box plot for HIF-1 activity. Evaluation of the AMPK and HIF-1 activities of luminal A (n = 531), luminal B (n = 135), HER2+ (n = 36), and TN (n = 141) breast carcinoma (E), and HCC (n = 373), CCRCC (n = 534), and PAC (n = 179) from TCGA (F; see Supplementary Information Section S7 for details of methods).
Next, we applied the AMPK and HIF-1 signatures to evaluate the metabolic activities of different subtypes of breast carcinoma and different types of tumors across TCGA. HER2+ and TN breast tumors clearly show higher HIF-1 activity and lower AMPK activity than luminal tumors (Fig. 5E), consistent with the experimental observation that TN cells are more dependent on glycolysis than luminal cells (35). PAC and CCRCC show higher HIF-1 activity and lower AMPK activity compared with HCC (Fig. 5F; Supplementary Fig. S6), consistent with experimental observations (see Supplementary Fig. S6 for details).
For each of the three cancer types shown in Fig. 5A–C, we further evaluated the activities of the oncogenic pathways in samples with various metabolic phenotypes. First, each cancer sample was assigned a metabolic state (“W,” “W/O” and “O”) by k-mean clustering analysis of the whole dataset (Fig. 5, left). The activities of the oncogenic pathways (MYC, c-SRC, and AKT) were quantified by the oncogenic scores derived experimentally in previous studies [see Supplementary Information Section S8 for details regarding the calculation of the oncogene score (36)]. The classification and scoring methods allow us to evaluate the enrichment of samples with high oncogenic activity for each metabolic phenotype (Fig. 5A–C). Interestingly, high MYC activity was observed to be enriched in the “W” phenotype, whereas high c-SRC and AKT (representing RAS) activities are enriched in both the “W” and the “W/O” phenotypes. Low oncogenic activity is associated with the “O” phenotype and low RS group (better prognosis; ref. 34) of single LUAD cells (Supplementary Figs. S7 and S8). These results are consistent with the modeling results and experimental evidence discussed in the previous section. Moreover, we evaluated the 5-year overall survival rates for different metabolism phenotypes and found that patient samples in “O” phenotype are significantly correlated with better survival results compared with patients' samples in either “W” or “W/O” (Fig. 5A–C). Consistently, single LUAD cells with higher risk scores (worse prognosis; ref. 34) maintain higher HIF-1 activity (Fig. 5D). Thus, we believe that our HIF-1 and AMPK signatures provide an easy and reliable scoring method to quantify the metabolic states of glycolysis and OXPHOS based on RNA-seq data from both bulk tumors and single cells.
The advantages of the hybrid metabolic phenotype and experimental support
We argue that cancer cells in the “W/O” hybrid state can have a significant advantage in supporting their survival, proliferation, and metastasis. First, cells in the hybrid state have higher metabolic plasticity to better adapt to various microenvironments, such as hypoxia and acidic conditions, because these cells are flexible in their use of available nutrients, such as glucose and fatty acids, to produce energy. Second, the hybrid state allows cancer cells to efficiently produce energy by both OXPHOS and glycolysis while at the same time using lactate and pyruvate, the byproducts of glycolysis, to generate biomass for cell proliferation (7). Third, cells in the hybrid state can modulate ROS at a moderate level, so that cells can take advantage of ROS signaling (20, 37) to promote metastasis meanwhile avoid excessive DNA damage (38). Here, ROS scavenging is achieved by pyruvate and NADPH, the byproducts of glycolysis and fatty acid oxidation, respectively. Fourth, the hybrid state might be specifically associated with metastasis, as suggested by some experimental evidence listed below. In summary, cancer cells at the hybrid metabolic phenotype have advantage in multiple aspects over cells with either the “W” or the “O” phenotypes.
Our metabolic circuit model is consistent with many recent observations that, although glycolysis is a typical feature of cancer cells, mitochondrial respiration plays a crucial role in tumor invasion and metastasis (9, 16–21, 39). It has been observed that some aggressive tumor cells, such as SiHa and HeLa, have not only oxidative respiration but also glycolysis due to stronger HIF-1 activation from lactate (21). Moreover, the inhibition of the mitochondrial respiratory chain contributes to the reduction of multidrug resistance of slow cycling melanoma cells (14). It has also been shown that the surviving pancreatic cancer cells after doxycycline withdrawal depend on OXPHOS and are highly sensitive to OXPHOS inhibitors (13). Also, the super-metastatic tumor cells obtained by experimental selection in vitro (SiHa-F3 cells) and in vivo (B16F10 and B16-M1 to M5 tumor cells) have increased OXPHOS with higher ROS production (20). TGFβ1 treatment of non–small lung carcinoma A459 cells induces metastasis, and they were found to have decreased fatty acid synthesis and increased oxygen consumption (11). The inhibition of OXPHOS by graphene for multiple cancer cells, including breast cancer cells and hepatocellular carcinoma cells, can effectively inhibit tumor migration and invasion (39). In the context of 4T1 mammary epithelial cancer cells, PGC-1α promotes metastasis through activating mitochondrial biogenesis and OXPHOS (10). The breast cancer cells can also shift from glycolysis to mitochondrial OXPHOS after radiation exposure to generate more ATP for survival (16). A detailed list of these experimental findings is presented in Supplementary Table S1. These data taken together indicates a critical role of OXPHOS in tumorigenesis, and supports our prediction of the hybrid “W/O” state. Our model and its experimental support suggest that the hybrid metabolic phenotype could be a good target for metabolic drugs.
Modeling therapeutic strategies targeting cancer metabolism
It has been shown that targeting cellular metabolism is a promising strategy for fighting against cancer (details in Supplementary Table S10; ref. 30). Certain metabolic drugs have been shown to be effective in treating cancers in some cases. These drugs, for example, S-trans, trans-farnesylthiosalicylic acid (FTS; ref. 40), 2-deoxy-d-glucose (2DG; ref. 41), 3-bromopyruvate (3BP; ref. 42), metformin (43), and AICAR (44), typically target glycolysis or mitochondrial OXPHOS. It has also been suggested that combinations of different drugs could prove more effective (40, 45, 46). For example, administrating metformin and 2DG together can induce massive ATP depletion in cancer cells and further trigger cellular processes to induce cell death, such as autophagy or p53-dependent apoptosis (45). Often, metabolic drugs work well in some cases but not in others, and the underlying mechanisms are not completely understood. Here we propose that these therapeutic strategies might be effective partly because they target the hybrid metabolic phenotype of cancer cells. An effective treatment outcome could be achieved by shifting the metabolic phenotype from the hybrid state to the other allowed states of the cancer cells. By doing so, the drugs might sensitize the cancer cells, therefore, other therapies, such as chemotherapy, can be more effective in killing the cancer cells.
We evaluated the dynamic response of the AMPK:HIF-1:ROS circuit to various treatment strategies, each of which exerts its own regulatory mechanism on the circuit (See Supplementary Information Section S9 for details). In particular, hyperbaric oxygen therapy can effectively reduce hypoxia and accelerate degradation of HIF-1; 3BP, similar to 2DG, targets the glycolytic enzymes such as glucose transporters (GLUT) and hexokinase (HK), thus effectively reducing cellular glycolysis; metformin activates AMPK, inhibits ETC Complex-1, and inhibits mTOR in an AMPK-independent manner, which further inhibits HIF-1. The effective circuit diagrams for each case are shown in Fig. 6, and details of the modeling procedures are shown in Supplementary Information Section S5 with the relevant parameters given in Supplementary Table S11. The model allows us to calculate the steady states and simulate the time course of the levels of pAMPK, HIF-1, and ROS. We are especially interested in how different treatments drive phenotypic transitions among different metabolic states.
The bifurcation diagrams of the levels of pAMPK in response to the oxygen level in the hyperbaric oxygen therapy (A), 3BP level (B), metformin level (C), and combined 3BP and metformin therapy (D). In A, the normal condition of oxygen level is 5%. In B–D, the drug (1/2/3) represents a hypothetic signal with opposite effect of the corresponding drug(s). In A–D, blue, green, and magenta solid lines represent phenotype O, W/O, and W, respectively. The red dashed lines represent the unstable steady states. The dashed lines with arrows indicate transitions between two phenotypes. Different background colors represent different phases.
A comparison of possible metabolic therapies
Starting from the model where cancer cells have the “W,” the “W/O,” and the “O” states, we evaluated the effects of various treatment strategies on allowing the cell to escape from the tristability of the metabolic circuit. Each of the panel of Fig. 6A–D shows the bifurcation diagram of the levels of pAMPK and HIF-1 with respect to the level of the drugs, for each treatment. The results show that these treatment strategies can shift the cancer metabolism from the tristable phase to the monostable phase (“O”) by increasing the dose level.
From the bifurcation analysis, metformin is less effective compared with the hyperbaric oxygen therapy and 3BP. The inefficiency of metformin is caused by the decrease of mitochondrial potential and increase of mtROS production by the drug (46). Various treatment strategies were further evaluated by simulating the time course of the levels of pAMPK, HIF-1, and ROS (Supplementary Figs. S12–S15). We assume that an effective therapy needs to drive cancer cells away from the hybrid states. We found that the administration of metformin alone is less effective in avoiding the hybrid state, while the combined therapies are more effective. The prediction is consistent with the finding that, in human gastric and esophageal cell lines, the administration of both metformin and 2DG (whose effect is similar to 3BP) is more effective than the administration of each drug alone (46).
Discussion
In this study, we established a theoretical framework for modeling genetic regulation of cancer metabolism. By integrating existent data, a network was constructed that features the regulation of mitochondrial OXPHOS and glycolysis by both regulatory protein and metabolites. We further coarse-grained the network into a core regulatory circuit that is composed of HIF-1, AMPK, and ROS. Although we made several simplifications, the reduced circuit was still sufficient to capture regulation of both glycolysis and OXPHOS, as well as the major differences in metabolism between normal and cancer cells. As the regulatory links in the core circuit are supported by multiple experiments (Supplementary Table S1), the circuit model is expected to be robust for studying the behavior of glycolysis and OXPHOS.
Normal cells typically use either OXPHOS or glycolysis at a fixed time depending on the availability of energy sources, because of tight regulation of the metabolic circuit and the competition between different metabolic modes. Cancer cells, however, kidnap the same gene regulatory circuit of metabolism for their own advantage. Because of either high energetic and oxidative stresses or the activation of specific oncogenes, cancer cells can be in a hybrid metabolic phenotype utilizing both glycolysis and OXPHOS. As cells in the hybrid “W/O” state have increased plasticity, they could have an advantage in survival over cells with the other phenotypes. The modeling results are supported by the TCGA patient data from multiple cancer types. The model also explains the phenomenon of oxygen shock (47). Cells that largely rely on glycolysis are exposed to large amount of oxygen when they reach the blood vessel. These cells could switch the metabolic phenotype to the “W/O” state, which further induces metastasis (10). In future work, the landscape approach (48, 49) could be utilized to quantify the transition processes among glycolysis, OXPHOS, and hybrid metabolism phenotypes. Note that, in addition to the use of glycolysis and OXPHOS, cancer cells can also utilize glutamine as a nitrogen source (5), while cancer-associated fibroblasts can perform glycolysis to fuel cancer cells by producing lactate, known as reverse Warburg effects (50). These aspects are out of the scope of this current study, but they are worth further investigation by using modeling approach.
To help test the basic motions of our core circuit model, we developed two metabolic signatures, one for glycolysis based on the expressions of HIF-1 downstream genes and the other for OXPHOS based on the expression of AMPK downstream genes. These metabolic signatures were applied to distinguish patient samples from multiple cancer types or single cells with different metabolic phenotypes, to identify the anticorrelations between the AMPK and HIF-1 activities, and to elucidate the association of oncogenic pathways with the metabolic states. In general, we expect these metabolic signatures to be a powerful tool to predict the metabolic phenotypes of cancer cells directly from gene expression data.
We proposed to design metabolic therapies by considering the hybrid metabolic phenotype (W/O). A putative strategy would be to drive cancer cells away from the hybrid states. Guided by this assumption, we evaluated the effectiveness of several therapeutic strategies by model simulation. We found that therapy using metformin is more likely to drive cells into the hybrid metabolic phenotype during the treatment time-course than the therapies using glycolysis blockers, such as 3BP and 2DG. Our model also provides an explanation of why some combined therapies are more effective (45, 46), because these therapies are efficient in driving cells out of the hybrid state.
Disclosure of Potential Conflicts of Interest
No potential conflicts of interest were disclosed.
Authors' Contributions
Conception and design: L. Yu, M. Lu, E. Ben-Jacob, B.A. Kaipparettu, J.N. Onuchic
Development of methodology: L. Yu, M. Lu, D. Jia, E. Ben-Jacob, H. Levine, B.A. Kaipparettu, J.N. Onuchic
Analysis and interpretation of data (e.g., statistical analysis, biostatistics, computational analysis): L. Yu, M. Lu, D. Jia, E. Ben-Jacob, H. Levine, B.A. Kaipparettu
Writing, review, and/or revision of the manuscript: L. Yu, M. Lu, D. Jia, H. Levine, B.A. Kaipparettu, J.N. Onuchic
Administrative, technical, or material support (i.e., reporting or organizing data, constructing databases): L. Yu, B.A. Kaipparettu
Study supervision: M. Lu, J. Ma, E. Ben-Jacob, B.A. Kaipparettu, J.N. Onuchic
Grant Support
J.N. Onuchic and H. Levine were awarded the Physics Frontiers Center NSF grant PHY-1427654. H. Levine was awarded the NSF grant DMS-1361411. J.N. Onuchic and H. Levine were awarded the Cancer Prevention and Research Institute of Texas (CPRIT) grants R1110 and R1111, respectively. M. Lu had a training fellowship from the Keck Center for Interdisciplinary Bioscience Training of the Gulf Coast Consortia (CPRIT Grant RP140113). J. Ma is thankful for support from the National Institutes of Health (R01-GM067801, R01-GM116280) and the Welch Foundation (Q-1512). E. Ben-Jacob was also supported by the Tauber Family Funds and the Maguy-Glass Chair in Physics of Complex Systems. B.A. Kaipparettu was awarded by the National Cancer Institute (NCI) grants R21CA173150 and R21CA179720.
The costs of publication of this article were defrayed in part by the payment of page charges. This article must therefore be hereby marked advertisement in accordance with 18 U.S.C. Section 1734 solely to indicate this fact.
Footnotes
Note: Supplementary data for this article are available at Cancer Research Online (http://cancerres.aacrjournals.org/).
This work is a joint collaboration of the two laboratories of B.A. Kaipparettu and J.N. Onuchic.
↵†Deceased
- Received August 5, 2016.
- Revision received December 2, 2016.
- Accepted December 19, 2016.
- ©2017 American Association for Cancer Research.
References
Volume 77, Issue 7
