Metabolic Modeling of the C3-CAM Continuum Revealed the Establishment of a Starch/Sugar-Malate Cycle in CAM Evolution

The evolution of Crassulacean acid metabolism (CAM) is thought to be along a C3-CAM continuum including multiple variations of CAM such as CAM cycling and CAM idling. Here, we applied large-scale constraint-based modeling to investigate the metabolism and energetics of plants operating in C3, CAM, CAM cycling, and CAM idling. Our modeling results suggested that CAM cycling and CAM idling could be potential evolutionary intermediates in CAM evolution by establishing a starch/sugar-malate cycle. Our model analysis showed that by varying CO2 exchange during the light period, as a proxy of stomatal conductance, there exists a C3-CAM continuum with gradual metabolic changes, supporting the notion that evolution of CAM from C3 could occur solely through incremental changes in metabolic fluxes. Along the C3-CAM continuum, our model predicted changes in metabolic fluxes not only through the starch/sugar-malate cycle that is involved in CAM photosynthetic CO2 fixation but also other metabolic processes including the mitochondrial electron transport chain and the tricarboxylate acid cycle at night. These predictions could guide engineering efforts in introducing CAM into C3 crops for improved water use efficiency.


INTRODUCTION
Crassulacean acid metabolism (CAM) photosynthetic CO 2 fixation is an evolutionary descendant of C 3 photosynthesis. CAM photosynthesis is known to have evolved independently multiple times in at least 35 plant families comprising about 6% of flowering plant species (Winter and Smith, 1996a;Silvera et al., 2010). CAM is an adaptation of photosynthetic CO 2 fixation typically associated to limited water availability (Cushman and Borland, 2002). By closing their stomata during the light period and fixing atmospheric and/or respiratory carbon dioxide (CO 2 ) exclusively in the dark period, CAM allows plants to use water more efficiently while fixing carbon for growth (Figure 1). The engineering of CAM into C 3 crops has been suggested as a possible strategy to meet the demands on agriculture for food, feed, fiber, and fuels, without exacerbating the pressures on arable land area due to climate change (Borland et al., 2014). However, as a carbon-concentrating mechanism, CAM is thought to be more metabolically expensive than C 3 FIGURE 1 | Schematic diagram of carbon flow in C 3 , CAM, CAM cycling, and CAM idling. Gray and yellow backgrounds represent the dark and light phases, respectively. The recycling of carbon skeleton from carbohydrate to organic acids in the dark phase was not illustrated for simplicity. (Winter and Smith, 1996b), which suggests that transferring a CAM pathway into C 3 crops would incur a crop yield penalty. To investigate the energetics of C 3 and CAM, large-scale metabolic models were applied which showed that engineering CAM into C 3 plants does not impose a significant energetic penalty given the reduction in photorespiration from the carbon-concentrating mechanism (Cheung et al., 2014;Shameer et al., 2018).
Besides the phylogenetic and ecological diversity of CAM plants, there is remarkable plasticity in its metabolism with multiple defined variations of CAM including, CAM cycling, and CAM idling (Ting, 1985;Cushman, 2001;Winter, 2019). CAM cycling primarily fixes CO 2 in the light period with refixing respiratory CO 2 behind closed stomata at night, leading to a small diel organic acid flux; CAM idling lacks diel gaseous exchange with closed stomata across the 24 h light/dark cycle and has a small continued diel fluctuation in organic acids level (Figure 1; Sipes and Ting, 1985;Cushman, 2001;Winter, 2019). Various variations of CAM, including CAM cycling, and CAM idling, were grouped into a "C 3 + CAM" phenotype, which is suggested to be an evolutionary accessible and stable precursor to CAM (Edwards, 2019). Silvera et al. (2010) generalized the idea of plasticity of CAM into a continuum of CAM levels, due to the differences in the degree of nocturnal and daytime net CO 2 uptake. Bräutigam et al. (2017) took the idea further to include C 3 as part of the C 3 -CAM continuum and suggested that the evolution of C 3 to CAM only required incremental increases in metabolic fluxes.
The changes in diel gaseous exchange and fluctuation in organic acids level along the C 3 -CAM continuum are often accompanied by a change in diel fluctuation in carbohydrate, mostly starch and/or sugars. These changes in turn have a knock-on effect on other parts of the metabolic system. An increasing popular approach for studying systems-level behavior of metabolism is constraint-based modeling, which has been successfully applied to study the metabolism and energetics of C 3 and CAM plants (Cheung et al., 2014;Shameer et al., 2018). Constraint-based modeling involves the application of biologically relevant constraints, such as the export of sucrose and amino-acids into the phloem for modeling leaf metabolism, to a stoichiometric model of metabolism. An objective function, typically related to energetic efficiency such as minimizing the amount of nutrients required or maximizing the production of biomolecules, is applied to obtain biologically relevant predictions of metabolic fluxes. One of the strengths of constraint-based modeling is that it is a predictive method, which is particular useful for investigating systems that are difficult to study experimentally such as the evolution of metabolic systems. In this study, large-scale constraint-based metabolic modeling was applied to investigate how CAM cycling and CAM idling fit into the continuum of CAM evolution and to identify the changes in metabolic fluxes along the C 3 -CAM continuum. The results from our modeling study provide novel insights into the energetics and metabolic alterations from C 3 to CAM, which could guide engineering efforts aimed at introducing CAM into C 3 plants.

