C4 photosynthesis in C3 rice: a theoretical analysis of biochemical and anatomical factors.

Engineering C4 photosynthesis into rice has been considered a promising strategy to increase photosynthesis and yield. A question that remains to be answered is whether expressing a C4 metabolic cycle into a C3 leaf structure and without removing the C3 background metabolism improves photosynthetic efficiency. To explore this question, we developed a 3D reaction diffusion model of bundle-sheath and connected mesophyll cells in a C3 rice leaf. Our results show that integrating a C4 metabolic pathway into rice leaves with a C3 metabolism and mesophyll structure may lead to an improved photosynthesis under current ambient CO2 concentration. We analysed a number of physiological factors that influence the CO2 uptake rate, which include the chloroplast surface area exposed to intercellular air space, bundle-sheath cell wall thickness, bundle-sheath chloroplast envelope permeability, Rubisco concentration and the energy partitioning between C3 and C4 cycles. Among these, partitioning of energy between C3 and C4 photosynthesis and the partitioning of Rubisco between mesophyll and bundle-sheath cells are decisive factors controlling photosynthetic efficiency in an engineered C3 -C4 leaf. The implications of the results for the sequence of C4 evolution are also discussed.

Introduction

In C3 species, photosynthesis mainly occurs in mesophyll cells of leaves. Carbon dioxide (CO2) diffuses through stomata and intercellular spaces before it enters these cells and reacts with Ribulose bisphosphate (RuBP) in the chloroplasts to form the 3‐carbon compound 3‐phosphoglycerate. This reaction is catalysed by RuBP carboxylase–oxygenase (Rubisco) (Lorimer 1981; Hall & Rao 1999). Oxygen can also be fixed by Rubisco and leads to the production of 2‐phosphoglycolate. This compound is subsequently recycled to RuBP by the photosynthetic carbon oxidation (PCO) cycle. The latter process is associated with additional energy cost and results in the release of some CO2, decreasing photosynthetic efficiency by about one third (Ehleringer & Monson 1993).

In C4 species, the initial fixation reaction is not catalysed by Rubisco, but instead, CO2 is converted to bicarbonate (HCO3 −) in the cytosol of mesophyll cells and subsequently fixed into C4 acids by phosphoenolpyruvate carboxylase (PEPC). The C4 acids are transported from mesophyll cell into bundle‐sheath cells through plasmodesmata and decarboxylated back to CO2. This CO2 is subsequently refixed by Rubisco located in bundle‐sheath chloroplasts. The high affinity of PEPC to bicarbonate and the low permeability of the bundle‐sheath‐to‐mesophyll cell interface leads to a high CO2 concentration in bundle‐sheath cells. This high CO2 concentration allows for much lower Rubisco oxygenation and photorespiration rates (Furbank et al. 2004). Three types of C4‐photosynthesis have been historically distinguished based on the enzymes catalysing the C4‐acid decarboxylation: an NADP‐ME type, an NAD‐ME type and a PCK type (von Caemmerer & Furbank 2003). In many plants, these different decarboxylation pathways occur simultaneously to varying degrees (Lorimer 1981; Sommer et al. 2012; Muhaidat & McKown 2013; Wang et al. 2014a). In plants relying predominantly on the NADP‐ME pathway, such as Sorghum biocolor, Miscanthus and Zea mays (maize), decarboxylation of C4 acids in bundle‐sheath chloroplasts is mainly catalysed by NADP‐malic enzyme (NADP‐ME). C4 plants usually have an inherent higher photosynthetic CO2 uptake rate and higher conversion efficiency of solar energy compared to C3 plants (Zhu et al. 2008, 2010). Besides higher light use efficiency, C4 plants have higher use efficiency of water and nitrogen and have higher yields under warmer temperatures (Hibberd et al. 2008). These features suggest that expressing a C4 pathway in C3 crop species may be a useful strategy to improve crop yields (Leegood 2002).

In addition to differences in metabolism, C3 and C4 plants have different anatomical features. In rice, irregularly arranged, heavily lobed mesophyll cells are located between the bundle‐sheath cells. Chloroplasts in these cells occupy about 66% of the protoplast volume and cover about 97% of cell periphery, which is thought to maximize the diffusive conductance of CO2 into the stroma (Sage & Sage 2009). Rice bundle‐sheath cells contain fewer chloroplasts than mesophyll cells, and the chloroplasts do not form a continuous boundary around the periphery of the cell, with only 21% to 52% of the cell surface covered by chloroplasts (Sheehy et al. 2008). The leaves of many C4 plants are characterized by a so‐called Kranz anatomy, that is, each vascular bundle is surrounded by an inner ring of large bundle‐sheath cells and an outer ring of mesophyll cells (Furbank et al. 2004). C4 plants have fewer mesophyll cells between neighbouring bundle‐sheath cells and the interval distance between neighbouring bundle‐sheath cells is shorter (Dengler et al. 1994). The mesophyll cells in, for example, maize are lobed, but not so extensively as in rice (Giannoutsou et al. 2013; Warner et al. 2014). There are fewer chloroplasts in C4 mesophyll cells compared to those of related C3 species, and they do not cover the complete cell periphery (Stata et al. 2014). In addition, maize bundle‐sheath cells are larger than in rice and contain large chloroplasts, which are centrifugally (towards the mesophyll cells) arranged (Maai et al. 2011). Generally, low inter‐veinal distances and a high ratio between bundle‐sheath and mesophyll cell volume are characteristic for C4 lineages (Griffiths et al. 2013), but it remains unclear whether the rice leaf anatomy precludes an efficient C4 photosynthetic cycle.

Most of current work related to C4 engineering focuses on expressing a complete NADP‐ME type C4 metabolism in a C3 crop (Kajala et al. 2011; Miyao et al. 2011; Sage & Zhu 2011). Biochemical factors that affect this efficiency have been systematically evaluated earlier (Laisk & Edwards 2000; Wang et al. 2014b). However, it is not clear whether it is advantageous to engineer such a C4 metabolic cycle into a typical C3 leaf without removing the original C3 metabolic processes and by utilizing a C3 leaf anatomy. This is an important question to answer because it will determine whether the current C4 engineering work needs to remove the existing C3 metabolic cycle from C3 leaves and whether anatomical changes to the leaf are necessary. von Caemmerer (2003) analysed the efficiency of an engineered single‐cell C4‐type concentrating mechanism in rice. They found that a single‐cell approach limits the energy efficiency of C4 photosynthesis, and suggested that compartmentation of CO2 decarboxylation in the bundle‐sheath may be a far more successful strategy. Therefore, in this work we developed a reaction diffusion model that accounts for a two‐compartment C4 metabolism being expressed in a C3 background, while also accounting for the typical leaf anatomy observed in rice.

Materials and Methods

A number of different approaches were used in this study to describe the rate of photosynthesis in C3 leaves and in engineered leaves which contain elements of both C3 and C4 photosynthesis. Details of the reaction diffusion models, that is, the 3D model structure, the mass balance equations describing the rates of concentration changes of diffusible substrates and also the rate equations, are shown below and in Supplemental Files S1–S3. We also summarized the differences among these models in Table 1.

Table 1

Different model approaches used in this study

