Decomposition analysis on soybean productivity increase under elevated CO 2 using 3-D canopy model reveals synergestic effects of CO 2 and light in photosynthesis PART OF A SPECIAL ISSUE ON FUNCTIONAL-STRUCTURAL PLANT GROWTH MODELLING

• Background and Aims Understanding how climate change influences crop productivity helps in identifying new options to increase crop productivity. Soybean is the most important dicotyledonous seed crop in terms of planting area. Although the impacts of elevated atmospheric [CO 2 ] on soybean physiology, growth and biomass ac- cumulation have been studied extensively, the contribution of different factors to changes in season-long whole crop photosynthetic CO 2 uptake [gross primary productivity (GPP)] under elevated [CO 2 ] have not been fully quantified. • Methods A 3-D canopy model combining canopy 3-D architecture, ray tracing and leaf photosynthesis was built to: (1) study the impacts of elevated [CO 2 ] on soybean GPP across a whole growing season; (2) dissect the contri- bution of different factors to changes in GPP; and (3) determine the extent, if any, of synergism between [CO 2 ] and light on changes in GPP. The model was parameterized from measurements of leaf physiology and canopy architec- tural parameters at the soybean Free Air CO 2 Enrichment (SoyFACE) facility in Champaign, Illinois. • Key Results Using this model, we showed that both a CO 2 fertilization effect and changes in canopy architecture contributed to the large increase in GPP while acclimation in photosynthetic physiological parameters to elevated [CO 2 ] and altered leaf temperature played only a minor role in the changes in GPP. Furthermore, at early developmental stages, elevated [CO 2 ] increased leaf area index which led to increased canopy light absorption and canopy photosynthesis. At later developmental stages, on days with high ambient light levels, the proportion of leaves in a canopy limited by Rubisco carboxylation increased from 12.2 % to 35.6 %, which led to a greater enhancement of elevated [CO 2 ] to GPP. • Conclusions This study develops a new method to dissect the contribution of different factors to responses of crops under climate change. We showed that there is a synergestic effect of CO 2 and light on crop growth under elevated CO 2 conditions.


INTRODUCTION
Soybean is the most important dicotyledonous seed crop in terms of area planted and mass produced at the global scale. It is also the largest single source of vegetable protein for food and feed in the world. Understanding the response of soybean growth and productivity to global atmospheric change will be critical to an accurate projection of future soybean production and global food security (Parry et al., 2004;Long et al., 2006). The impacts of elevated [CO 2 ] on soybean physiology, growth, development and yield have been extensively documented (Kimball, 1983;Ainsworth et al., 2002;Bernacchi et al., 2005Bernacchi et al., , 2006Bernacchi et al., , 2007Morgan et al., 2005). Soybean grown under free air [CO 2 ] enrichment shows a 37 % increase of biomass dry weight and 24 % increase in yield with a [CO 2 ] increase to about 700 ppm (Ainsworth et al., 2002). When [CO 2 ] is increased to about 550 ppm, above-ground biomass increases by 17-18 % and yield increases by 15 % . The yield increase is significant, but less than the expected increase under elevated [CO 2 ] . Similar observations of lower than expected increase in crop yields under elevated [CO 2 ] have been made in other species as well .
To a first approximation, the yield of a crop under optimal growth conditions is the product of the solar radiation receipt for the growing season and the efficiencies with which the crop intercepts that radiation (ε i ), converts it into biomass (ε c ) and partitions the biomass into the harvested organ (ε p ), i.e. seed in the case of soybean (Monteith, 1972). Under elevated [CO 2 ], the partitioning efficiency (also known as harvest index) of soybean is only slightly reduced, compared to under ambient [CO 2 ]. Measured harvest index is about 0.55 and 0.53 for