Core Metabolic Model for Modeling C 3 , CAM, CAM Cycling, and CAM Idling
The mass-and charge-balance core plant metabolic model in Shameer et al. (2018) was used for modeling the metabolism of leaves operating in C 3 , CAM, CAM cycling, and CAM idling. A number of minor modifications were made to the Shameer et al. (2018) model to more accurately model metabolism of C 3 , CAM, CAM cycling, CAM idling, and the C 3 -CAM continuum. Firstly, a reaction ("OXYGEN_MOLECULE_e_dielTransfer") for the accumulation of oxygen, which is produced from water splitting in the photosynthetic light reactions, was added to the model such that we could run the simulations for CAM and CAM idling as oxygen exchange was constrained to zero during the day in these two scenarios. Another modification from Shameer et al. (2018) was how the acidification of the vacuole was modeled. Instead of directly setting different fixed pH for the vacuole in C 3 and CAM, protons were allowed to accumulate in the vacuole in our model ("PROTON_v_dielTransfer"). Reactions allowing protons to freely flow in and out of the cytosol ("unlProtHYPO_c1" and "unlProtHYPO_c2") were blocked such that pH homeostasis can be modeled through the accumulation of protons in the vacuole. For linking the cytosolic and mitochondrial proton pools, proton transport between the cytosol and the mitochondrial intermembrane space ("H_ic1" and "H_ic2") was set to be reversible. The modification of this constraint only led to a very minor change in the flux predictions (data not shown). Irreversible proton transporters were added from the vacuole to the cytosol ("H_vc1" and "H_vc2") and from extracellular to the cytosol ("H_ec1" and "H_ec2") to allow leakage of protons down the electrochemical gradients. Lastly, the compartmentation of metabolites in the reaction "HEXOKINASE_RXN_MANNOSE_c" was corrected to be in the cytosol. The modified core model can be found in SBML and Excel formats (Supplementary File 1). The scripts for model modification can be found in Supplementary File 2.

Model Simulations With Flux Balance Analysis
Based on the constraints and objective function stated in the "Results" section, parsimonious flux balance analysis (pFBA) was performed using scobra 1 , an extension of cobrapy (Ebrahim et al., 2013). The scripts for running the simulations in this study can be found in Supplementary File 2. In this study, we primarily reported the results from the pFBA simulations ( Supplementary Tables 1-3). The conclusions made based on the pFBA results for C 3 , CAM, CAM cycling, and CAM idling were confirmed using flux variability analysis (Mahadevan and Schilling, 2003) applied on the primary objective function, which identifies the flux range of each reaction in the optimal solution space (Supplementary Table 4).