Name	Diffusion limitations	Enzyme limited reaction metabolites	Location
C3 reaction diffusion model	3D reaction diffusion modela	CO2, HCO3 −	Materials and Methods
C3–C4 reaction diffusion model	3D reaction diffusion modela	CO2, HCO3 −	Supplemental File S1
Extended C3–C4 reaction diffusion model	3D reaction diffusion modela	CO2, HCO3 −, OAA, malate, PEP, pyruvate	Materials and Methods, Supplemental File S2, S3
C3 biochemical model	Resistance modelb	CO2	Supplemental File S4
C3–C4 biochemical model	Resistance modelb	CO2	Supplemental File S4

The facilitating effect of CA was considered by accounting for the rates of CO2 hydration and HCO3 − dehydration (Tholen & Zhu 2011).

The facilitating effect of CA was considered by assuming full equilibrium between CO2 and HCO3 − (Evans et al. 2009).

3D structure

We constructed a two‐cell reaction diffusion model of C3 rice photosynthesis, which consists of a mesophyll cell connected to a bundle‐sheath cell. The CO2 concentration at the mesophyll cell boundaries was assumed to be in equilibrium with the intercellular air space of a typical leaf. Under this assumption, we modelled a typical rice mesophyll cell (Sage & Sage 2009) containing six lobes, and in each lobe we assumed there is a cluster of mitochondria and a layer of chloroplasts that nearly completely covers the cell wall adjacent to intercellular spaces. The lobed mesophyll cell was modelled as a combination of spheres: a central sphere was fused with six peripheral spheres of the same dimensions and at regular distances from each other, representing the cell lobes. The vacuole of the mesophyll cell was located in the middle of the central sphere. Each bundle‐sheath cell contained a layer of chloroplasts proximal to the mesophyll cell, followed by two clusters of mitochondria. The vacuole in the bundle‐sheath cell was located distal to the mesophyll cell. The 3D geometry of the mesophyll and bundle‐sheath cell, including all organelles is shown in Fig. 1.

Figure 1

The 3D model geometry representing a rice mesophyll (upper‐left) cell connected to a bundle‐sheath cell (lower‐right). Layers of chloroplasts are indicated in green, clusters of mitochondria in red, vacuoles are blue and the cytosol is gray.

The radius of each lobe in these cells was 4 μm. The distance from the centre of the central sphere to the centre of each lobe was 5.77 μm, and there was considerable overlap resulting in a total cell volume of 1550.6 μm3. The mesophyll chloroplast surface exposed to intercellular spaces was about 91.7%. The bundle‐sheath cell was built as an ellipsoid with length, width and depth being 10 μm, 8 μm and 8 μm, respectively. The surface area of bundle‐sheath chloroplast facing the mesophyll cell wall was about 22%. Additional parameters related to organelles in the model are listed in Table S1.

Mass balance equations

C3 reaction diffusion model

Biochemical reactions in the C3 rice model are based on Tholen & Zhu (2011). CO2 and HCO3 − and related enzymes (carbonic anhydrase (CA) and Rubisco) were included in this the model. We assumed that all enzymes are uniformly distributed throughout their respective organelles. Metabolites were allowed to freely diffuse within an organelle, and permeabilities were defined for transport between different compartments.

Extending the models presented in Cowan (1986) and Tholen & Zhu (2011), the equation of CO2 during the process of diffusion and reactions was described by: where ∇2[CO 2] = ∂2[CO 2]/∂x2 + ∂2[CO 2]/∂y2 + ∂2[CO 2]/∂z2, D (m2 s−1) is the liquid phase diffusion coefficient for CO2, η is a dimensionless factor which represents the relative viscosity of the compartment in which the diffusion takes place, [CO] (mol m−3) is the CO2 concentration, v (mol m−3 s−1) and v (mol m−3 s−1) are the rates of CO2 fixation by Rubisco in mesophyll and bundle‐sheath chloroplasts, respectively, h (mol m−3 s−1) is the net CO2 hydration rate, r (mol m−3 s−1) is the respiration rate, and v (mol m−3 s−1) and v (mol m−3 s−1) are the photorespiration rates in mesophyll and bundle‐sheath mitochondria, respectively. For Eqn 1 in the C3 reaction diffusion model, v ≠ 0 in bundle‐sheath chloroplast, v ≠ 0 in mesophyll chloroplast, r ≠ 0 in mesophyll and bundle‐sheath mitochondria, and v ≠ 0 and v ≠ 0 in bundle‐sheath mitochondria and mesophyll mitochondria, respectively.

Similarly, the equation for HCO3 − diffusion and reactions was: where D (m2 s−1) is the liquid‐phase diffusion coefficient for bicarbonate, and [HCO] (mol m−3) is the bicarbonate concentration.

C3–C4 reaction diffusion model

To model C4 photosynthesis in a C3 background, the model described above was extended with a rate equation describing the fixation of bicarbonate by PEPC. Strictly speaking, bicarbonate is converted into OAA by PEPC in the mesophyll cytosol. For better comparison with current biochemical models, we followed the assumption that the conversion into HCO3 − is not limiting the C4 cycle (von Caemmerer 2000). This means that the rate of carboxylation by PEPC (v in mol m−3 s−1) is accounted for directly in Eqn 1. A more realistic approach will be taken below under ‘Extended C’. Thus, the equation for CO2 becomes: where v (mol m−3 s−1) is the rate of carboxylation by PEPC, and v (mol m−3 s−1) is the decarboxylation rate in the bundle‐sheath chloroplast. The total rate of decarboxylation in the bundle‐sheath chloroplast was constrained to be equal to the rate of carboxylation by PEPC in the mesophyll cytosol. To incorporate the facilitating effect of CA on diffusion in the model (Cowan 1986; Evans et al. 2009; Tholen & Zhu 2011), we described HCO3 − diffusion using Eqn 2. Detailed reaction rate equations are given in Supplemental File S1.

Extended C3–C4 reaction diffusion model

Current models of C4 photosynthesis, including the reaction diffusion model described in the previous paragraph, assume that the enzyme‐limited rate of C4 photosynthesis is limited by Rubisco or by PEP carboxylation. We developed an extended C3–C4 reaction diffusion model by including the hydration reaction for CO2 in the C4 cycle and also adding C4 related metabolites: oxaloacetic acid (OAA), malate, pyruvate, phosphoenolpyruvate (PEP) and the C4 related enzymes: PEPC, malate dehydrogenase (NADP‐MDH), NADP‐ME, pyruvate and phosphate dikinase (PPDK). Biochemical reactions and CO2 flows are shown in (Fig 2). The equation for CO2 reactions and diffusion through the liquid phase is then: where v (mol m−3 s−1) is the decarboxylation rate of C4 acids in bundle‐sheath chloroplast. The equation for HCO3 − was: where v (mol m−3 s−1) is the rate of the reaction from bicarbonate to OAA catalysed by PEPC. Additional mass balance equations of C4 metabolites are given in Supplementary File S2.

Figure 2

Schematic overview of the biochemical reactions in a rice plant expressing a C4 metabolism. OAA: oxaloacetic acid; PEP: phosphoenolpyruvate; PEPC: phosphoenolpyruvate carboxylase; NADP‐MDH: malate dehydrogenase; NADP‐ME: NADP‐malic enzyme; PPDK: pyruvate and phosphate dikinase. Metabolites, reactions and enzymes are indicated in black. All metabolites shown in the diagram can diffuse between different compartments in mesophyll and bundle‐sheath cells. Blue arrows: CO2 flux.