PART OF A SPECIAL ISSUE ON FUNCTIONAL-STRUCTURAL PLANT GROWTH MODELLING
soybean grown under atmospheric [CO 2 ] of 370 and 550 ppm, respectively . Furthermore, given that the root : shoot ratio changes relatively little for soybean under elevated [CO 2 ] (Ainsworth et al., 2002), the increase in yield and biomass can mainly be attributed to an increase in total CO 2 uptake of the whole canopy, i.e. gross primary productivity (GPP). GPP includes the photosynthetic CO 2 uptake of all leaves in the canopy, including not only those leaves at the top of the canopy which are mostly light-saturated, but also those at lower layers which are mostly light-limited.
GPP is affected by many environmental factors besides light, including [CO 2 ], temperature, humidity and plant water status. It is affected also by canopy size, canopy architecture, leaf chlorophyll content and leaf photosynthetic parameters. When soybean is grown under elevated [CO 2 ] many physiological parameters change Ainsworth & Long, 2005), such as leaf area index (LAI) (Dermody et al., 2006), the maximal rate of carboxylation under RuBP and CO 2 saturation (V cmax ) Morgan et al., 2004), stomatal conductance  and microclimatic parameters especially leaf temperature . These changes differentially influence GPP. Increasing LAI at early crop developmental stages leads to increased GPP, but at later developmental stages an increase in LAI may decrease GPP (Srinivasan et al., 2016). Decreasing V cmax can monotonically decrease GPP, decreasing stomatal conductance (g s ) can similarly monotonically decrease GPP, while changes in leaf temperature (T leaf ) influence GPP in an non-linear manner dependent on ambient temperature. Given the complexity and non-linearity in the impacts of these factors on GPP, it is difficult to dissect the contribution of these individual factors and their interactions on changes in GPP and correspondingly observed increases in biomass and seed production.
In this study, we used a mathematical modelling approach to dissect the contribution of different factors to change in GPP. In fact, many mathematical models of canopy photosynthesis have been developed and used to study GPP. The sunlit-shaded canopy photosynthesis model, which divides leaves in a canopy dynamically into those that are sunlit and those that are shaded, has been used widely in predicting canopy photosynthesis of plants under different conditions (Norman, 1980;Givnish, 1988;DePury & Farquhar, 1997). This method has been used to evaluate the consequences of altering Rubisco kinetic properties and changing the speed of relaxation of photoprotection on GPP (Zhu et al., 2004a, b). It has also been used to evaluate the potential changes in GPP under elevated [CO 2 ] and to dissect the factors controlling changes in GPP (Wittig et al., 2005). The sunlit-shaded model uses aggregated parameters to represent canopy architectural features. Recently, 3-D canopy models were developed which use canopy architectural features directly and can be used to predict fine details of canopy light distribuiton (Song et al., 2013;Kim et al., 2016). Such models can be used to study the impacts of changes in architecture on canopy photosynthesis. For example, these models have been used to evaluate the optimal architectures that maximize canopy photosynthetic CO 2 uptake rate under different ambient [CO 2 ] (Song et al., 2013), and to predict the impacts of assymetric row spacing and row orientation on GPP in sugarcane agronomy (Wang et al., 2017).
In this study, we use a 3-D canopy photosynthesis model to dissect the contribution of different environmental, architectural and physiological parameters to the changes in GPP under elevated CO 2 . Using the model, we further show that there is a synergistic effect between [CO 2 ] and light for soybean grown under elevated CO 2 .

SoyFACE facility
The data used for model parameterization in this study were collected from the SoyFACE facility (www.soyface.uiuc.edu), which employs a free air concentration enrichment (FACE) technology. SoyFACE is a 32-ha facility at the University of Illinois at Urbana-Champaign (40°03′21.3″N, 88°12′3.4″W, 230 m elevation). The soil at SoyFACE is a Drummer-Flanagan series (fine-silty, mixed, mesic Typic Endoaquoll; Morgan et al., 2005), typically very deep and formed from loess and silt parent material deposited on the till and outwash plains. The average ground surface slope is <1 % at this site, with tile drains at a depth of 1-2 m below the ground surface. This rain-fed field site has been following a continuous crop rotation practice typical for the US Midwestern corn belt of soybean and maize. Extended descriptions of site, including micrometeorology and climate, have been previously described (Leakey et al., 2004;Rogers et al., 2004). SoyFACE elevated [CO 2 ] treatments consist of four blocks each with two octagonal plots of 20 m diameter within the 16 ha planted with soybean. Each block contains one control plot at an ambient [CO 2 ] of 370 ppm and one fumigated plot to a target [CO 2 ] of 550 ppm, using the FACE technology of Miglietta et al. (2001). Fumigation was performed between sunrise and sunset and began 3 d after planting and operated over the remainder of the growing season until crop harvest.
Leaf photosynthetic parameters V cmax and J max of canopy top leaves for ambient and elevated [CO 2 ] were taken from Morgan et al. (2005). An exponential distribution model was used to predict V cmax and J max for leaves at different depth (z) of the canopy (eqn 1) following previous observations of V cmax , J max and leaf nitrogen content in different layers of tghe canopy (Morgan et al., 2004;Srinivasan et al., 2016). The relationship between V cmax and cumulative leaf area index (cLAI) from the top of the canopy was based on the vertical distribution of leaf nitrogen content in the canopy. V cmax,top is the V cmax of leaves in the top layer of a canopy. (1)

Meteorological parameters
Air temperature (°C), relative humidity, photosynthetic photon flux density (PPFD) and diffuse PPFD were recorded with 10-min intervals at the weather station located at SoyFACE.

Calculation of GPP
The overall workflow used to calculate GPP was as follows: (1) development of a 3-D soybean canopy photosynthesis model; (2) simulation of the light environments inside a canopy; (3) calculation of the photosynthetic rate of each leaf; and (4) calculation of the canopy photosynthetic rate by integrating the total photosynthetic rate for all leaves in a canopy. These steps are described in detail in the following sections.