Predicted Metabolic Fluxes of C 3 , CAM, CAM Cycling, and CAM Idling
In this study, we simulated the metabolism of leaves undergoing C 3 , CAM, CAM cycling, and CAM idling using a recently published core plant metabolic model which was used to model C 3 and CAM plants (Shameer et al., 2018). It is a constraintbased stoichiometric model capable of simultaneously simulating leaf metabolism over the diel cycle with a light phase and a dark phase. The two phases interact by allowing metabolites to accumulate in one phase that are utilized in the other phase. Minor modifications of the model were outlined in the "Materials and Methods" section. The constraints for simulating the core metabolic functions of mature leaves, namely export of 1 https://github.com/mauriceccy/scobra sucrose and amino acids into the phloem and ATP and NADPH costs for cellular maintenance, were set based on the values in Shameer et al. (2018). All simulations, except CAM idling, were constrained to have a phloem export rate of 0.259 µmol m −2 s −1 based on the value of C 3 plants in Shameer et al. (2018). The set of constraints for modeling the four different modes of photosynthesis are summarized in Table 1. The primary objective function of minimizing the total photon demand was used throughout this study, which allows us to study the metabolic efficiencies of the different modes of photosynthesis. From an evolutionary perspective, this primary objective function assumes maximum efficiency in energy conversation by the plants. pFBA, i.e., minimization of absolute sum of fluxes, was applied as a secondary objective to eliminate substrate cycles. Flux variability analysis performed on the primary objective (Mahadevan and Schilling, 2003) was used to confirm that the conclusions made based on the pFBA results were not affected by the presence of alternative optimal solutions (Supplementary Table 4).
The model predictions of C 3 and CAM were very similar to that in Shameer et al. (2018) given the similarities in the constraints used. Without any constraints on malate decarboxylation enzyme and carbohydrate storage, the model predicted net carbon fixation during the light period in the C 3 flux prediction, whereas in CAM carbon was fixed in the dark period with phosphoenolpyruvate carboxykinase (PEPCK) being the main predicted route for malate decarboxylation. Starch was predicted to be the main carbohydrate storage in both C 3 and CAM. These results are consistent with the findings in Shameer et al. (2018) where starch-storing PEPCK subtype were predicted to be the most energy efficient. The effect of the choice of decarboxylation enzymes (PEPCK vs malic enzyme) on the model predictions was explored by constraining other decarboxylating enzymes to carry zero flux. It was found that the choice of decarboxylation enzymes makes little qualitative difference with respect to the results presented (Supplementary Table 5). From here on, the results presented were model predictions with no constraints on the decarboxylation enzymes. As for carbohydrate storage, simulations were performed with starch, sucrose or fructan as the sole storage compound. Except for reactions involved in the synthesis, accumulation and degradation of carbohydrate storage, the predicted fluxes in central carbon metabolism were largely similar between the three carbohydrate storages tested (Supplementary Table 6). In this study, we mostly presented the results from simulations with starch as the carbohydrate storage. Similar conclusions can be made for using sugar as the carbohydrate storage. The core set of metabolic fluxes for C 3 , CAM, CAM cycling, and CAM idling with starch as the carbohydrate storage is depicted in Figure 2.