Rate equations

Rate equations for the light reactions

Electron transport rate (J mol m−2 s−1) was calculated following Ögren & Evans (1993) and von Caemmerer (2000): where I (mol m−2 s−1) is the light absorbed by photosystem II (PS II) in the chloroplasts, J (mol m−2 s−1) is the maximum capacity of electron transport chain and θ is an empirical curvature factor assumed to be around 0.7 (Evans 1989; von Caemmerer 2000). I was calculated as: where I (mol m−2 s−1) is the incident irradiance, α is the absorptance of leaves (assumed to be 0.85) and f is the fraction of absorbed photons that do not drive electron generation and was set at 0.15 (Evans 1987). The factor 1/2 indicates that 50% of the energy is assumed to be absorbed by PS II.

Rate equations for metabolic processes

C3 reaction diffusion model

It has been suggested (von Caemmerer 2000; Leegood 2008) that the rice bundle‐sheath contributes to photosynthesis. Moreover, in barley, the concentration of Rubisco in bundle‐sheath cells is similar to that in mesophyll cells (Koroleva et al, 2000). We therefore assumed that the concentrations and kinetic properties of Rubisco between these two cell types in rice are the same. Rubisco carboxylation rate in mesophyll and bundle‐sheath cells were distributed on the basis of their relative chloroplast volume. The maximum carboxylation capacity per unit leaf area by Rubisco was assumed to be 80 μmol m−2 s−1. We assumed a linear electron transport chain operates in both mesophyll and bundle‐sheath chloroplast and photosynthesis is limited by the amount of ATP available. Assuming a Q cycle operates in the photosystem (Sacksteder et al, 2000; Kramer & Evans 2011), the ratio between proton transport cross thylakoid membrane and electron flow (H+/e−) is 3. The H+/ATP ratio varies (Kramer & Evans 2011), and here we assumed a ratio of 4 (Sheehy et al. 2000). Thus, the e−/ATP ratio for linear electron transport flow is 3/4. We further assumed that the ratio of electron transport rate between bundle‐sheath and mesophyll cells equals the ratio between the total volume of mesophyll and bundle‐sheath chloroplasts. The CO2 fixation rates in mesophyll (v mol m−3 s−1) and bundle‐sheath (v mol m−3 s−1) chloroplasts can therefore be described as: where K  = K (1 + [O 2]/K ), [CO] (mol m−3) and [CO] (mol m−3) are the CO2 concentrations in mesophyll and bundle‐sheath chloroplasts, V (μmol m−2 s−1) is the maximum carboxylation rate of Rubisco per unit leaf area, f and f are the fraction of bundle‐sheath chloroplasts and mesophyll chloroplasts volume relative to the total chloroplast volume respectively, K (mol m−3) is the effective Michaelis–Menten constant for Rubisco in mesophyll cell, Γ* (mol m−3) is the CO2 compensation point in the absence of mitochondrial respiration, K (mol m−3) is the Michaelis–Menten constant for Rubisco carboxylase, K (mol m−3) is the Michaelis–Menten constant for Rubisco oxygenase, [O] (mol m−3) is the O2 concentration, f and f are the fraction of energy partitioned for photosynthetic carbon reduction (PCR) cycle and PCO cycle in mesophyll and bundle‐sheath cells respectively, S (m2 m−2) is the mesophyll surface exposed to intercellular spaces area per unit leaf area, and V (m3 m−2) and V (m3 m−2) are mesophyll and bundle‐sheath chloroplast volumes per unit mesophyll surface exposed to intercellular spaces area.

The volumetric respiration rate, r (mol m−3 s−1), was assumed to be equal in mesophyll and bundle‐sheath mitochondria, and was described as: where R (mol m−2 s−1) is the respiration rate per unit leaf area, and V (m3 m−2) and V (m3 m−2) are mesophyll and bundle‐sheath mitochondria volumes per unit mesophyll surface exposed to intercellular spaces area.

In the C3 biochemical model, the rate of photorespiration is calculated as Rubisco carboxylation rate multiplied by the CO2 compensation point in the absence of respiration, and divided by the chloroplastic CO2 concentration (von Caemmerer 2000). Assuming that no photorespiratory intermediates are transported between mesophyll and bundle‐sheath cells, we represented the local volumetric photorespiration rate (v (mol m−3 s−1) in mesophyll mitochondria and v (mol m−3 s−1) in bundle‐sheath mitochondria) as the integral of Rubisco carboxylation rate over the chloroplasts volume multiplied by Γ* and dividing by the volume of mitochondria (see also Tholen & Zhu (2011)):

The net CO2 hydration rate h (mol m−3 s−1) was calculated as (modified from Spalding & Portis (1985)): where [HCO] (mol m−3) is the HCO3 − concentration, k (s−1) is the catalytic rate of CA, X (mol m−3) is the CA active site concentration per unit leaf area, K (mol m−3) is the equilibrium constant, and K and K (mol m−3) are Michaelis–Menten constants for hydration and dehydration, respectively, pH is the pH value. The values of X and pH differ among different organelles (Table 2).

Table 2

Default anatomical and biochemical parameters and constants used in the C3 and C3–C4 reaction diffusion model (at 25 °C)