3-D soybean canopy model
A 3-D soybean canopy model was developed based on the 3-D canopy modelling algorithm of Song et al. (2013) (the MATLAB-based program, mCanopy-soybean, is available from the authors upon request).
The different measurements required for the 3-D canopy model development are summarized in Supplementary Data Tables S1-S8. Leaf lengths, leaf widths, petiole lengths and angles for the left, middle and right trifoliate leaves (Tables S1-S3) were obtained using digital photography and image processing. An in-situ non-destructive leaf photo scanner was used to obtain a picture of the soybean trifoliate leaves at all nodes for five plants per plot. The scanner holds a leaf between two flat sheets, with the top sheet being transparent, to obtain a photo of the leaf. A standard length was placed within the frame of the picture for reference. The camera was mounted normal to the leaf plane and the leaf image was recorded using a digital camera. Leaf widths, lengths and the angles were digitally processed using the javabased image processing and analysis software ImageJ (1.48j 11 December 2013, https://imagej.nih.gov/ij/).
The internode lengths and petiole lengths (Supplementary Data Table S4) for each node in the main stem and the branches were measured for the same five plants in each plot using a ruler. Petiole angle, branch angle and leaf angle were measured using a protractor (iGaging digital protractor, http://www. igaging.com/) at midday. These measurements were made for the same five plants in each plot.
Measurements of leaf senescence were made on the same five plants in each plot, recorded twice a week for each node. Leaves that turned pale yellow to brown were considered nonphotosynthetic and senesced. Note that there was a finite duration of time between post-yellowing of leaves and complete litterfall (a few days to a week). In our 3-D model, we assume all leaves that turned yellow have senesced. During the growing season, ten plants from each plot were identified to measure development, and the presence of branches and their node lengths were measured. Using these data, a branching probability was computed (Supplementary Data Tables S5-S7).
A soybean architecture model (features in Fig. 1A-C) was developed based on these measured architecture data. The description of nodes and branches is given in Table 1. In the main stem, the base node before the first node is the VU node, which is the node growing unifoliate leaves (Fig. 1A), and the first trifoliate leaf is on the first node V1. The second trifoliate leaf grows on the second node V2 and so on. The first branch (Br1) grows from the V1 node and the second branch (Br2) from the V2 node and so on. Internode length is the distance between two nodes (e.g. V3 internode length was the length between nodes V3 and V2). Branch angle is the angle between the main stem and branch (Fig. 1B), petiole angle 1 is the angle between the main stem (or branch) and the common petiole of a trifoliate leaf (Fig. 1C), and petiole angle 2 (Fig. 1B) is the angle between a common petiole and the petiole of the mid-leaf of the trifoliate leaf. Mid-leaf angle is the angle between the petioles of the mid-leaf and the main vein of the mid-leaf (Fig. 1b). Leaf length and leaf width are the maximal length and width of a leaf. Leaf angle L, leaf angle R and leaf angle M (Fig. 1C) are the angles between the common petiole of a trifoliate leaf and the main veins of left, right and middle leaves when these leaves are laid on a horizontal plane.

3-D soybean canopy models applied to different days in the growing season
Row spacing was 38 cm and the planting distance was 5 cm under both ambient and elevated CO 2 conditions. The integrated model was run for the year 2002 between days 168 (V1) and 267 (V16) every 3 d (Fig. 1D, E).
To simulate canopy photosynthesis throughout a growing season ( Fig. 1), canopy architectural parameters (i.e. leaf length, leaf width, leaf angle, internode length, etc.) were measured at different stages of plant growth (Supplementary Data Tables S1-S4). To model the variation of node number among different plants, we used a randomization algorithm to determine the node number (V x ) for each main stem and branch as follows.
First, the maximal node number (V x_max ) for the main stem and branches for different days of the year (DOY) were determined based on previous measurements (Castro et al., 2009). Second, the probability p(n) (probability for node number = n) of the main stem or branches on different days were calculated based on measurement data collected in the field (Supplementary Data Table S5). A random value i between 0 and 1 (uniform distribution) was then generated and if i > sum(p(n<N)) and i < sum(p(n>N+1)), N is used for V x_r , which is the randomized node number. Lastly, if V x_r is less than V x_max , V x_r is used as the node number, and if V x_r is larger than V x_max , V x_max is used as the node number (Tables S5-S7). Pseudocode for the above process is given in the Supplementary data Methods.
Senescenced leaves were excluded from the model. The number of senescenced leaves was counted during the growing season and the number of senescenced leaves every 3 d was based on measurement data (Supplementary Data Table S8). The total number of senescenced leaves (N s ) on a DOY is the sum of the senescenced leaf number from the start day (168 DOY) to the current day. If N s is not an integer, a randomization algorithm is used to determine the total senescenced leaf number. For example, if N s = 2.3, then leaf number being 2 is used with a probability of 0.7 and leaf number being 3 is used with a probability of 0.3. In the model, newly formed leaves were smaller than mature leaves. The sizes of the top three newest formed leaves (from top to bottom) were assumed to be 25 %, 50 % and 75 % of their mature sizes.

Ray tracing algorithm
The light environment in the soybean canopy was simulated using a ray tracing algorithm, fastTracer (Song et al., 2013).
The program for fastTracer is available upon request from the corresponding author. Measured ambient direct and diffuse photosynthetic photon flux density (PPFD) data were used as input to fastTracer. The ray tracing algorithm simulates the direct PPFD, diffuse PPFD and leaf scattering PPFD absorbed by every leaf. The simulations were conducted 15 times a day, from 0500 to 1900 h with intervals of 1 h. Using ray tracing, the calculation of PPFD distribution in the canopy was accurate for evaluating the impacts of small differences in canopy architecture, such as leaf size and leaf number and the PPFD distribution under different weather conditions, both of which were important for this study.

Leaf photosynthesis calculation
Leaf photosynthetic CO 2 uptake rate was calculated with the steady-state biochemical model of C 3 leaf photosynthesis (Farquhar et al., 1980).
Leaf photosynthesis rate at any given CO 2 , light, O 2 and temperature conditions was calculated from eqn (2): where: Γ* is the CO 2 compensation point in the absence of dark respiration, C i is the leaf intercellular CO 2 concentration, W c is   The first node Fully developed trifoliate leaf at the node above the unifoliate node V2 The second node Two nodes on the main stem with fully developed trifoliate leaves V(n) The nth node n nodes on the main stem with fully developed trifoliate leaves BV1 The first node of a branch First fully developed trifoliate leaf at the first node of a branch BV2 The second node of a branch Two nodes on a branch with fully developed trifoliate leaves BV(n) The nth node of a branch n nodes on a branch with fully developed trifoliate leaves Br1 The first branch Branch developed on the main stem at V1 Br2 The second branch Branch developed on the main stem at V2 Br(n) The nth branch Branch developed on the main stem at V(n) the Rubisco-limited rate of carboxylation and W j is the RuBP regeneration-limited rate of carboxylation, which were calculated from eqns (3) and (4): (3) where K c is the Michaelis constant for CO 2 (404.9 μmol mol −1 ), K o is the Michaelis constant for O 2 (278.4 mmol mol −1 ) and J is the potential photosynthetic electron transport rate. The parameters describing impacts of temperature on Rubisco kinetics follow Bernacchi et al. (2001). Leaf temperature was assumed to be equal to air temperature T for the ambient condition and leaf temperature was assumed to be 1.5 °C higher under elevated CO 2 than under ambient conditions based on the data from Long et al. (2006). The potential electron transport rate, J, was calculated as: where Θ PSII (0.864 was used) is the convexity of the nonrectangular curve, I 2 is the PPFD absorbed by photosystem II (PSII), and J max is the maximal electron transport rate (Chen, Zhu & Long, 2008). I 2 was calculated from: where I sun is the PPFD incident upon a facet in a leaf simulated by the ray tracing algorithm described by Song et al. (2013), α l is leaf absorbance (0.85 was used), Φ PSII,max is the maximal quantum yield of PSII (0.85 was used) and β (0.5 was used) is the maximal fraction of quanta that reaches PSII (Chen et al., 2008). The parameters used in describing the temperature response were taken from .

Estimation of Ci and gs
During the calculation of leaf photosynthesis, C i and g s were estimated using eqns (7)-(12) as in previous studies (Humphries & Long, 1995;Song et al., 2013) based on Ball et al. (1987). Equation (12) shows the calculation of intercellular CO 2 partial pressure (C i , µbar) based on the CO 2 partial pressure on the leaf surface (C s , µbar), photosynthetic CO 2 assimilation rate (A, µmol m −2 s −1 ), stomatal conductance (g s , mmol m −2 s −1 ) and air pressure (P a , bar). Equation (11) shows the calculation of CO 2 partial pressure on the leaf surface (C s ) based on the ambient CO 2 partial pressure (C a , µbar), photosynthetic CO 2 assimilation rate (A), leaf boundary conductance (g b , mmol m −2 s −1 ) and air pressure (P a ). The parameters a, b and c calculated in eqns (8)-(10) are those used in eqn (7), which is an empirical equation used to calculate g s based on the CO 2 partial pressure on the leaf surface (C s ), photosynthetic CO 2 assimilation rate (A), leaf boundary conductance (g b ), relative humidity (RH), partial pressure of the saturated water vapour for the air temperature (e air , mbar) and for leaf temperature (e leaf , mbar), stomatal coefficient g 0 (20) and stomatal coefficient g 1 (11.35). The equations used in this study to calculate stomatal conductance are suitable to model plants without water stress (Ball et al., 1987). To adapt the Ball-Berry model to simulate plants under water stress conditions, additional factors associated with leaf water potential or soil water content are needed (Buckley et al., 2003;Li et al., 2012).
Temperature response of photosynthetic parameters Leaf temperature affects photosynthetic parameters and eqns (13)-(21) were used to describe the temperature responses of photosynthetic parameters.
Iterative calculation of Ci, gs and P under constant leaf temperature Equations (1)-(21) together form a system of equations, which was solved by using the Newton-Raphson method. During this computation, the initial C i was calculated (eqn 15) and then used to calculate P (eqn 2) and g s (eqn 7); C i was updated using eqn (12). The calculation was terminated when either C i reached a stable value or the maximum number of iterations was reached.

Calculation of GPP
GPP was calculated by integrating the photosynthesis rates of all leaves in a canopy (eqn 22).
where P i is leaf photosynthesis rate for the ith leaf facet (the 'i' here is the sequence ID of a facet in the data file, 'leaf ID' is used to label a leaf for parameterization of leaf photosynthetic parameters) and S i is the corresponding leaf area of the leaf facet. S ground is the ground area occupied by the canopy.
To study GPP on different days during the growing season and to compare the whole season GPP of soybean under ambient and elevated CO 2 conditions, the GPP per day or per season was calculated by integrating GPP during a day or throughout the growing season (eqn 23).
where GPP day is the GPP per day, P c,t is total canopy photosynthetic CO 2 uptake for a unit ground area during a particular hour; P c,u is total canopy photosynthetic CO 2 uptake at the midpoint of an hour, i.e. at the 30th minute in each hour.
Cumulative GPP (cGPP) from the start day (168 DOY) to different days (n) between 168 DOY and 267 DOY was calculated by adding GPP day from the start day to day n (eqn 24).
Calculation of GPP under ambient and elevated CO 2 The model was parameterized for both ambient and elevated CO 2 conditions ( Table 2). The parameters used include [CO 2 ], V cmax , J max , air temperature and canopy structure on different days from 168 DOY to 267 DOY (for each day, the model was run three times).
Estimation of above-ground biomass from cGPP Model-estimated above-ground biomass was calculated based on cGPP, the harvest index, carbon content in different organs and root : shoot ratio. Briefly, we first calculated the total carbon in biomass (W c ) by subtracting respiration from cGPP assuming that respiration is 57.5 % of total photosynthetic CO 2 uptake (Amthor, 2000). Second, the proportion of carbon for the whole plant was calculated based on the harvest index, the proportion of carbon in pod and seed, and the proportion of carbon in other parts of the soybean, which was assumed as (C 6 H 10 O 5 ) n . We assumed that the composition of soybean pod and seed was carbohydrate (29 %), protein (37 %), lipid (18 %), lignin (6 %), organic acid (5 %) and mineral (5 %); this results in a proportion of carbon (C pod+seed ) of 53 % in pod and seed (Amthor, 2000). We also assumed that the composition of root, stem and leaf biomass was (C 6 H 10 O 5 ) n , which results in a proportion of carbon (C other ) being 44.4 %. We further assumed a harvest index η of 0.57 (Pedersen & Lauer, 2004;Spaeth et al., 2010), and then calculated the proportion of carbon in a whole soybean plant (C plant ) as: Whole plant biomass was given by: Finally, assuming a ratio of root to total biomass (p root ) as 18.7 % (Clough & Peet, 1981), we calculated above-ground biomass: Dissecting the factors contributing to changes in GPP under elevated CO 2 To evaluate the relative contribution of physiological parameters (i.e. V cmax and J max ), architectural parameters and environmental factors (i.e. air temperature and CO 2 ) to the changed GPP under elevated CO 2 , we simulated GPP under different scenarios. The method was adapted from the sensitivity analysis of a model that is commonly used in previous studies (Zhu et al., 2007;Wu & Cournède, 2010). The different scenarios are listed in Table 3. The contribution of all four factors, i.e.  was calculated from eqns (29) Statistical analysis Pearson's correlation coefficient was calculated with the R software function cor. Student's t-test was calculated with the R software function t.test.

Canopy architectural and physiological data used to develop soybean canopy models
Soybean grown under elevated [CO 2 ] revealed about a 10-90 % increase in leaf lengths and leaf widths for different leaves (Supplementary Data Tables S1 & S2). Leaf angle distributions were assumed to be the same for both CO 2 treatments based on field observations (Table S3). Internode distances were obtained from direct measurements under ambient [CO 2 ] conditions (Table S4). Internode distances for soybean grown under elevated [CO 2 ] were 6 % longer than under ambient [CO 2 ] conditions (Ainsworth et al., 2002). The probabilities of node numbers of the main stem and branches at the final stage were assumed to be the same between the two [CO 2 ] conditions and were measured (Table S5). For the main stem, all measured plants had more than ten nodes; 95 % of plants had more than 11 nodes; with an increase in node number, the percentage of plants having more than the node number gradually decreased (Table S5). Node probabilities differed dramatically between branches. For example, the probability of having more than one node was 10 % for Br1, 50 % for Br2 and 30 % for Br3 (Table S5). Maximal node numbers of the main stem and different branches during a growing season under ambient and elevated [CO 2 ] conditions are given in Tables S6 and S7, which was compiled based on a previous study (Castro et al., 2009 (Table S8). Leaf photosynthetic parameters, V cmax and J max , collected based on Bernacchi et al. (2005), are shown in Table S9. GPP during a growing season under elevated and ambient CO 2 conditions Simulated GPP day was higher under elevated [CO 2 ] compared to ambient [CO 2 ] throughout the growing season ( Fig. 2A), with the relative increase in simulated GPP day being higher in the early growth season (DOY<200) (Fig. 2B). The seasonal trends in simulated daily GPP (GPP day ) mimicked the behaviour of measured LAI with an initial fast increase, followed by a peak in GPP day near the date of peak LAI and a final decline towards the end of the growing season ( Fig. 2A, C). Simulated cGPP showed a high degree of correlation (Pearson coefficient > 0.99) with measured aBM , under both ambient and elevated [CO 2 ] conditions ( Fig. 2D) Table S10). To compare the model simulation with measured data, we further estimated the above-ground biomass with cGPP. The model-estimated above-ground biomass was linearly correlated with measured values, with slope being 0.96 and 0.97 under elevated and ambient CO 2 , respectively (Fig.  2E). The R 2 values for these two relationships were 0.995 and 0.992 for elevated and ambient CO 2 respectively (Fig. 2E).

(Supplementary Data
Contributions of CO 2 concentration, canopy structure, leaf temperature and V cmax and J max to the increase of GPP under elevated [CO 2 ] Between soybean plants grown under elevated [CO 2 ] and ambient [CO 2 ], four factors (i.e. CO 2 concentration, temperature, V cmax and J max , and canopy structure, which includes leaf size and leaf number) differed (Table 2). Here we used canopy photosynthesis models to dissect the contribution of each of these Table 3

. Scenarios used to calculate GPP, which is used to dissect the contributions of individual factors ([CO 2 ], V cmax and J max , T, canopy structure) and their interactions to changes in GPP
Scenarios Elevated CO 2 (X), Ambient CO 2 (-)  (Fig. 3). Figure 3A shows the averaged contributions over the growing season of CO 2 (76.7 %), canopy structure (17.2 %) and their interaction (2.6 %) to increases in ΔGPP (Fig. 3A). The increase in air temperature under elevated CO 2 showed a mild negative impact (−0.1 %) on ΔGPP, but the interaction between CO 2 and temperature showed a substantial positive impact (6.9 %) on ΔGPP (Fig. 3A). V cmax and J max negatively influenced ΔGPP, i.e. decreasing GPP by −6.7 % of ΔGPP, but their interactions with CO 2 had a positive impact of 3.3 % on ΔGPP (Fig. 3A). The contributions of other interactions on ΔGPP are not significant (<1 %) (Fig. 3A).
To further investigate the impacts of these factors on ΔGPP in different growth stages and weather conditions, we chose a Sunny day in the Early developmental stage (SE), a Cloudy day in the Late developmental stage (CL), and a Sunny day at a Later developmental stage (SL) for analysis ( Fig. 3B-D). The contribution of CO 2 concentration was smaller in SE (62.5 %) than in SL (99.8 %) and the contribution of canopy structure was larger in SE (20 %) than in SL (4.8 %) (Fig. 3B, D). The contribution of [CO 2 ] was smaller in CL (93.8 %) than in SL (99.8 %) and the contribution of canopy structure was greater in CL (6 %) than in SL (4.8 %) (Fig. 3C, D).  (Dermody et al., 2006) for soybean grown under elevated [CO 2 ] and ambient [CO 2 ] conditions. (D) Correlation between above-ground biomass  and calculated cumulative GPP (cGPP) at different stages. The R 2 of linear fitting was > 0.99 for both ambient [CO 2 ] and elevated [CO 2 ] conditions. (E) Correlation between the measured above-ground biomass and model-estimated above-ground biomass from cGPP.

Canopy absorbance and PPFD distribution in different stages and weather conditions
To explore the greater relative contribution of canopy structure to ΔGPP in the earlier growth stages, we analysed the canopy absorbance for SE, SL and CL days (Fig. 4). In the earlier stages, canopy absorbance decreased with time in the morning and increased in the afternoon because the canopy was not closed to fully cover the ground (LAI = 2.6 for ambient and 2.9 for elevated [CO 2 ]), resulting in more light penetrating through the canopy and reaching the soil surface when solar elevation angle increased (Fig. 4B). In later growth stages, canopy absorbance was relatively constant over the course of a day (Fig. 4E, H). The difference in daily averaged canopy absorbance between the two [CO 2 ] conditions was larger in SE (6.2 %) than SL (1.2 %) and CL (1.0 %) (Fig. 4C, F, I).
The PPFD distribution in SE was more uniform than in SL (Fig. 5A, D). The canopy in SL experienced more scattered PPFD or sunflecks in the botttom layers than that in CL (Fig.  5D, G). The light extinction coefficients for canopies developed under elevated [CO 2 ] and ambient [CO 2 ] were almost the same in either SL, CL or SE days (Fig. 5B, E, H). With increasing cLAI, canopy absorbance increased but did not saturate in SE, increased and approached saturation at about cLAI = 6 in SL, and increased and approached saturation at about cLAI = 3 in CL. The marginal canopy absorbance, defined as d(Abs)/d(cLAI), was much higher (0.24) in SE, comparing to those at the later stages (0.04 for SL and 0.02 for CL) (Fig. 5C, F, I).
Synergistic effect of PPFD on the contribution of [CO 2 ] to ΔGPP For mature canopies (DOY > 207), daily ΔGPP CO 2 was positively correlated (R 2 = 0.727) with the daily averaged ambient solar PPFD (Fig. 6A). In Fig. 6A, diurnal average PPFD is the diurnal averaged ambient solar PPFD, while points represent PPFD of different days of the later developmental stages. For these days, we also calculated ΔGPP contributed by elevated [CO 2 ], and found that on days with high ambient PPFD (i.e. sunny days), 'ΔGPP contributed by elevated CO 2 ' was higher; on days with low ambient PPFD (i.e. cloudy days), the 'ΔGPP contributed by elevated CO 2 ' was lower. Although daily GPP is a function of daily total intercepted solar radiation, here we used the average ambient PPFD as the x-axis to ease comparison with Fig. 6B, which shows the simulated light response curves of leaf photosynthetic CO 2 uptake rate (AQ curve) for both ambient and elevated CO 2 conditions on 219 DOY (Fig. 6B) (data of AQ curves for 210-252 DOY are given in Supplementary Data Table S11). When PPFD was lower than about 800 µmol m −2 s −1 and photosynthesis was limted by RuBP regeneration, the increase in leaf photosynthesis rate (ΔP) was about 11 % under elevated CO 2 compared to ambient CO 2 condition (Fig. 6B). When PPFD was greater than 800 µmol m −2 s −1  and photosynthesis was limited by Rubisco, ΔP was as high as 24 % (Fig. 6B).

DISCUSSION
This study presents a new integrative framework that coupled an explicit 3-D soybean architecture model with a ray tracing algorithm (Song et al., 2013) and a leaf photosynthesis model (Farquhar et al., 1980;Ball et al., 1987;Monteith & Unsworth, 2007) to compute whole canopy photosynthetic response under different environments. In addition, the integrated model also incorporated the responses of photosynthetic parameters to temperature (Bernacchi et al., 2001. The integrated model was employed in this study to dissect the contribution of different factors to the changes in GPP of soybean grown under elevated [CO 2 ]. Model simulations over the entire growing season demonstrated that CO 2 fertilization and structural acclimation significantly increased whole canopy ΔGPP, while photosynthetic acclimation and leaf temperature changes played a minor effect (Fig. 3). Furthermore, we show that the impacts of these different factors on the observed ΔGPP varied with different canopy architecture parameterizations and weather conditions (Fig. 3). Finally, we demonstrate a synergetic effect of CO 2 and light on ΔGPP. Specifically, an increase in LAI under elevated [CO 2 ] during the early developmental stage dramatically increased light absorption and hence canopy photosynthesis (Figs 4 and 5); at later developmental stages, on days with high ambient light (i.e. on bright sunny days), the contribution of elevated CO 2 to GPP was larger than that under lower ambient light (i.e. on cloudy days), because the proportion of leaves in a canopy undertaking Rubisco-limited photosynthesis increased when ambient PPFD increased (Fig. 6).

A new method to dissect the contribution of different factors and their interactions to ΔGPP
This study used the canopy architectural parameters for each devleopmental stage to construct 3-D canopy models along the growing season (Fig. 1). This enabled quantification of GPP for plants at each developmental stage ( Fig. 2A, B). Previous studies have used the sunlit-shaded model to simulate light environments inside the canopy and canopy photosynthesis (DePury & Farquhar, 1997;Wang & Leuning, 1998). In this study, we applied a ray tracing algorithm coupled with 3-D canopy architecture model to simulate the fine details or heterogeneities of light environments in the soybean canopies (Fig. 5A, D, G).
Another major feature of the model was the 3-D canopy photosythesis model along the growing season. In some previous studies, 3-D canopy photosynthesis models have been built for a   particular developmental stage (Zheng et al., 2008;Song et al., 2013;Pound et al., 2014). Here we developed a routine to extrapolate the architectural parameters along the whole growth cycle using architectural parameters measured in a representative developmental stage. Our method uses probability tables to randomize node numbers on the main stem and branches (Supplementary Data Tables S5 and S6). The senescence of soybean leaves begins from the bottom of a canopy, with senescent leaves dropping to the ground causing a decrease in LAI at later developmental stages (Setiyono et al., 2008). Such leaf senescence has been modelled by decreasing the LAI (Yin, 2000). In this study, we modelled the number of senescent leaves with a probability table and those senescenced leaves was removed directly from the 3-D canopy model (Table S8; Fig. 1). The model-calculated GPP was linearly correlated with the measured above-ground biomass  along the growing season, and the obtained coefficient of determination, R 2 , was higher than 0.99 (Fig. 2D), which suggests that the model can be uesd to study factors influencing variations of GPP and hence biomass production. Assuming a fixed root : shoot ratio (Ainsworth et al., 2002) and a fixed fraction of dark respiration (Amthor et al., 2001), the model predicted a 21.4 % increase in above-ground biomass when increasing [CO 2 ] from 370 to 550 ppm ( Fig. 2; Supplementary Data Table  S10). Our model prediction is largely consistent with results from earlier studies. Previous studies showed that when [CO 2 ] is increased to about 700 ppm, total biomass increases by about 37 % (Ainsworth et al., 2002); when [CO 2 ] is increased to about 550 ppm, above-ground biomass increases by about 17-18 % . The difference between the measured and predicted increase in biomass can be attributed to a number of factors: difference between the predicted vs. measured V cmax and J max for leaves at different layers of a canopy, and potential heterogeneity of microclimatic factors other than light (e.g. CO 2 , humidity) inside canopies. These areas need to be improved in future canopy photosynthesis modelling studies.

Dissection of the contribution of different factors to ΔGPP
This model provides a unique opportunity to dissect the contribution of different factors to changes in GPP (ΔGPP) under elevated [CO 2 ] as compared to that under ambient [CO 2 ]. The CO 2 fertilization effect showed the greatest contribution (76.7 %), followed by canopy architecture (17.2 %), to ΔGPP (Fig. 3A). This dominant role of elevated [CO 2 ] to ΔGPP is consistent with an earlier study which showed that changes in the canopy photosynthetic energy conversion efficiency contributed 80 % and the interception efficiency contributed 20 % to the increase in soybean yield under elevated [CO 2 ] (Dermody et al., 2008). The contribution of canopy structure, as shown in the present study, can be further divided into the contribution of LAI and contribution from architecture. LAI directly determines the total photosynthesis leaf area and greatly influences canopy photosynthesis (Song et al., 2013). In this study, the leaf size of soybean was 1. Changes in leaf temperature showed a minor influence (−0.1 %) on ΔGPP. Temperature can influence a number of factors related to ΔGPP. First, it influenced leaf photosynthetic rates. In SoyFACE, the soybean growth temperature was around the optimal temperature for soybean leaf photosynthesis ( Supplementary Data Fig. S1), which was reported to be at about 28 °C based on the temperature response to photosynthesis (Bernacchi et al., 2001Medlyn et al., 2002). A further increase in leaf temperature in SoyFACE decreased canopy photosynthesis. Second, temperature can alter leaf respiration. In a previous study, the negative impact of temperature on rice yield has been attributed to increased respiration at night (Peng et al., 2004). Increased leaf temperature can increase leaf respiration rate, which is usually modelled via the Q10 parameter. Note that total canopy respiration is usually positively correlated with LAI and canopy total nitrogen content. LAI increased by about 19 % under elevated [CO 2 ] (cLAI in Fig. 5E, H) which caused and increase in total canopy respiration; however, the leaf nitrogen content on a leaf area basis decreases by an average of 3.9 % in crops under elevated [CO 2 ] as result of Rubisco acclimation (Leakey et al., 2009).
Under elevated [CO 2 ] conditions, while Rubisco content decreases by 19 %, V cmax at 25 °C decreases only by 15 % due to increased Rubisco activation state (Ainsworth & Long, 2005). The photosynthetic acclimation to growth under elevated [CO 2 ] observed in later developmental stages under the elevated [CO 2 ] condition increases the nitrogen use efficiency of the photosynthesis system  and contributed to ΔGPP by about −6.7 % (Fig. 3).

Synergetic effect of CO 2 and light on ΔGPP
This study reveals a synergetic effect of CO 2 and light on ΔGPP. On the one hand, elevated CO 2 promoted photosynthesis and growth of the plant canopy, which resulted in a higher LAI. The higher LAI led to more light absorption, which increased canopy photosynthesis. At early developmental stages, when LAI was relatively low, canopy absorbance under elevated [CO 2 ] can be 6.2 % higher than that under ambient [CO 2 ] (Fig.  4C). This positive effect of LAI on light absorbance decreased with increasing LAI; at later developmental stages, the difference in canopy absorbance between canopies grown under elevated vs. ambient [CO 2 ] was only about 1-1.2 % because the canopy intercepted about 90 % of total PPFD when LAI > 6 ( Fig. 5F, I). The increase in leaf area under elevated CO 2 can increase canopy respiration, which decreased the positive impact of increasing LAI on net canopy photosynthesis. In fact, greater than optimal LAI can even lead to decreased canopy photosynthesis (McCree & Troughton, 1966;Anten et al., 1995;Song et al., 2013).
The synergistic effect between CO 2 and light on ΔGPP is shown by the correlation between daily averaged PPFD and ΔGPP (Fig. 6A), which shows that the impact of elevated [CO 2 ] on ΔGPP was positively correlated with ambient PPFD. Canopy photosynthetic rate is the integral of photosynthetic rates of all leaves in a canopy. Every facet of the every leaf in a canopy is under either Rubisco limitation or RuBPregeneration limitation depending on its photosynthetic parameters and the absorbed PPFD. The proportion of leaves in which photosynthesis is limited by Rubisco is higher on sunny days than on cloudy days. CO 2 can influence photosynthesis through two mechanisms, i.e. suppressing photorespiration and increasing substrate availability to Rubisco . When photosynthesis was limited by RuBP regeneration, leaf photosynthetic rate under a [CO 2 ] of 550 ppm was 11 % higher than that under a [CO 2 ] of 370 ppm, mainly due to suppression of photorespiration (at low light in Fig. 6B); in contrast, when photosynthesis was limited by Rubisco, a potential 24 % increase in leaf photosynthetic rate was predicted as a result of both suppressing photorespiration and increased substrate availablity for Rubisco (at high light in Fig. 6B). Therefore, with an increase in ambient PPFD, the proportion of leaves under Rubisco-limited photosynthesis increases from 12.2 % under low PPFD to 35.6 % under high PPFD (Supplementary Data Fig. S2), which leads to an increased contribution of elevated CO 2 on canopy photosynthesis or ΔGPP.

CONCLUSION
We integrated an existing 3-D canopy photosynthesis model (Song et al., 2013) with data describing the architectural and physiological changes of soybean grown under ambient and elevated [CO 2 ] to build a soybean canopy model throughout the whole growing season. The integrated 3-D soybean models were then used to dissect the contribution of the different acclimation responses, i.e. [CO 2 ], V cmax , J max , canopy architecture, leaf temperature and their interactions, on whole-canopy GPP. The study demostrates the synergetic effect of [CO 2 ] and light on GPP under elevated [CO 2 ].

SUPPLEMENTARY DATA
Supplementary data are available online at https://academic. oup.com/aob and consist of the following. Table S1: Leaf length of soybean in the mature stage under both ambient and elevated CO 2 conditions. Table S2: Leaf widths of soybean in the mature stage under both ambient and elevated CO 2 conditions. Table  S3: Branch angle, petiole angle, mid-leaf angle and leaf angles of leaf, middle and right leaves in a trifoliate leaf in the mature stage used for parameterization of the soybean architecture models under both ambient and elevated CO 2 conditions. Table S4: Internode and petiole lengths for each node in the main stem and those for branches used for parameterization of models under ambient CO 2 conditions. Table S5: Probabilities of node numbers for the main stem and the branches used for parameterization of the canopy architecture model under both ambient and elevated CO 2 conditions. Table S6: Maximal node numbers of the main stem and the branches for Br1 to Br2 of soybean under ambient and elevated CO 2 conditions. Table S7: Maximal node numbers of the main stem and branches for Br3 to Br6 of soybean under ambient and elevated CO 2 conditions. Table S8: Senescence leaf numbers in every 3 d for soybean grown under ambient and elevated CO 2 conditions. Table S9: V cmax and J max for top leaves of soybean grown under ambient and elevated CO 2 conditions during a growing season. Table  S10: Biomass, model-calculated cumulative GPP. Table S11: Leaf photosynthesis at different PPFD under both ambient and elevated [CO 2 ] conditions. Figure S1: Air temperature, relative humidity and photosynthetic photon flux density during the growing season from 168 DOY to 267 DOY used for modelling canopy photosynthesis. Figure S2: The proportion of Rubiscolimited photosynthesis under high light and under low light. Methods: Pseudocode for the randomized process.