CAM Cycling
Similar to C 3 plants, CAM cycling fixes carbon in the light period. CAM cycling is characterized by closed stomata in the dark period with refixation of respiratory CO 2 and a small diel organic acid flux (Sipes and Ting, 1985;Cushman, 2001;Winter, 2019). To model CAM cycling, we applied the C 3 constraints with an additional constraint of setting CO 2 and O 2 exchange at night to zero to simulate the closure of the Phloem export rate was set based on the predicted value of C 3 plants in Shameer et al. (2018). RuBisCO carboxylase:oxygenase ratio was set to 3:1 when stomata is opened [a typical value for C 3 plants (Gutteridge and Pierce, 2006)], and 5.15:1 when stomata is closed [the median value from experimental estimates of CAM plants (Lüttge, 2010)] similar to Shameer et al. (2018). The photorespiratory pathway is shown in chloroplast for simplicity, which in reality spans multiple compartments. Flux from 3-phosphoglycerate to PEP was taken as the flux for glycolysis and gluconeogenesis. Flux for succinate dehydrogenase was taken as the TCA cycle flux. Note that only cyclic mode of the TCA cycle was shown. RuBisCO carboxylase flux was taken as the flux through the Calvin-Benson cycle.
stomata (Table 1). This resulted in a flux distribution that resembled a weak version of CAM, with nocturnal malate accumulation and increased light period starch accumulation ( Figure 2C). Phosphoenolpyruvate carboxylase (PEPC) was predicted to be active only at night in CAM cycling for CO 2 refixation, in contrast to C 3 where PEPC was only active during the light period (Supplementary Table 2). Another major difference between CAM cycling and C 3 is malate accumulation. While C 3 was predicted to have a very small amount of malate accumulation during the light period, CAM cycling was predicted to have substantial amount of nocturnal malate accumulation (∼20% of the amount of malate accumulation in CAM) (Figure 2; Supplementary Table 2), which is consistent with known behavior of CAM cycling (Ting, 1985;Cushman, 2001). The nocturnal malate accumulation and respiratory CO 2 refixation via PEPC under the CAM cycling scenario were accompanied by changes in fluxes in other parts of metabolism. Malate decarboxylation during the light period was predicted to be active in CAM cycling but not in C 3 (Figure 2). In CAM cycling, flux through gluconeogenesis was predicted to convert malate into starch in the light period, which led to 2.68-fold more starch accumulation during the light period in CAM cycling compared to C 3 (Figure 2; Supplementary Table 2). Given that CAM cycling has a higher starch accumulation in the light period, it was predicted to have a 6.53-fold larger glycolytic flux in the dark to convert starch into phosphoenolpyruvate (PEP) for CO 2 refixation, compared to C 3 (Figure 2; Supplementary Table 2). The activities of most of the other reactions at night were similar in CAM cycling and in C 3 , with CAM cycling having a slightly higher flux (15-17%) through the tricarboxylic acid (TCA) cycle and the mitochondrial electron transport chain (ETC) than C 3 , presumably to produce extra ATP for transporting malate into the vacuole for storage at night.

CAM Idling
Crassulacean acid metabolism idling is characterized by the lack of diel gaseous exchange and a small continued diel fluctuation in the organic acids level because of internally recycled CO 2 (Sipes and Ting, 1985;Winter, 2019). It is usually an adaptation in water-stressed plants, which results in the closure of stomata for the whole 24-h cycle. To model this, the CO 2 and O 2 exchange during the light and the dark periods were constrained to carry zero flux ( Table 1). Given that there is no CO 2 exchange, we assumed that there is no net carbon fixation, hence phloem export was constrained to zero for CAM idling. The primary metabolic demand for plants in CAM idling is cellular maintenance. The model predicted a starch-malate cycle where starch accumulated in the light period is metabolized in the dark period mainly through glycolysis and the oxidative pentose phosphate pathway (OPPP) to produce ATP and NADPH for maintenance processes (Figure 2D). While the majority of PEP was used as precursor for carbon refixation by PEPC, a significant proportion (17%) of PEP was predicted to be metabolized further through the TCA cycle to feed the mitochondrial ETC for ATP synthesis (Figure 2D). Given that it is a closed system with respect to carbon, CO 2 produced in the OPPP and the TCA cycle is refixed by PEPC, which ultimately leads to the accumulation of malate in the dark. In the light period, PEP from malate decarboxylation was recycled to produce starch via gluconeogenesis, while the CO 2 produced from malate decarboxylation was refixed via the Calvin-Benson cycle similar to the scenario for CAM (Figure 2). With no net carbon import or export, the amount of carbon stored in starch in the light period was predicted to be equaled to the amount of carbon storage in malate at night. The starch-malate cycle was primarily driven by the energy from the light reactions of photosynthesis, and it acted as a carbon neutral way of storing and transferring energy from the light period to the dark period. Similar results were obtained when sucrose or fructan was used as the sole carbohydrate storage instead of starch (Supplementary Table 6), meaning that a sugar-malate cycle can serve the same function as the starch-malate cycle in sugarstoring plants.