Name	Symbol	Default value	Units	Notes and references
Oxygen concentration	[O2]	0.21	bar	Assuming 21% oxygen concentration
CO2 concentration in intercellular air space	Ci	9.18 × 10−3	mol m−3
Diffusion constant of HCO3 −	Db	9.52 × 10−10	m−2 s−1	(Hoofd et al. 1986)
Cell wall thickness of bundle‐sheath cells	dw_bs	1.5 × 10−7	m	Assumed
Diffusion constant of CO2	Dc	1.83 × 10−9	m−2 s−1	(Hoofd et al. 1986)
Diffusion constant of malate	Dmal	1.22 × 10−9	m−2 s−1	Assumed
Cell wall thickness of mesophyll cells	dw_ms	1.5 × 10−7	m	(Scafaro et al. 2011)
Diffusion constant of OAA	Doaa	1.22 × 10−9	m−2 s−1	(Yaws 1995)
Diffusion constant of PEP	Dpep	1.12 × 10−9	m−2 s−1	Assumed
Diffusion constant of pyruvate	Dpyr	1.12 × 10−9	m−2 s−1	(Yaws 1995)
Fraction of energy partitioning for PCR and PCO cycle in bundle‐sheath	fJ_bs	0.185		Assumed
Fraction of energy partitioning for PCR and PCO cycle in mesophyll	fJ,c_ms	0.63		Assumed
Fraction of energy partitioning for C4 cycle regeneration by PPDK in mesophyll chloroplast	fJ,c4_ms	0.185		Assumed
Fraction of Rubisco partitioning in bundle‐sheath chloroplast	fV_bs	0.185		Assumed
Fraction of Rubisco partitioning in mesophyll chloroplast	fV_ms	0.815		Assumed
Mesophyll cell wall and plasmalemma conductance	Gwall	0.1	mol m−2s−1	Assumed
Maximum electron transport rate per unit leaf area	Jmax	1.6 × 10−4	mol m−2 s−1	(Gu et al. 2012)
Carbonic anhydrase turnover rate	kca	3 × 105	s−1	(Pocker & Ng 1973)
Michaelis–Menten constant for Rubisco carboxylase	Kc	239	μbar	(von Caemmerer 2000)
Equilibrium constant for NADP‐MDH	Ke_mdh	4.45 × 103		(Laisk & Edwards 2000)
Equilibrium constant for NADP‐ME	Ke_me	0.051	mol m−3	(Harary et al. 1953)
Equilibrium constant for hydration	Ke_ca	5.6 × 107		(Pocker & Miksch 1978)
Inhibition constant of malate for PEPC	Ki,mal_pepc	0.5	mol m−3	(Gao & Woo 1996)
Inhibition constant of PEP for PPDK	Ki,pep_ppdk	0.16	mol m−3	(Kanai & Edwards 1999)
Effective Michaelis–Menten constant for Rubisco	Km	14.05 × 10−3	mol m−3	Calculated
Michaelis–Menten constant of NADP‐MDH for malate	Km,mal_mdh	32	mol m−3	(Kagawa & Bruno 1988)
Michaelis–Menten constant of CA for bicarbonate	Km,b_ca	34	mol m−3	(Pocker & Miksch 1978)
Michaelis–Menten constant of CA for CO2	Km,c_ca	1.5	mol m−3	(Pocker & Ng 1973)
Michaelis–Menten constant of NADP‐ME for CO2	Km,c_me	1.1	mol m−3	(Jenkins et al. 1987)
Michaelis–Menten constant of PEPC for HCO3 −	Km,b_pepc	0.02	mol m−3	(Uedan & Sugiyama 1976)
Michaelis–Menten constant of NADP‐ME for malate	Km,mal_me	0.23	mol m−3	(Detarsio et al. 2003)
Michaelis–Menten constant of NADP‐MDH for OAA	Km,oaa_mdh	0.056	mol m−3	(Kagawa & Bruno 1988)
Michaelis–Menten constant of PEPC for PEP	Km,pep	0.1	mol m−3	(Mukerji 1977)
Michaelis–Menten constant of NADP‐ME for pyruvate	Km,pyr_me	3	mol m−3	(Detarsio et al. 2003)
Michaelis–Menten constant of PPDK for pyruvate	Km,pyr_ppdk	0.082	mol m−3	(Jenkins & Hatch 1985)
Effective Michaelis–Menten constant of PEPC for CO2	Kp	2.6 × 10−3	mol m−3	(von Caemmerer 2000)
Michaelis–Menten constant of Rubisco oxygenase	Ko	266	mbar	(von Caemmerer 2000)
Length of plasmodesmata	Lpd	0.2	μm	Assumed
Chloroplast viscosity	ƞch	10		(Tholen and Zhu, 2011)
Cytosol viscosity	ƞcy	2		(Tholen & Zhu 2011)
Mitochondria viscosity	ƞmi	10		(Tholen & Zhu 2011)
Vacuole viscosity	ƞv	1		Assumed
Air pressure	P	105	Pa	Assumed
The fraction of plasmodesmata surface area relative to the total bundle‐sheath cell/mesophyll cell interface area	ɸ	0.03		Assumed
CO2 permeability in chloroplast membranes	Pc,ch	0.0035	m s−1	(Evans et al. 2009)
CO2 permeability in mitochondria membranes	Pc,mi	0.0035	m s−1	(Evans et al. 2009)
Cytosol pH	pHcy	7.3		(Tholen & Zhu 2011)
HCO3 − permeability chloroplast membranes	Pb,ch	5 × 10−7	m s−1	(Felle & Bertl 1986)
HCO3 − permeability mitochondria membranes	Pb,mi	5 × 10−7	m s−1	(Felle & Bertl 1986)
Mitochondria pH	pHmi	8.0		(Tholen & Zhu 2011)
Stroma pH	pHch	8.0		(Tholen & Zhu 2011)
Effective bundle‐sheath cell wall porosity	ptbs	0.1		(Evans et al, 2009)
Effective mesophyll cell wall porosity	ptms	0.2		(Evans et al. 2009)
Dark respiration	Rd	4 × 10−7	mol m−2 s−1	(von Caemmerer 2000)
CO2 solubility	sc	3.29 × 10−4	mol m−3 Pa−1
Mesophyll surface exposed to intercellular spaces area per unit leaf area	Smes	10.04	m2 m−2	(Giuliani et al. 2013), (Hanba et al. 2004)
The fraction of bundle‐sheath and mesophyll cell interface cell wall area per unit leaf area	Sw/Sl	0.8		Assumed
Maximum carboxylation rate of Rubisco per unit leaf area	Vcmax	80	μmol m−2 s−1	(Gu et al. 2012)
Maximum NADP‐MDH catalysed activity per unit leaf area	Vmdh	90	μmol m−2 s−1	(Kanai & Edwards 1999)
Maximum NADP‐ME catalysed activity per unit leaf area	Vme	90	μmol m−2 s−1	(Jenkins et al. 1987), (Kanai and Edwards, 1999) with modification
Maximum PEP carboxylase activity per unit leaf area	Vpepc	120	μmol m−2 s−1	(von Caemmerer 2000)
Maximum PPDK catalysed activity per unit leaf area	Vppdk	90	μmol m−2 s−1	(Kanai & Edwards 1999)
CA concentration in cytosol	Xca,,cy	0.5	mol m−3	(Rumeau et al. 1996)
CA concentration in stroma	Xca,,ch	0.3	mol m−3	(Atkins et al. 1972a,1972b)
CO2 compensation point in the absence of respiration	Γ *	1.56 × 10−3	mol m−3	(von Caemmerer 2000)
Empirical curvature factor	θ	0.7		(Evans 1989; von Caemmerer 2000)
Absorptance of leaves	α	0.85		Assumed
Fraction of absorbed photons that do not drive electron generation	f	0.15		(Evans 1987)

Extended C3–C4 reaction diffusion model

In the extended C3–C4 reaction diffusion model, to represent the scenario of engineering a C4 pathway into a complete C3 leaf without altering the leaf anatomy and biochemical properties, the calculation of electron transport rate based on light intensity and maximum electron transport rate was assumed to be the same as those in C3 reaction diffusion model. Thus, the Calvin–Benson–Bassham cycle was assumed to be ATP limited. Equations for CO2 fixation in chloroplasts limited by ATP were assumed to be the same as those in C3 reaction diffusion model. The C4‐cycle has an additional demand for ATP by the reaction catalysed by PPDK, and the rate of the reaction from pyruvate to PEP (v mol m−3 s−1) limited by ATP was described as follows: where f is fraction of the energy partitioned to C4 cycle. The calculation for CO2 fixation, photorespiration, respiration and hydration are the same as those in C3 reaction diffusion model. Further equations for C4 cycle are described in Supplemental File S2.

Model parameterization, algorithm for solving the models

Model parameterization

The range of mesophyll surface exposed to intercellular spaces area per unit leaf area (S) in rice varies between 10 and 24 m2 m−2 (Hanba et al. 2004; Giuliani et al. 2013). Our model represents a rice leaf with an S of 10 (m2 m−2). Parameters that are related to biochemical and physical processes were taken from Tholen & Zhu (2011). C4 photosynthesis related enzyme kinetic properties were based on Wang et al. (2014b). Plasmodesmata properties and C4 acid diffusion coefficients are given in Table 2. For the C3 reaction‐diffusion model, we assumed that the partitioning of the ATP between mesophyll and bundle‐sheath is linked to the relative volume of chloroplasts in these cells (4.4:1), and the sum of energy fractions allocated to the mesophyll and the bundle‐sheath equals 1 (f + f = 1). For the extended C3–C4 reaction‐diffusion model, we assumed that the mesophyll chloroplasts had the same amount of energy available as in the C3 case (4.4), but a fraction of this energy has to be allocated to PEP regeneration; we assumed that this fraction was balanced with the demands of the C4 cycle (1:1). This leads to an energy partitioning of 3.4:1:1 for the PCR and PCO cycle in the mesophyll, PEPC regeneration in the mesophyll and the PCR and PCO cycle in the bundle‐sheath (f + f + f = 1, f, f, f defaults are given in Table 2).

Algorithm for solving reaction diffusion models

Each reaction diffusion model was discretized into 812 424 elements. Biochemical reactions and boundary conditions (Supplemental File S3) were set up for various metabolites and applied to each subdomain. The models were solved using the finite element method by COMSOL Multiphysics, version 4.3b by a time‐dependent solver (start at time was zero and end time was 106 s). Solving the model resulted in an estimate of the steady‐state concentrations of metabolites at each discrete element of the model structure. Calculation of photosynthetic rate, mesophyll conductance and bundle‐sheath conductance are described in Supplemental File S3.

Results

A comparison of reaction diffusion models with classical biochemical models

We simulated photosynthetic CO2 uptake rate (A) at different intercellular CO2 concentrations (A–C response curve) for a C3 rice leaf and a C3 rice leaf expressing a C4 metabolic cycle using our C3, C3–C4 and the extended C3–C4 reaction diffusion models, and compared the results with the commonly used biochemical models of photosynthesis (Farquhar et al. 1980; von Caemmerer 2000). The simulation results from the C3 reaction diffusion model were compared with the classical C3 biochemical model. The photosynthetic rate in C3 biochemical model was calculated as the sum of photosynthesis in mesophyll and bundle‐sheath cell (Supplemental File S4). The C3 biochemical model predicted slightly higher rates of photosynthesis compared to the C3 reaction diffusion model (Fig. 3a) using the default parameterization for a typical rice leaf (Table 2, Table S1).

Figure 3

Comparison of reaction diffusion models with commonly used biochemical models. (a) Predicted photosynthetic CO2 uptake rate (A) versus intercellular CO2 partial pressure (C) in the C3 reaction diffusion model and the classical C3 biochemical model by Farquhar et al. (1980) (Supplemental file S1) under saturating light. (b) A–C response curves for the C3–C4 biochemical model (Supplemental file S1), the C3–C4 reaction diffusion model and the extended C3–C4 reaction diffusion model under saturating light. (c) A–C response curves for different models under saturating light. (d) Predicted photosynthetic CO2 uptake rate (A) versus light intensity (PPFD) in different models. Intercellular CO2 partial pressure equaled 28 Pa. (e) Predicted net CO2 fixation rates in the mesophyll cell (MSC) and bundle‐sheath cell (BSC) at different intercellular CO2 partial pressures (C) for the C3 reaction diffusion model under saturating light. (f) Predicted net CO2 fixation rates by C4 photosynthesis and C3 photosynthesis at different intercellular CO2 partial pressures (C) for the extended C3–C4 reaction diffusion model under saturating light.

We also compared the C3–C4 and the extended C3–C4 reaction diffusion model to a C3–C4 biochemical model. Because a rice leaf expressing a C4 metabolism may still perform some C3 photosynthesis in mesophyll cells, this C3–C4 biochemical model represented an engineered leaf as the sum of C3 photosynthesis in the mesophyll and C4 photosynthesis in both the mesophyll cell and the bundle‐sheath cells (Supplemental File S4). Similar to the C3 model, the C3–C4 reaction diffusion model predicted slightly lower rates of photosynthesis compared to the biochemical model (Fig. 3b). When a more detailed description of C4 cycle enzymes was included in the extended version of the model, photosynthesis rates increased somewhat (Fig. 3b). A comparison at different light intensities also showed only little differences between the biochemical model and reaction diffusion models (Fig. S1).

The models predict that adding C4 reactions to the C3 model increased the predicted rate of photosynthesis when C is lower than about 45 Pa or PPFD is higher than 240 μmol m−2 s−1 (Fig. 3c, d). Bundle‐sheath photosynthesis only contributed little to the total rate of photosynthesis in the C3 model (Fig. 3e). But the net CO2 fixation rate in bundle‐sheath cells by C4 photosynthesis contributed significantly to the total assimilation rate (Fig. 3f). The presence of a working CO2 concentrating mechanism is further supported by the observations that in the extended C3–C4 reaction diffusion model, the CO2 concentration in the bundle‐sheath cells was higher than in the mesophyll (Fig. S2b). In addition, the flux of malate into bundle‐sheath cells and the net CO2 fixation rate in bundle‐sheath cells reached a maximum at low C (Fig. S2b). The later result contrasts with that for the C3 reaction diffusion model, where the CO2 fixation rate in mesophyll and bundle sheath cells increased gradually with increasing C (Fig. S2a). The extended reaction diffusion model also enables predictions of the concentrations of OAA, malate, PEP and pyruvate in different organelles (Fig. S3). In addition, lowering CA activity had a minor effect on the shape of the A–C curve (Fig. S4).

Factors influencing the photosynthetic CO2 uptake rate of an engineered C3–C4 leaf

We further tested how the predicted rates of photosynthesis in the extended C3–C4 reaction diffusion model change in response to several anatomical and biochemical factors. With respect to anatomy, we investigated the effect of the amount of mesophyll surface that was covered by chloroplasts, the effect of bundle‐sheath cell wall thickness and the effect of permeability of the bundle‐sheath chloroplast envelopes. When the amount of mesophyll surface area covered by chloroplasts was modified, total chloroplast and cytosol volumes were kept constant per unit leaf area by increasing the thickness of the chloroplast, and adjusting the vacuole size. In addition, the concentration of all enzymes in these compartments were kept constant with different coverage. We found that decreasing the surface area of the mesophyll covered by chloroplast from 91.7% to 72.8% had a negligible impact on A in the C3–C4 reaction diffusion model, even at low C, although it corresponded to a decrease in mesophyll conductance of nearly 28% (Fig. 4a). A thickening of the bundle‐sheath cell wall, or decreased permeability of bundle‐sheath chloroplast membrane, had only small effects on A, even although the bundle‐sheath conductance significantly decreased (Fig. 4b, c, d).

Figure 4

The effect of the surface area of chloroplasts exposed to intercellular spaces relative to the surface area of mesophyll cells (S), bundle‐sheath cell wall thickness (d) and bundle‐sheath chloroplast envelope permeability to CO2 (P) and bicarbonate P on photosynthesis and conductance. (a) Net photosynthetic CO2 uptake rate (A, continuous line) and mesophyll conductance (g) between intercellular airspaces and the site of initial CO2 fixation (dashed line) versus the mesophyll chloroplast coverage adjacent to intercellular spaces for the extended C3–C4 model. Intercellular CO2 partial pressures were 28 Pa. Five geometries representing different coverages (indicated by the points) were analysed under saturating light. The dashed vertical line indicates the default coverage given in Table 2. (b) The predicted A (continuous line) and bundle‐sheath conductance (dashed line) versus bundle‐sheath cell wall thickness under an intercellular CO2 partial pressure of 28 Pa in the extended C3–C4 reaction diffusion model under saturating light. The dashed vertical line indicates the default wall thickness given in Table 2. (c, d) Predicted A (continuous line) and bundle‐sheath conductance (dashed line) for different bundle‐sheath chloroplast membrane permeabilities to CO2 (c) and to bicarbonate (d) at an intercellular CO2 partial pressure of 28 Pa under saturating light. The dashed vertical lines indicate default permeabilities of the model (Table 2).