Energetics and Metabolite Accumulation in C 3 , CAM, CAM Cycling, and CAM Idling
The metabolic flux predictions of C 3 , CAM, CAM cycling, and CAM idling were compared to see how CAM cycling and CAM idling fit into the evolution of CAM from C 3 . Table 2 summarizes the predicted fluxes related to energetics and metabolic accumulation in the four simulations. CAM idling was predicted to use the fewest photons, which was expected given that it does not have the metabolic demand for exporting sucrose and amino acids into the phloem. For the same metabolic demand, CAM requires more photons than C 3 , as expected. It is interesting to see that the photon demand for CAM cycling  Frontiers in Plant Science | www.frontiersin.org falls between C 3 and CAM. A similar trend was observed for other fluxes related to energy metabolism including the ATP and NADPH production by the photosynthetic light reactions and the ATP production by the mitochondrial ATP synthase ( Table 2). The same trend was also reflected in the energetic demands of the Calvin-Benson cycle in terms of ATP and NADPH consumption (Supplementary Table 2). Metabolite accumulation showed a different pattern compared to the energetics ( Table 2). C 3 had the lowest daytime starch accumulation, followed by CAM cycling and CAM idling which had about 2-3 times more starch accumulation than C 3 . CAM had the highest light period starch accumulation with more than nine times the amount associated with C 3 . This suggested that CAM cycling and CAM idling could potentially be intermediate steps in CAM evolution with respect to the regulation of starch accumulation. A similar pattern can be observed for malate accumulation. A very small amount of malate was predicted to accumulate during the day for C 3 plants, whereas a large nocturnal malate accumulation was predicted for CAM as part of CAM photosynthesis. CAM cycling and CAM idling had intermediate level of nocturnal malate accumulation (∼20% of that in CAM), which was related to the refixation of nocturnal CO 2 by PEPC. Reactions related to the starch/sugar-malate cycle, including glycolysis and PEPC flux in the dark period, and gluconeogenesis and malate decarboxylation during the light period, showed a similar trend (Supplementary

Predicting the Metabolic Transitions During C 3 -CAM Evolution
The behavior of diel CO 2 exchange is the main diagnostic indicator between C 3 and CAM (Silvera et al., 2010). To model the potential metabolic transitions that could happen during the evolution of CAM from C 3 , we varied the CO 2 uptake rate during the light period from 13.12 µmol m −2 s −1 (the predicted value for C 3 ) to 0 µmol m −2 s −1 (which had the same effect as gradually increasing nocturnal CO 2 uptake given the assumption of a fixed carbon export to the phloem). This simulates the decrease in gaseous exchange during the light period by stomatal closure, hence a similar constraint was set for light period oxygen exchange. As the stomata closes in the light period, i.e., light period CO 2 uptake decreases, the proportion of ribulose-1,5-bisphosphate carboxylase/oxygenase (RuBisCO) flux going through the carboxylase reaction was assumed to increase linearly from 75% [carboxylase to oxygenase ratio of 3:1; a typical value for C 3 plants (Gutteridge and Pierce, 2006)] to 83.74% [carboxylase to oxygenase ratio of 5.15:1; the median value from experimental estimates of CAM plants (Lüttge, 2010)] to account for the reduction of photorespiration from the increase in CAMlike behavior with increasing internal CO 2 concentration from malate decarboxylation. All other constraints remained the same as the C 3 and CAM simulations. This analysis aimed to explore the metabolic transitions during the evolution of CAM from C 3 assuming the existence of an increasing evolutionary pressure to reduce diurnal CO 2 uptake, e.g., water conservation. The full results from this simulation can be found in Supplementary Table 3.
Given that the metabolic demands remained constant throughout the analysis, a decrease in CO 2 uptake in the light period led to a shift from C 3 to CAM photosynthesis with an increase in flux through the starch-malate cycle including starch degradation, glycolysis, PEPC, and malate accumulation at night, and malate decarboxylation and starch accumulation during the light period (Figures 3A,B and Supplementary Table 3). Note that dark period CO 2 uptake increased as light period CO 2 uptake decreased due to the carbon balance of the model in exporting a fixed amount of sucrose and amino acids into the phloem. CAM cycling occurs at the point when dark period CO 2 uptake is zero.
A similar analysis on a rice model with varying RuBisCO carboxylase to oxygenase ratio but without a change in light period CO 2 uptake showed that decreasing RuBisCO oxygenase contribution led to a decrease in photon influx (Chatterjee et al., 2017). In our analysis, despite the constrained decrease in RuBisCO oxygenase contribution as light period CO 2 uptake decreased, the amount of energy (in terms of photons) required to sustain the same metabolic demand increased by about 7% from C 3 to CAM (Figure 3C) as extra energy is needed to run the starch-malate cycle. This is correlated with the increase in flux through the photosynthetic light reactions. Besides plastidial ATP synthesis, there was also an increase in ATP synthesis by the mitochondrial ETC in the light period as the simulation shifted from C 3 to CAM (Figure 3D). The contribution of mitochondrial ATP synthesis increased from 18.2% in C 3 to 35.6% in CAM (Figure 3E), which is likely to be related to the increase in NADH produced during malate decarboxylation. In our simulations, the RuBisCO carboxylase flux was predicted to be remain relatively constant while the total RuBisCO flux (carboxylase + oxygenase) decreased from C 3 to CAM due to the decrease in RuBisCO oxygenase activity ( Figure 3F). There were two major factors affecting RuBisCO carboxylase flux, (i) refixation of photorespiratory CO 2 , and (ii) starch accumulation to support energy demand in the dark period. In this case, the two factors counteract each other throughout the simulation where photorespiration decreases and the energy demand for running the starch-malate cycle (mostly for pumping malate into the vacuole) increases from C 3 to CAM. For the simulations with sucrose or fructan as the sole carbohydrate storage, the model predicted an increase in RuBisCO carboxylase flux from C 3 to CAM as the energy required for running the sugar-malate cycle is higher than the starch-malate cycle (due to the cost of pumping sugars into the vacuole for storage).
During the night, other than the increase in glycolytic flux as part of the starch-malate cycle from C 3 to CAM, the model predicted an 87% increase in flux through the TCA cycle and an 83% increase in flux through the mitochondrial ETC ( Figure 3G). This increase in mitochondrial ATP synthesis was mostly used to support the ATP-dependent tonoplast proton pump for the increasing nocturnal vacuolar malate accumulation. The cytosolic OPPP flux was predicted to decrease by 30% in the night from C 3 to CAM (Figure 3H). This could be explained by the increase in the TCA cycle flux which contributed to the production of NADPH in the mitochondrion by the NADP-isocitrate dehydrogenase. This lessened the demand for the production of cytosolic NADPH required to be shuttled into the mitochondrion for maintenance processes.

CAM Cycling and CAM Idling as Viable Evolutionary Steps for Establishing the Starch-Malate Cycle
CAM cycling is considered as a weak form of CAM in which stomata are open during the day and are closed at night (Lüttge, 2004;Silvera et al., 2010;Winter, 2019). With these constraints, our model predicted the known features of CAM cycling including the refixation of respiratory CO 2 in the dark period, and a small amount of nocturnal malate accumulation (Cushman, 2001;Winter, 2019). To support these metabolic behaviors, our model predicted the establishment of a starchmalate cycle in CAM cycling, which included increased flux through malate decarboxylation, gluconeogenesis and starch synthesis and accumulation during the light period, and starch degradation and glycolysis during the dark period, when compared to C 3 plants. Based on our modeling results, the main metabolic advantage of CAM cycling over C 3 is its higher carbon conversion efficiency when photosynthesis is limited by stomatal conductance in the light period, i.e., carbon limited. Given the same metabolic outputs, CAM cycling was predicted to require 20% less external CO 2 compared to C 3 due to the refixation of nocturnal respiratory CO 2 . This comes with a minor cost of 4.8% more photons and 1.6% more RuBisCO activity required, assuming that there is no reduction in photorespiration, which could be affected by limiting stomatal conductance and internal CO 2 generation from malate decarboxylation. One of the proposed hypotheses for the adaptive benefits of CAM cycling was the conservation of respired CO 2 at night that would otherwise be lost (Martin and Zee, 1983;Martin and Jackson, 1986). While the hypothesis was thought to be untenable based on energetic considerations (Martin et al., 1988), our modeling results suggest that conservation of respired CO 2 could be a potential adaptive benefit for CAM cycling under carbon-limited conditions. Given an environment that limits stomatal conductance in the light period, e.g., high temperature and drought, the evolution of CAM cycling, together with the establishment of the starch/sugar-malate cycle, was predicted to be advantageous in maximizing carbon conversion efficiency. The metabolic activities of all reactions in the starch-malate cycle in CAM cycling were predicted to be at an intermediate level between C 3 and CAM. The same applies to other supporting reactions such as the TCA cycle in the dark and the mitochondrial ETC during the light and dark periods. Our findings agree with studies in the Crassulaceae (Monson, 1989) and the Portulacaceae (Guralnick and Jackson, 2001), suggesting that CAM cycling is likely to be a possible evolutionary step along the path to the evolution of CAM.
Crassulacean acid metabolism idling is considered to be a survival mechanism of severely droughted CAM plants (Winter, 2019). In CAM idling, stomata remain closed throughout the day and night with small, sustained diel fluctuations in organic acids (Cushman, 2001;Silvera et al., 2010;Winter, 2019). By constraining our model with closed stomata in both the light and dark periods, the model predicted the operation of the starch/sugar-malate cycle as the most energy efficient way to sustain cellular activities. This is consistent with the experimental measurements of diel fluctuation in titratable acidity in plants operating in CAM-idling (Rayder and Ting, 1983;Sipes and Ting, 1985). With the lack of experimental studies on the metabolism on CAM-idling, our modeling results provide a valuable resource in understanding the metabolism of plants operation in CAMidling. CAM-idling has been suggested to be an important mechanism from survival (Hanscom and Ting, 1978). From an evolutionary perspective, if a plant often experiences conditions that require the closure of stomata throughout day and night, such as long periods of severe drought, the evolution of CAM idling would be advantageous for the plant to stay alive. While the evolution of CAM through CAM cycling seems more likely given its similarities to C 3 , it is not impossible that some lineages could establish the starch/sugar-malate cycle through CAM idling.

Stomatal Conductance as a Determinant Along the C 3 -CAM Continuum
It has been proposed that CAM evolution occurs along a continuum from C 3 to CAM (Silvera et al., 2010;Bräutigam et al., 2017). Our model analysis supports this proposal. By varying the CO 2 exchange in the light period, as a proxy for stomatal conductance, our model predicted a C 3 -CAM continuum with gradual metabolic changes along the continuum (Figure 3). A recent modeling study demonstrated that CAM could be a result of a trade-off between water saving and leaf productivity (Töpfer et al., 2020), which is in line with our results that CAM could emerge under an evolutionary pressure for reduced light period stomatal conductance such as water saving. The key metabolic changes along the C 3 -CAM continuum included the processes in the starch/sugar-malate cycle, the TCA cycle at night, and the chloroplastic and mitochondrial ETCs, which are consistent with what we know about metabolic fluxes in CAM (Silvera et al., 2010). The fact that a gradual continuum was predicted to be the most energetically favorable way to adapt to a change in stomatal conductance suggests that the fitness landscape between C 3 and CAM with regards to metabolic changes is a smooth one. Given our results, it is not surprising to see many facultative CAM plants which can switch between C 3 and CAM triggered by environmental changes (Winter, 2019). Based on our model predictions, it is hypothesized that we could find plants operating anywhere on the C 3 -CAM continuum. For example, the early development of a Kalanchoe pinnata leaf showed an ontogenetic progression from C 3 to CAM, which went through a stage like CAM cycling (Winter, 2019). Given the flexibility shown in facultative CAM plants and our results on the C 3 -CAM continuum, it could be possible to find existing plants or engineer new plants that can switch not only between C 3 and CAM but also at different points on the continuum depending on the environmental conditions.

CONCLUSION
Using a core plant metabolic model, we were able to model the metabolic behaviors of CAM, CAM cycling, and CAM idling by changing a few simple constraints on gaseous exchange and phloem export. Our modeling results suggest that CAM cycling and CAM idling could potentially be evolutionary intermediates on the path to CAM evolution by establishing an intermediate flux through the starch/sugar-malate cycle. By varying the light period CO 2 exchange as a proxy for stomatal conductance, the model predicted a continuum from C 3 to CAM with gradual metabolic changes, suggesting a smooth fitness landscape between C 3 and CAM in terms of changes in metabolism. Besides the insights gained in CAM evolution, the results from this study are informative to guide engineering efforts aiming to introduce CAM into C 3 crops by identifying the metabolic changes required to convert C 3 to CAM. In additional to the starch/sugar-malate cycle involved in CAM photosynthesis, our model showed that the fluxes of other metabolic processes, including the TCA cycle and the mitochondrial ETC, need to be altered from C 3 to optimize the operation of CAM.

DATA AVAILABILITY STATEMENT
The original contributions presented in the study are included in the article/Supplementary Material, further inquiries can be directed to the corresponding author/s.