In a C3 leaf with engineered C4 metabolism, ATP is used by three processes: the PCR and PCO cycle in the mesophyll (f), the PCR and PCO cycle in the bundle‐sheath (f) and for the regeneration of PEP by PPDK in mesophyll chloroplasts (f). The energy partitioning between the existing C3 metabolic pathway and the introduced C4 cycle depends on the relative amount of C3 and C4 photosynthetic machinery allocated in mesophyll and bundle‐sheath cells. In the engineered C3–C4 leaf, the default energy partitioning f = 63%, f = 18.5% and f = 18.5% (see details in Materials and Methods). We tested the effect of the energy distribution on photosynthetic rate under the default settings given in Table 2. Maximum photosynthetic rates could be achieved when f = 60%, f = 25% and f = 15% (Fig. 5).

Figure 5

The predicted photosynthetic CO2 uptake rate (A) at an intercellular CO2 partial pressure of 28 Pa versus the fractions of energy partitioned to the PCR and PCO cycle in the mesophyll (f) and the bundle‐sheath (f) in the extended C3–C4 reaction diffusion model under saturated light. The interval used during the sensitivity analysis was 0.05 for each parameter ranging from 0 to 1. The sum of f, f and the fraction of energy partitioning for PEP regeneration by PPDK (f) is 1.

In this study, we assumed that the default Rubisco content in mesophyll and bundle‐sheath chloroplasts scales with the volume of the chloroplast (f, f in Table 2), and therefore only 18.5% of the total amount of Rubisco was present in the bundle‐sheath cells in the default C3 leaf. Increasing the proportion of Rubisco partitioned to the bundle‐sheath cells did not increase photosynthetic rates when this energy partitioning was not modified, and even resulted in a decreased photosynthesis when the fraction of Rubisco partitioned to bundle‐sheath chloroplasts exceeded about 40% (Fig. 6).

Figure 6

The predicted photosynthetic CO2 uptake rate (A) versus the fraction of Rubisco partitioned to bundle‐sheath chloroplast (f) at an intercellular CO2 partial pressure of 28 Pa for the extended C3–C4 reaction diffusion model under saturating light. The dashed vertical line indicates the default value of f given in Table 2.

A sensitivity analysis of the energy partitioning was conducted to determine the maximum CO2 assimilation rate under different fractions of Rubisco partitioning to bundle‐sheath chloroplasts (f) (Fig. 7). As expected, increasing the amount of Rubisco in the bundle‐sheath allows for less energy to be allocated to the C3‐cycle in the mesophyll and overall higher A (Fig. 7). Interestingly, when more than half of Rubisco was allocated to bundle‐sheath chloroplasts, maximum photosynthesis was achieved when no energy was partitioned to C3 photosynthesis in mesophyll chloroplasts. In this case, 60% of the energy was partitioned to PCR and PCO cycle in the bundle‐sheath and 40% was partitioned to PEP regeneration by PPDK in mesophyll chloroplasts (Fig. 7).

Figure 7

The maximal photosynthetic rate (A, red line) achievable with an optimal energy partitioning (determined from an analysis as shown in Fig. 5) for different fractions of Rubisco partitioning to bundle‐sheath chloroplasts (f) in the extended C3–C4 reaction diffusion model under saturating light. The energy partitioning (PCR and PCO cycle in the mesophyll (f), PEP regeneration in the mesophyll (f) and PCR and PCO cycle in the bundle‐sheath (f)) required for this optimal rate is also shown.

Discussion

A new modeling framework to explore the consequences of engineering C4 photosynthetic metabolism in a C3 leaf

This study presents a new modeling framework that couples biochemical reactions with cellular structural features and related gas‐diffusion processes in bundle‐sheath and mesophyll cells. Compared to the earlier modeling efforts (Laisk & Edwards 2000; Wang et al. 2014b), our new framework explicitly considers the anatomical features, diffusional processes, in addition to the biochemical processes. A recent two‐dimensional reaction diffusion model of a maize leaf was built to explore the effect of diffusion and biochemistry in C4 plant (Retta et al. 2016). Here, we presented a full 3D model that not only enables study of the intricate interaction between different biochemical and anatomical features, but also enables evaluation of the impact of modifying these different features on the photosynthetic efficiency.

We compared the results from the reaction diffusion model with the classical biochemical model of photosynthesis (Farquhar et al. 1980; von Caemmerer 2000). When a normal C3 rice metabolism was simulated by the models, there were very few differences in the photosynthetic rate between both models (Fig. 3a). These differences were attributed to the more realistic representation of the 3D leaf structure in a reaction diffusion model, as compared to the biochemical model where the resistances of all components (the cell wall, cytosol, chloroplast envelope, etc.) were considered as serial and simply added (Supplemental File S4). In addition, (photo)respiratory CO2 release may result in a variable effective mesophyll conductance (Tholen & Zhu 2011; Tholen et al. 2012), which was not accounted for in the classical biochemical model. Similarly, when adding PEP carboxylation and decarboxylation to the models, a difference between the C3–C4 reaction diffusion model and the classical biochemical model was observed (Fig. 3b). The explanation for these observations was again the difficulty in accurately estimating and describing the diffusion resistances between the different compartments (i.e. the resistance to CO2 leakage) using the biochemical model.

To enable a comprehensive evaluation of the potential factors controlling photosynthetic efficiency in a leaf where C4 photosynthetic metabolism is engineered into a C3 metabolic background, we developed an extended C3–C4 reaction diffusion model where both the 3D anatomical features and the key enzymes involved in the C4 cycle. Although the influence of engineering any particular C4 component on a leaf can be studied directly through a transgenic approach, a mechanistic model can be used to quickly study the expected consequences of genetic engineering. Our extended C3–C4 reaction diffusion model reached a stable solution and can predict commonly observed A–C curve and light response curve (Fig. 3). It is worth noting that the mesophyll conductance in the extended C3–C4 reaction diffusion model (Fig. 4) was in the range of recent estimates for mesophyll conductance in C4 plants (Barbour et al. 2016).

To illustrate the added capacity of this new reaction diffusion model, we first used it to predict the impact of modifying CA on the rate of photosynthesis (Fig. S4). The concentration of cytosolic CA in the model was based on current estimates for C3 plants (see Tholen & Zhu (2011) for a discussion). In our simulation, a reduction of the CA concentration by 98% resulted in only a small decrease (11% at ambient CO2 conditions) in photosynthesis (Fig. S4). This is consistent with a recent report showing that decreasing the activity of CA in maize by about 97% does not influence photosynthesis much under current or elevated CO2 concentrations (Studer et al. 2014). However, it is worth noting as well that photosynthetic rate was drastically decreased when CA is less than 5% of wild type in C4 dicot plants Flaveria bidentis, which has 10‐fold higher CA activity than maize (von Caemmerer et al. 2004; Cousins et al. 2008). This difference in the control of CA over photosynthetic efficiency might be the result of alternative mechanism in C4 monocots evolution (see the discussion in (Studer et al. 2014)).

Engineering C4 metabolism into a C3 leaf can lead to increased photosynthetic rates

When a two‐cell C4‐cycle was added to an existing C3 metabolism, the model predicted that photosynthesis under saturating light conditions would be enhanced until the CO2 partial pressure in the intercellular airspaces rises above 45 Pa (which is far above the levels corresponding to current atmospheric conditions) (Fig. 3c). Furthermore, under ambient CO2 levels, that is, C = 28 Pa, photosynthetic rate was enhanced when light intensity is above 240 μmol m−2 s−1 (Fig. 3d). The higher predicted A was because of the CO2 concentrating mechanism, as demonstrated by the elevated CO2 concentrations in the bundle sheath cell and also a lower CO2 saturating point in the extended C3–C4 reaction diffusion model, as compared to the C3 reaction diffusion model (Fig. S2). Although adding a C4‐cycle into a C3 leaf decreased the CO2 assimilation in mesophyll cells, the total photosynthetic rate increased compared to the default C3 leaf (Fig. 3c, e, f). These results show that expressing a C4 metabolic cycle in a C3 leaf can increase photosynthesis, even without removing the original C3 metabolism or making extensive changes to the mesophyll anatomy. Therefore, although expressing a single‐cell C4‐cycle in the mesophyll is unlikely to improve photosynthesis (von Caemmerer 2003), our results show that compartmentation of the C4 metabolism in a mesophyll and bundle‐sheath part in a C3 leaf allows for higher A. However, it is important to note that here we assumed that every mesophyll cell was connected to a bundle‐sheath cell. The rice mesophyll has much larger numbers of mesophyll cells per bundle‐sheath (Smillie et al. 2012), and this anatomical feature would reduce the efficiency of such a C4 rice plant. Thus, increasing vein density (Tolley et al. 2012) is necessary to achieve higher rates of photosynthesis in such plants.

Factors influencing the photosynthetic efficiency of an engineered C3–C4 leaf

To test whether the above conclusions are robust, we performed a sensitivity analysis for a number of anatomical and biochemical features used in the reaction diffusion model. Firstly, we examined the influence of chloroplast coverage adjacent to intercellular spaces on the conductance to CO2 diffusion and on the rate of photosynthesis. In C3 plants, mesophyll conductance is thought to correlate well with the proportion of chloroplast area exposed to intercellular air spaces (Laisk et al. 1970; Terashima et al. 2006). Rice has a high proportion of chloroplast surface area covering cell wall (Sage & Sage 2009), which might have been the result of a strong selection pressure. In C4 plants, however, high chloroplast coverage may be counter‐productive as it would decrease the amount of cytosol (and thus PEPC) adjacent to the point of CO2 entry. Mesophyll cells of C4 plants commonly have fewer chloroplast, and these chloroplasts are located further from the cell walls (Stata et al. 2014). von Caemmerer (2003) attributed the inefficient C4 photosynthesis in a single C3 cell partly to the large chloroplast surface coverage. Our results show that the conductance to CO2 diffusion between the intercellular airspace and site of fixation increases with the amount of chloroplast coverage for C3 plants (Fig. S5), but indeed decreases with coverage in C4 photosynthesis in C3 plants. However, both effects on assimilation rates were only minor (Fig. 4, S5), suggesting that the chloroplast coverage in the mesophyll would not significantly limit CO2 uptake and might not need to be modified during C4 engineering.

In C4 plants, the efficiency of the C4 cycle is limited by leakage of CO2 or bicarbonate from the bundle‐sheath back into the mesophyll (Farquhar 1983; von Caemmerer 2000; Kromdijk et al. 2008; Bellasio & Griffiths 2014). If leakage is high, more energy needs to be spent to maintain a high CO2 concentration in the bundle‐sheath, resulting in less efficient photosynthesis. Given this importance of leakage to C4 efficiency, many C4 plants have anatomical features that prevent excessive leakage of carbon from the bundle‐sheath, such as suberized lamellae between the bundle‐sheath and mesophyll cell walls, or centripetally arranged chloroplasts in the bundle‐sheath cells (Leegood 2002). In addition, the diffusion across cell and chloroplast membranes is facilitated by aquaporins acting as a CO2‐pore (Kaldenhoff 2012), and this opens possibilities for differences in membrane permeability or its regulation between C3 and C4 plants. Single‐cell C4 plants similarly have anatomical adaptations that spatially separate the initial CO2 fixation reaction from the decarboxylation reaction near Rubisco (Edwards et al. 2004).

To investigate whether features described above influence the efficiency of a C3 plant expressing a C4 cycle, we varied the thickness of the cell wall, and the permeability of the bundle‐sheath chloroplast membranes to CO2 and bicarbonate. The results show that although an increased permeability of the wall or membranes increased the conductance to carbon, the effect on photosynthesis was minor (Fig. 4). Thus, our results indicate that a mixed C3–C4 metabolism does not benefit much from a high resistance between bundle‐sheath and mesophyll. This is because in the engineered C3–C4 leaf, the CO2 uptake rate contributed by the C4 photosynthesis only accounts for minor (about one third) of the total photosynthetic rate when C equaled 28 Pa (Fig. 3f). Effects of the bundle‐sheath cell wall thickness and bundle‐sheath chloroplast permeability for bicarbonate on the photosynthesis are more significant when more Rubisco and energy are allocated to the bundle‐sheath (Fig. S6).

An efficient C4 metabolism must optimize the energy allocation between PEP regeneration in the mesophyll and the PCR and PCO cycle in the bundle‐sheath. We examined what would be the optimal energy distribution among PEP regeneration, the PCR and PCO cycle in the mesophyll and in the bundle‐sheath. Our results indicate that if 18.5% of the Rubisco is partitioned in bundle‐sheath cells (default value in Table 2), allocating about 50% of the available energy to the C4‐cycle leads to optimal A (Fig. 7). If the amount of Rubisco in the bundle‐sheath can be increased over that in mesophyll cells, correspondingly more energy is required in C4 photosynthesis, then there is no need to maintain C3 photosynthesis in the mesophyll (Fig. 7). These results suggest that a coordinated partitioning of Calvin‐cycle enzymes and energy is an important aspect for achieving improved rates of photosynthesis with C4 engineering. One caveat in this study is that the transfer of PGA/triose phosphate between bundle sheath cell and mesophyll cell, as discussed in a recent C4 metabolism model (Wang et al. 2014b), was not considered. Therefore, the actual energy partitioned into the bundle sheath cell would be higher than the value shown in Fig. 7 because part of the PGA generated in BSC will be phosphorylated and reduced in the mesophyll cells. We also analysed the optimal energy distribution under different light intensities. Our results indicate that under low light, more energy has to be partitioned to C3 photosynthesis to achieve maximal A (Fig. S7). This is a consequence of C3 photosynthesis needing less ATP compared to C4 photosynthesis.

For each CO2 fixed, C4 photosynthesis requires an additional 2 ATP to maintain the carbon concentrating mechanism. The extra ATP requirement can be met by an increased capacity of cyclic electron transport in NADP‐ME type C4 plants compared to C3 plants (Kanai & Edwards 1999; Nakamura et al. 2013). In our current analysis, we focused on examining the role of diffusion and enzyme limited biochemistry on the photosynthetic rates in an engineered C3–C4 leaf. An analysis of the effect of changes in the electron transport stoichiometry is beyond the scope of our work. Models have been developed that describe the required cyclic electron transport in C3 and C4 leaves (Zhu et al. 2005; Yin & Struik 2012; Walker et al. 2014). However, because of the current incomplete understanding of the regulation of electron transport and the interaction between cyclic and linear electron transfer capacity, a fully mechanistic model of this process is yet to be developed.

Implications for the evolution and engineering of C4 photosynthesis

C4 photosynthesis differs from C3 photosynthesis in many anatomical and biochemical aspects (Sage & Zhu 2011). The current notion is that these features were acquired in a stepwise manner during the C4 evolutionary processes (Sage et al. 2012). Mallmann et al. (2014) suggest that the re‐balancing the nitrogen metabolism between bundle sheath and mesophyll cells after the re‐localization of glycine decarboxylase from mesophyll to bundle sheath cells might have been a major evolutionary driving force for the establishment of the C4 metabolic cycle. The present study suggests that formation of a C4 metabolic cycle in a C3 photosynthetic background can also lead to increased photosynthetic efficiency and hence function as an evolutionary driving force (Fig. 3). This benefit is larger under low CO2 concentrations (Fig. 3), which is consistent with the report that C4 species emerged during the Oligocene period, a geological period that is feature by a low atmospheric CO2 concentration (Christin et al. 2008).

The results obtained from this simulation study are consistent with the current understanding of the evolutionary trajectories of C4 photosynthesis. The metabolic structure of the simulated C3–C4 species in this study mimics C4‐like species, where a C4 photosynthetic metabolism is incorporated into a C3 background without establishment of cell specific expression of photosynthetic enzymes, in particular Rubisco, between the bundle sheath and mesophyll cells (Sage et al. 2012). After the formation of C4‐like species, photosynthetic efficiency is further optimized through establishment of cell specific expression patterns for key photosynthetic enzymes, e.g. Rubisco (Sage et al. 2012). Our results show that further redistribution of Rubisco content between bundle sheath and mesophyll cells can lead to increased photosynthesis, and this is consistent with such an evolutionary sequence (Fig. 7). In addition, we found the percentage of energy required for maximal C4 photosynthesis in the C3–C4 rice leaf increased with more Rubisco distributed to bundle‐sheath cell. Thus, increased energy partitioning to the bundle‐sheath is an expected trend during C4 evolution.

Conclusions

We developed a 3D reaction‐diffusion model by incorporating a NADP‐ME‐type C4 metabolism and a C3 biochemical model in a realistic geometry representing part of a rice leaf. The model was then used to explore the influence of anatomical and biochemical features of an engineered C4 rice leaf. Our results suggest that expressing a two‐cell C4 metabolism in a C3 rice leaf may lead to an increased photosynthetic efficiency. Furthermore, we found that Rubisco allocation and energy partitioning are crucial to gain an increased photosynthesis in the engineered leaf.

Author Contributions

Xin‐Guang Zhu, Danny Tholen and Shuyue Wang designed the research. Shuyue Wang developed the model and performed the analysis. Shuyue Wang, Danny Tholen and Xin‐Guang Zhu wrote the paper.

Conflict of Interest

The authors declare that they have no conflict of interest.

Table S1. Additional anatomical and biochemical parameters in the reaction diffusion and biochemical models (at 25 °C).

Figure S1. Comparison of the predicted photosynthetic CO2 uptake rate (A) versus light intensity (PPFD) in the C3–C4 biochemical model, C3–C4 reaction diffusion model and the extended C3–C4 reaction diffusion model. Intercellular CO2 partial pressure equaled 28 Pa.

Figure S2. Comparison of CO2 concentration between C3 and extended C3–C4 reaction diffusion models. Predicted CO2 partial pressures in mesophyll chloroplast (continuous lines) and bundle‐sheath chloroplast (dashed lines) at different intercellular CO2 partial pressures (C) for C3 (a) and extended C3–C4 reaction diffusion (b) models under saturating light. For the extended C3–C4 model, the flux of malate into the bundle‐sheath is also shown.

Figure S3 Metabolite concentrations (OAA: oxaloacetic acid; PEP: phosphoenolpyruvate; malate and pyruvate) in specific organelles versus intercellular CO2 partial pressure (C) predicted by the extended C3–C4 reaction diffusion model under saturating light. MS: mesophyll cell; BS: bundle‐sheath cell.

Figure S4 Predicted photosynthetic CO2 uptake rate (A) at saturating light versus intercellular CO2 partial pressure (C) in the extended C3–C4 reaction diffusion model at default and at low (2% in mesophyll cell) CA levels.

Figure S5. Net CO2 fixation rate (continuous line) and mesophyll conductance between intercellular airspaces and the site of initial CO2 fixation (dashed line) versus the fraction of the mesophyll surface that is covered by chloroplasts (S/S) for the C3 reaction diffusion model under saturated light. Intercellular CO2 partial pressure was 28 Pa. Five different model geometries (indicated by the points) were analysed. The dashed vertical line indicates the default coverage.

Figure S6. The effect of anatomical features, including surface area of chloroplasts exposed to intercellular spaces relative to surface area of mesophyll cells (S), bundle‐sheath cell wall thickness (d), bundle‐sheath chloroplast envelope permeabilities to CO2 (P) and bicarbonate (P) on photosynthesis (A) and conductance for the extended C3–C4 reaction diffusion model with different fractions of Rubisco partitioned in bundle‐sheath chloroplast (f). (a, b) Net photosynthetic CO2 uptake rate (A) (a) and mesophyll conductance (g) (b) between intercellular airspaces and the site of initial CO2 fixation versus the mesophyll chloroplast coverage adjacent to intercellular spaces. (c, d) The predicted A (c) and bundle‐sheath conductance (d) versus bundle‐sheath cell wall thickness. (e, f) Predicted A (e) and bundle‐sheath conductance (f) for different bundle‐sheath chloroplast membrane permeabilities to CO2. (g, h) Predicted A (g) and bundle‐sheath conductance (h) for different bundle‐sheath chloroplast membrane permeabilities to bicarbonate. Simulations were done at an intercellular CO2 partial pressures as 28 Pa under saturating light. Corresponding energy partitioning with different f was from Fig. 7. Other parameters were taken from Table 2.

Figure S7. The maximal photosynthetic rate (A, red line) achievable with an optimal energy partitioning (determined from an analysis as shown in Fig. 7) for different light levels. The energy partitioning (PCR and PCO cycle in the mesophyll (f), PEP regeneration in the mesophyll (f) and PCR and PCO cycle in the bundle‐sheath (f)) required for achieving these rates is also shown.

File S1: Reaction rates in the C3–C4 reaction diffusion model.

File S2: Mass balance equations for C4 acids and rate reactions of the C4 cycle in the extended C3–C4 reaction diffusion model.

File S3: Boundary conditions and calculation of conductance for reaction diffusion models.

File S4: Biochemical models for C3 and C3–C4 photosynthesis

Supporting info item

Click here for additional data file.