# Functional-Structural Plant Modeling Highlights How Diversity in Leaf Dimensions and Tillering Capability Could Promote the Efficiency of Wheat Cultivar Mixtures

Emmanuelle Blanc, Pierre M Barbillon, Christian Fournier, Christophe Lecarpentier, Christophe Pradal, Jérôme Enjalbert

► To cite this version:

Emmanuelle Blanc, Pierre M Barbillon, Christian Fournier, Christophe Lecarpentier, Christophe Pradal, et al.. Functional-Structural Plant Modeling Highlights How Diversity in Leaf Dimensions and Tillering Capability Could Promote the Efficiency of Wheat Cultivar Mixtures. Frontiers in Plant Science, 2021, 12, pp.734056. 10.3389/fpls.2021.734056. hal-03358998

HAL Id: hal-03358998

https://inria.hal.science/hal-03358998

Submitted on 7 Oct 2021

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

# Functional-Structural Plant Modeling Highlights How Diversity in Leaf Dimensions and Tillering Capability Could Promote the Efficiency of Wheat Cultivar Mixtures

**OPEN ACCESS**

**Edited by:**
Marcello Mastronili,
Council for Agricultural and
Economics Research (CREA), Italy

**Reviewed by:**
Simon Pearson,
University of Lincoln, United Kingdom
Alessandro Persiani,
Council for Agricultural and
Economics Research (CREA), Italy

**\*Correspondence:**
Emmanuelle Blanc
emanuelle.blanc@
univeriste-paris-saclay.fr
Jérôme Enjalbert
jerome.enjalbert@inrae.fr

**Specialty section:**
This article was submitted to
Plant Breeding,
a section of the journal
Frontiers in Plant Science

**Received:** 30 June 2021
**Accepted:** 30 August 2021
**Published:** 29 September 2021

**Citation:**
Blanc E, Barbillon P, Fournier C, Lecarpentier C, Pradal C and Enjalbert J (2021) Functional-Structural Plant Modeling Highlights How Diversity in Leaf Dimensions and Tillering Capability Could Promote the Efficiency of Wheat Cultivar Mixtures. *Front. Plant Sci.* 12:734056.
doi: 10.3389/fpls.2021.734056

Emmanuelle Blanc¹*, Pierre Barbillon², Christian Fournier³, Christophe Lecarpentier⁴, Christophe Pradal⁵,⁶ and Jérôme Enjalbert¹*

¹ Université Paris-Saclay, INRAE, CNRS, AgroParisTech, GQE—Le Moulon, Gif-sur-Yvette, France, ² Université Paris-Saclay, AgroParisTech, INRAE, UMR MIA-Paris, 75005, Paris, France, ³ LEPSE, Univ Montpellier, INRAE, Institut Agro, Montpellier, France, ⁴ UMR 1095 GDEC, INRAE, Université Clermont Auvergne, Clermont-Ferrand, France, ⁵ CIRAD, UMR AGAP Institut, Montpellier, France, ⁶ INRIA and LIRMM, Univ Montpellier, CNRS, Montpellier, France

Increasing the cultivated diversity has been identified as a major leverage for the agroecological transition as it can help improve the resilience of low input cropping systems. For wheat, which is the most cultivated crop worldwide in terms of harvested area, the use of cultivar mixtures is spreading in several countries, but studies have seldom focused on establishing mixing rules based on plant architecture. Yet, the aerial architecture of plants and the overall canopy structure are critical for field performance as they greatly influence light interception, plant interactions and yield. The very high number of trait combinations in wheat mixtures makes it difficult to conduct experimentations on this issue, which is why a modeling approach appears to be an appropriate solution. In this study, we used WALTER, a functional structural plant model (FSPM), to simulate wheat cultivar mixtures and try to better understand how differences between cultivars in key traits of the aerial architecture influence mixture performance. We simulated balanced binary mixtures of cultivars differing for different critical plant traits: final height, leaf dimensions, leaf insertion angle and tillering capability. Our study highlights the impact of the leaf dimensions and the tillering capability on the performance of the simulated mixtures, which suggests that traits impacting the plants' leaf area index (LAI) have more influence on the performance of the stand than traits impacting the arrangement of the leaves. Our results show that the performance of mixtures is very variable depending on the values of the explored architectural traits. In particular, the best performances were achieved by mixing cultivars with different leaf dimensions and different tillering capability, which is in agreement with numerous studies linking the diversity of functional traits in plant communities to their productivity. However, some of the worst performances were also achieved by mixing varieties differing in their aerial architecture, which suggests that
diversity is not a sufficient criterion to design efficient mixtures. Overall, these results highlight the importance of simulation-based explorations for establishing assembly rules to design efficient mixtures.

**Keywords:** aerial architecture, competition for light, variety mixture, sensitivity analysis, functional traits, tillering plasticity

## INTRODUCTION

Increasing the cultivated diversity has been identified as a major leverage for the agroecological transition as it can help improve the resilience of low input cropping systems (Malézieux, 2012; Isbell, 2015). At the scale of the plot, diversity can be increased by mixing species or by increasing the genetic diversity within a species (varietal mixture for example). For wheat, which is the most cultivated crop worldwide in terms of harvested area (Food Agriculture Organization of the United Nations, 2018), the use of cultivar mixtures has been reported to present advantages for yield and quality as well as for diseases resistance, insect pests control, weed suppression, lodging limitation, exploitation of water and soil nutrients and yield stability (Borg et al., 2018). Moreover, the use of wheat cultivar mixtures is spreading in several countries (Faraji, 2011). In France, for example, wheat blends represented <1% of the total wheat surface in 2010 but they represented more than 11% of the total wheat surface in 2019 (FranceAgriMer, 2019). To support the spread of cultivar mixtures, assembly rules have been developed, but they're almost exclusively focused on diseases resistance (Borg et al., 2018). Studies have seldom focused on establishing mixing rules based on the plants architecture and the few existing recommendations usually advocate the association of homogeneous cultivars, but the demonstration of these guidelines seems to be lacking (Borg et al., 2018).

Yet, plant architecture and the overall canopy structure are critical for field performance, as they strongly influence light interception (Niinemets, 2010). Thus, many traits of the aerial architecture, such as height, leaf dimensions, leaf inclination and branching (i.e., tillering in cereals), have been identified for their impact on light interception. Indeed, leaf dimensions have a direct influence on the leaf area index (LAI) of the canopy, which is the most basic property affecting its light interception (Pugnaire and Valladares, 2007). Branching is also an important feature as it greatly impacts the total leaf area of the plant by changing its number of leaves (Valladares and Niinemets, 2008). On the other hand, plant height has a direct impact on the leaf area distribution in the canopy and on the overlap between neighboring leaves, which is critical for competition for light, especially in heterogeneous canopies (Givnish, 1982). As for leaf inclination, it determines the penetration of light into the canopy: the more upright the leaves are, the deeper light can penetrate into the lower layers of the canopy (Long et al., 2006). The impact of all of these traits has been considered for the definition of wheat ideotypes (i.e., ideal hypothetical plants for a given context) (Donald, 1968) and in breeding programs (Reynolds et al., 2012), but mostly in the context of single plants or pure stands. However, the impact of these architectural traits would

also be of interest for the establishment of assembly rules and the definition of mixture ideotypes (Litrico and Violle, 2015). To this end, it is necessary to better understand the interactions taking place in mixtures and the role played by architectural features in these interactions.

The very high number of trait combinations in wheat mixtures makes it difficult to conduct experimentations on this issue, which is why a modeling approach appears to be an appropriate solution. In particular, functional structural plant models (FSPM, Godin and Sinoquet, 2005; Vos et al., 2010) have the advantage of being individual-based, allowing for the explicit consideration of inter-individual variations in the stand (Zhang and De Angelis, 2020). Most importantly, FSPM explicitly represent the 3D architecture of plants and its interactions with the ecophysiological processes controlling plant development. These models are thus particularly adapted to study the impact of architectural traits on the performance of heterogeneous canopies (Evers et al., 2019; Gaudio et al., 2019). If past studies using FSPM were mainly focusing on sole crops (Sarlikioti et al., 2011; Chen et al., 2014; Da Silva et al., 2014; Streit et al., 2016; Perez et al., 2018), more and more work is now addressing mixtures of different species or cultivars (Barillot et al., 2014; Munz et al., 2018; Louarn et al., 2020). For wheat cultivar mixtures, Barillot et al. (2019) have studied the impact of differences in leaf inclination during the post-anthesis period; and Vidal et al. (2018) have explored how diversity in plant height could help control rain-borne diseases. However, no FSPM study has yet focused on the impact of more than one trait on the performance of wheat cultivar mixtures.

In this study, we used an FSPM to simulate balanced binary wheat cultivar mixtures differing in four traits of their aerial architecture. Using a sensitivity analysis approach, we were able to identify the major traits impacting the performance of the simulated stands. This study improves our understanding of how differences between cultivars in key traits of the aerial architecture influence mixture performance.

## MATERIALS AND METHODS

All the statistical analyses were done using the R software (R Core Team, 2018).

### WALTER: A 3D Wheat Model

WALTER (Lecarpentier et al., 2019) is an FSPM that simulates the development of the aerial architecture of winter wheat (*Triticum aestivum* L.) from sowing to maturity with a daily time step. In the model, the vegetative development of the plants (initiation, emergence, elongation and senescence of organs) follows a thermal time schedule and is based on a formalism
derived from ADEL-Wheat (Fournier et al., 2003). The geometry, size, developmental kinetics and the senescence of leaves and internodes are simulated thanks to deterministic equations. However, WALTER is not a completely deterministic model, as the position of the plants, the orientation of the organs, the duration before plant emergence and the emergence of tillers are partly random.

Most importantly, WALTER simulates the competition for light between plants within the field and the resulting plasticity of tillering (i.e., the branching ability of grasses), thanks to a radiative model (CARIBU: Chelle et al., 1998). The regulation of tillering in the model is based on three simple rules. (i) Empirically fixed probabilities control the emergence of new tillers. Then, (ii) the cessation of tillering is controlled by an early neighbor perception: plants stop emitting new tillers when the surrounding Green Area Index<fnref n="1" /> (GAI) reaches a critical value (GAIc). Finally, (iii) some of the tillers that were emitted regress: a tiller regresses if the amount of photosynthetically active radiation (PAR) it intercepts per unit area falls below a threshold (PARt).

Based on this formalism, WALTER produces meaningful outputs, such as the tillering dynamics of each plant (i.e., the number of axes on a plant for each day of the simulation) and the dynamics of light interception for each organ in the stand.

Since its publication in 2019 (Lecarpentier et al.), WALTER has undergone some changes aiming at improving the performances of the simulations, enhancing its realism and improving its ability to simulate mixtures of varieties, reusing several models distributed in the OpenAlea platform (Pradal et al., 2008, 2015). Thus, in the new version of WALTER, the representation of the leaves is based on the model by Fournier and Pradal (2012). The shape of the blades is based on the equations described in Dornbusch et al. (2009) and it is now possible to represent curved leaves based on the formalism described by Perez et al. (2016). The formalism for the dimensions of the leaves remains unchanged. In addition, the discretization of the sky was optimized by the use of a TURTLE sky (den Dulk, 1989) as implemented in the Python package Alinea.astk [package python]. (version 2.1.0, 2019). Furthermore, the option to consider an infinite periodic canopy, which is integrated in the radiative model (Chelle et al., 1998), was implemented in WALTER. Thus, it is no longer necessary to simulate additional plants to discard border effects. The computational cost of the simulations was further reduced by removing all dispensable (i.e., non-visible) organs from the 3D representation of plants. Moreover, a project manager, that serializes parameters, inputs, and outputs for each simulation, has been designed to distribute the computation on high-performance computing infrastructure. Modifications were also implemented to WALTER to allow the simulation of complex, heterogeneous canopies, such as varietal mixtures and populations, in which plants can differ from each other in many parameters. Finally, a fitness module has been implemented to allow the computation of the total fitness of the plot (Ftot) for each simulation. This new output is used as a proxy of the number of kernels produced by the plot and is

computed by WALTER as the sum of the fitness of all the axes in the plot. The computation of the fitness of an axis is based on the following equation:

$$F_{axis} = \frac{PAR_{i45}}{T_{m45}}$$

where Faxis is the fitness of the axis (expressed in molPAR·°Cd⁻<fnref n="1" />); PARi45 is the mean daily amount of PAR intercepted by the axis during the 45 days preceding flowering (in molPAR·d⁻<fnref n="1" />); and Tm45 is the mean temperature during the same period (in °C).

This equation is an adaptation of the photothermal quotient (Nix, 1976), which has been demonstrated to be linearly correlated to the number of grains produced in wheat (Midmore et al., 1984; Fischer, 1985; Abbate et al., 1995; Demotes-Mainard and Jeuffroy, 2001).

WALTER is available as an open source Python package on the OpenAlea platform (https://github.com/openalea/walter).

## Sensitivity Analysis

The impact of different architectural traits on the performance of wheat mixtures was evaluated through the simulation of balanced binary mixtures (mixtures of two varieties with 50% of each cultivar). To evaluate and hierarchize the influence of the different traits on the performance of the simulated mixtures, a sensitivity analysis was conducted with the variance-based Sobol method (Sobol, 1993). Because of the relatively high computational cost to run WALTER, a metamodeling approach was chosen to enable the computation of Sobol indices. The metamodel consists in a fast approximation of WALTER constructed from a limited number of model runs. We used a Gaussian process metamodel (Sacks et al., 1989; Currin et al., 1991), also known as Kriging metamodel, as in Marrel et al. (2009).

Three key outputs of the model were selected for the analyses to reflect the performance of each simulated plot. First, the mean number of ears per plant at maturity (Nears), which is a key component of the plot's yield, was selected. Second, the proportion of incident photosynthetically active radiation (PAR) intercepted by the stand (Lperc), was also taken into account. For each stand, the daily Lperc was averaged over the whole period of use of the radiative model (from around 100 days after sowing to the maturity of the stand). Third, the total fitness of the plot (Ftot) was used for the analyses.

To allow the comparison of the performance of a mixture with the performance of its pure components, the overyielding was also computed for each of these three outputs, using the following equation:

$$OY_{out} = \frac{out_{mix}}{\frac{out_{pure1} + out_{pure2}}{2}}$$

where OYout is the overyielding for the considered output (Nears, Lperc or Ftot); outmix is the value of the considered output for the mixture; and outpure1 and outpure2 are the values of the considered output for each of the two components of the mixture, when cultivated in pure stands.

<footnote n="1">¹The GAI is the ratio of photosynthetic surface of a plant to ground area.</footnote>
<table><caption>**TABLE 1** | List of parameters for the sensitivity analysis, definitions, traits impacted by the parameter, ranges of variation explored for the mean parameter value in the stand (Ref range) and for the difference of trait value between the two varieties of the mixture (Diff range), corresponding minimum and maximum values at the plant scale and units.</caption><thead><tr><th>Parameter</th><th>Description</th><th>Trait</th><th>Ref range</th><th>Diff range</th><th>Min value</th><th>Max value</th><th>Unit</th></tr></thead><tbody><tr><td>H<sub>MS</sub></td><td>Final height of the main stem</td><td>Plant height</td><td>65; 135</td><td>-50; 50</td><td>40</td><td>160</td><td>cm</td></tr><tr><td>Φ<sub>B</sub></td><td>Insertion angle between the blade and the stem<sup>a</sup></td><td>Leaf inclination</td><td>30; 50</td><td>-40; 40</td><td>10</td><td>70</td><td>degree</td></tr><tr><td></td><td>Final length of the longest blade of the main stem</td><td>Leaf dimensions</td><td>16; 27</td><td>-16; +16</td><td>8</td><td>35</td><td>cm</td></tr><tr><td>GAI<sub>c</sub></td><td>Green Area Index threshold above which the emission of tillers stops</td><td>Tillering</td><td>0.4; 0.7</td><td>-0.6; 0.6</td><td>0.1</td><td>1</td><td>-</td></tr><tr><td>PAR<sub>t</sub></td><td>PAR threshold below which a tiller does not survive</td><td>Tillering</td><td>0.3; 0.5</td><td>-0.2; 0.2</td><td>0.2</td><td>0.6</td><td>mol·cm<sup>-2</sup>·°Cd<sup>-1</sup></td></tr></tbody></table>

Five parameters associated with four key traits were considered (Table 1). For each of the five parameters, the impact of the differential between the two varieties in the mixture was considered (diff), as well as the impact of the mean value in the plot (ref). As a result, there were two input factors considered for each parameter, for a total of 10 input factors for the sensitivity analysis. The ranges of variation for all the input factors were selected to generate biologically realistic varieties with contrasted architectures. For example, the selected ranges allowed for the simulation of varieties with a final height between 40 and 160 cm (Table 1).

To build the Kriging metamodel, we used a set of WALTER runs following a 10 dimensions maximin latin hypercube sampling (LHS; McKay et al., 1979; Johnson et al., 1990) of 5,000 simulations generated with the DiceDesign R package (Dupuy et al., 2015) with the ranges detailed in Table 1. For each output (N<sub>ears</sub>, L<sub>perc</sub> and F<sub>tot</sub>), a Kriging approximation was generated using the DiceKriging R package (Roustant et al., 2012). The sensitivity R package (Iooss et al., 2017) was then used to estimate the Sobol indices. To be able to compute the overyielding for every mixture, another set of WALTER runs was generated from a five dimensions maximin LHS of 1,000 simulations. This second LHS only includes simulations of pure stands (no trait differential in the stand) and the ranges explored for each input parameter corresponded to the minimum and maximum values detailed in Table 1. Based on this second set of simulations, a Kriging model was constructed independently for each of the three outputs (N<sub>ears</sub>, L<sub>perc</sub> and F<sub>tot</sub>). Thus, the first design allowed the approximation of the behavior of mixtures while the second design allowed the approximation of the behavior of the corresponding pure components. Using the Kriging metamodels from both designs, it was possible to compute sensitivity indices for the overyielding of each output.

The simulated plots were composed of 110 plants sown at a density of 200 plants/m². An intermediate sowing density was chosen to allow for competition between plants while allowing plants to potentially emit and maintain several tillers. The number of plants to simulate was the object of preliminary work taking into account the impact of the stochasticity of WALTER on the model outputs (data not shown). The plants

were arranged in a grid pattern with plants equidistant from each other and the different cultivars were distributed randomly in the plot (Figure 1). Apart from the input factors of the sensitivity analysis described above and in Table 1, all genotypic parameters were identical for both varieties in the mixture and had values given in Lecarpentier et al. (2019) for the cultivar Maxwell. The climatic sequence used for the simulations was obtained by averaging the climatic sequences of five French locations over 10 years (from 2007 to 2017). The five locations, listed in Supplementary Material, were selected to represent contrasted climatic conditions. For the computation of light interception, only diffuse radiations were considered, according to the standard overcast sky (Moon and Spencer, 1942). Diffuse radiation was approximated using a set of 16 light sources according to a TURTLE sky (den Dulk, 1989).

## Optimization

From the Kriging metamodels, it is possible to seek for combinations of architectural traits maximizing the performance of the stand. However, the uncertainty of the metamodels can affect the results. Therefore, we decided to enrich the initial LHS design with new runs of WALTER, based on an Efficient Global Optimization (EGO) algorithm (Jones et al., 1998). We decided to focus on the maximization of F<sub>tot</sub>, as this output encompasses information about both N<sub>ears</sub> and L<sub>perc</sub>. The initial LHS design of 5,000 simulations was iteratively enriched with a new run of WALTER until 50 new runs were added. At each iteration, the parameters' values for the new run were selected based on the Expected Improvement (EI) criterion among a set of 10,000 parameter combinations sampled from a 10 dimensions latin hypercube. The EI criterion was computed with the DiceOptim R package (Roustant et al., 2012). The EI criterion is balancing the need to maximize F<sub>tot</sub> with the need to reduce the uncertainty of the Kriging metamodel.

## RESULTS

Outputs from the 5,000 initial simulations are highly variable (Supplementary Material), which indicates that at least some of the parameters varying in the LHS design have a strong impact on

<footnote><sup>a</sup>L<sup>B</sup><sub>max</sub> also has an indirect impact on the final length of the other blades and on their width. Thus, at the scale of the cultivar, the higher the value of L<sup>B</sup><sub>max</sub>, the greater the blade area. See **Supplementary Material** for details.</footnote>
<figure><img src="image_3.png" /><figcaption>**FIGURE 1** | 3D representation of a plot simulated with WALTER, at three different stages of development (from left to right: 100, 200 and 300 days after sowing). The simulated plot is a binary mixture composed of 110 plants sown at a density of 200 plants/m². The simulation considers the plot as the pattern of an infinite periodic canopy, removing border effects. Non-senecent organs are represented in green, except for ears that are represented in brown. Regressing tillers are represented in blue and dead organs are represented in yellow.</figcaption></figure>

the performance of the simulated stands. However, some outputs are much more variable than others. Thus N<sub>ears</sub> is the most variable output with values ranging from 1.0 to 4.3 ears per plant, even though 98% of the simulations resulted in N<sub>ears</sub> values below 3.2 ears per plant. On the other hand, L<sub>perc</sub> is much less variable, with values ranging from 0.76 to 0.96. As for F<sub>tot</sub>, the initial simulations gave results ranging from 0.73 to 1.5 mol·°Cd<sup>-1</sup>. For the OY indicators, simulations gave values between 0.40 and 1.6 for N<sub>ears</sub>, and values between 0.67 and 1.9 for F<sub>tot</sub>. Finally, the OY for L<sub>perc</sub> was less variable, with values ranging from 0.97 to 1.4.

## Sensitivity Analysis

The first-order sensitivity indices (MSI) provide information on the mean influence of the different inputs on each output and the total-order indices (TSI) aggregate mean influence and influence through interactions. The most important architectural traits regarding the performance of the stands are dependent on the considered output (Figure 2).

The mean number of ears per plant (N<sub>ears</sub>) is mainly influenced by first-order effects. For this output the dimensions of the leaves is the most influential trait in the simulated plots as it explains more than 68% of the total variance of N<sub>ears</sub> (taking into account the effect of the interactions). The mean value of L<sub>max</sub><sup>B</sup> in the stand has more impact on N<sub>ears</sub> than the difference of L<sub>max</sub><sup>B</sup> between varieties (TSI = 0.41 and 0.28, respectively). The tillering capability of the plants in the stand also has an important impact on N<sub>ears</sub>: GAI<sub>c</sub>-diff explains around 13% of the total variance of the output and PAR<sub>t</sub>-ref accounts for about 12%.

For L<sub>perc</sub>, as for N<sub>ears</sub>, the input factor with the most influence is L<sub>max</sub><sup>B</sup>-ref. Indeed, the variations of L<sub>max</sub><sup>B</sup>-ref explain almost 63% of the total variance of L<sub>perc</sub>. GAI<sub>c</sub> is also an important parameter for light interception with GAI<sub>c</sub>-ref accounting for about 18% of the total variance of L<sub>perc</sub> and GAI<sub>c</sub>-diff explaining more than 10% of its total variance. However, the impact of GAI<sub>c</sub> on L<sub>perc</sub> is mainly explained by interaction effects. Finally, the mean leaf insertion angle in the canopy has an important impact on L<sub>perc</sub> with a total-order sensitivity index of 0.15 for φ<sub>B</sub>-ref.

All the other input factors have very little influence on the light interception efficiency.

The input factor with the most impact on F<sub>tot</sub> is L<sub>max</sub><sup>B</sup>-diff (TSI = 0.62). The mean value of PAR<sub>t</sub> in the mixture (PAR<sub>t</sub>-ref) and the height differential between varieties (H<sub>MS</sub>-diff) also have a lot of influence on F<sub>tot</sub>, as they explain 12 and 10% of the total variance of this output respectively.

For the OY outputs, the two most influential input factors are always L<sub>max</sub><sup>B</sup>-diff and GAI<sub>c</sub>-diff, regardless of the output considered. Together, these two inputs always explain more than 77% of each OY output. Interestingly, interactions between input factors have an important impact on these outputs, particularly for the OY of N<sub>ears</sub> and L<sub>perc</sub>. Nevertheless, differences exist between the different OY outputs. For example, for the OY of N<sub>ears</sub>, GAI<sub>c</sub>-diff has more impact than L<sub>max</sub><sup>B</sup>-diff, whereas it is L<sub>max</sub><sup>B</sup>-diff that has the most influence on the OY of L<sub>perc</sub> and F<sub>tot</sub>.

To a lesser extent, other input factors also have a substantial influence on the OY outputs. PAR<sub>t</sub>-diff, GAI<sub>c</sub>-ref and H<sub>MS</sub>-diff each explain more than 10% of the total variance of the OY of N<sub>ears</sub>. The TSI of PAR<sub>t</sub>-diff for the OY of N<sub>ears</sub> even reaches 0.21, but this effect is mainly explained by interactions as its MSI is only around 0.067. For the OY of L<sub>perc</sub>, GAI<sub>c</sub>-ref, φ<sub>B</sub>-diff and L<sub>max</sub><sup>B</sup>-ref are important factors, each with a TSI >0.1. Finally, PAR<sub>t</sub>-diff explains almost 16% of the total variance of the OY of F<sub>tot</sub>.

## Optimization

The sensitivity indices (Figure 2) only provide information about the relative importance of the different input parameters for the determination of the performance of the stand. Heatmaps representing the outputs of the Kriging metamodels as a function of their most influential parameters (Figures 3–5) give more information on how each parameter affects each output. For these figures however, it is important to keep in mind that the less influential parameters were set to the mean value of their range of variation, thus discarding most of the interaction effect that they might have with the other inputs. Another critical point when
<figure><img src="image_4.png" /><figcaption>**FIGURE 2** | Mean estimates for the total-order (blue triangles) and first-order (red circles) Sobol indices computed for Nears, Lperc, Ftot and the OY of these three outputs, for the 10 input factors (x-axis). The standard-deviation of each index is represented with a vertical line.</figcaption></figure>

reading these heatmaps is the fact that, when the input variable is a parameter difference, it generates a plot that is symmetric on the axis diff = 0. And if the heatmap is a function of two parameters-diff, there is a central symmetry of the plot.

**Figure 3** shows that the mean value of $L_{max}^B$ in the stand ($L_{max}^B$-ref) is negatively correlated to Nears, meaning that canopies with small blades tend to produce more spikes than stands with larger leaves. This heatmap also reveals that a differential of $L_{max}^B$ between the cultivars of a mixture is unfavorable for the ear yield in the simulated plots. Unlike Nears, Lperc is positively correlated to $L_{max}^B$-ref: larger leaves allow for a higher light interception efficiency. Lperc is also positively correlated to GAIc-ref, but mainly for values of GAIc-ref lower than 0.5: above this value, the impact of GAIc-ref on Lperc is very low. This means that Lperc is higher for stands in which plants have a low sensitivity to the early competition (high GAIc values) and can thus emit many tillers.

As for Ftot (Figure 4), unlike the results observed for Nears (Figure 3), a differential of $L_{max}^B$ between the cultivars of a mixture is favorable for the performance of the stand. An absolute differential of around 10 cm between the two varieties in the mixture leads to the highest values of Ftot. Moreover, the highest values of Ftot can only be reached with low values of PARt-ref.

Indeed, there is a clear negative correlation between PARt-ref and Ftot, which means that stands tend to be more productive when plants in the canopy do not need much light to prevent the regression of their tillers. Finally, even though the impact of HMS-diff on Ftot is low for the combinations of parameters presented in Figure 4, the results show that homogeneous plant height in the canopies lead to higher values of Ftot than the mixture of cultivars with a height differential. Interestingly, for mixtures of varieties with a height differential, Ftot is higher if the tallest variety is the one with the smallest leaves.

For the OY outputs, $L_{max}^B$-diff and GAIc-diff are always the most influential inputs (Figure 2). However, the effect of these two inputs, when the other inputs are fixed to their mean values, are very different depending on the output considered (Nears, Lperc or Ftot) (Figure 5). Firstly, the combinations of parameters presented in Figure 5 lead almost exclusively to OY values ≥1 for Ftot. This means that, for these combinations of parameters, mixtures are almost always at least as performant as the mean of their pure components. Conversely, OY values for Nears may be as low as 0.6 and are rarely >1. Finally, for Lperc, the OY is mostly around 1, with very little variation in its value. For this output, the highest values of OY are obtained for the extreme values of $L_{max}^B$-diff and GAIc-diff by mixing plants with large leaves and
<figure><img src="image_5.png" /><figcaption>**FIGURE 3** | Mean outputs predicted by the Kriging metamodels as a function of their two most influential input parameters. Predictions were made with all other parameters set to their mean values. Left: $N_{\text{ears}}$ as a function of $L_{\text{max}}^{\text{B-diff}}$ and $L_{\text{max}}^{\text{B-ref}}$. Right: $L_{\text{perc}}$ as a function of $GAI_{\text{c}}$-ref and $L_{\text{max}}^{\text{B-ref}}$.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>**FIGURE 4** | $F_{\text{tot}}$ predicted by the Kriging metamodel as a function of the 3 parameters that have the most impact on the variance of this output ($L_{\text{max}}^{\text{B-diff}}$, $H_{\text{MS}}$-diff and $\text{PAR}_t$-ref). Predictions were made with all other parameters set to their mean values.</figcaption></figure>

high values of $GAI_c$ with plants with small leaves and low $GAI_c$ values. For $N_{\text{ears}}$, the highest values of OY are also obtained for the extreme values of $GAI_c$-diff, but for values of $L_{\text{max}}^{\text{B-diff}}$ close to 0. For mixtures with a differential in both $L_{\text{max}}^{\text{B}}$ and $GAI_c$, the OY of $N_{\text{ears}}$ is more important when the cultivar with the largest leaves is also the one with the highest value of $GAI_c$, as is the case for the OY of $L_{\text{perc}}$. In contrast, the OY of $F_{\text{tot}}$ is the highest when mixing plants with large leaves and low $GAI_c$ values with plants with small leaves and high $GAI_c$ values. Interestingly, despite the observed differences between the combinations of parameters leading to the highest OY for $N_{\text{ears}}$ and $F_{\text{tot}}$, $N_{\text{ears}}$ is positively correlated to $F_{\text{tot}}$ in the overall dataset (data not shown). Finally, it is worth mentioning that the metamodels predict an OY of 1 for mixtures with no trait differential between cultivars (Figure 5), which means that they perform exactly like the corresponding pure stands. This is expected as all the plants in such canopies have identical parameter values. Thus, this result highlights the consistency of the metamodels' predictions.

In the initial design, different combinations of parameters can describe identical mixtures (as seen in the heatmaps). For example, considering only $\phi_B$-diff and $H_{\text{MS}}$-diff, a stand with the minimum values of these two parameters is identical to a stand with the maximum values of the two parameters. Indeed both combinations of parameters describe a mixture of tall planophile plants with short erectophile plants. To avoid this redundancy and simplify the analysis of the results, Figures 6, 7, 9 have been represented only with positive values of $L_{\text{max}}^{\text{B-diff}}$ (all combinations are considered, but the ones with negative values of $L_{\text{max}}^{\text{B-diff}}$ have been reversed).

Focusing on the 50 runs of the initial design with the highest $F_{\text{tot}}$ values (Figures 6, 7), it is possible to identify the trait combinations leading to the best performances. Simulations with high $F_{\text{tot}}$ values usually combine high $L_{\text{max}}^{\text{B-diff}}$ values with low (negative) values of $\text{PAR}_t$-diff. This combination of parameters describes mixtures in which the genotype with the largest leaves is the one with the lowest value of $\text{PAR}_t$, i.e., the one
<figure><img src="image_7.png" /><figcaption>**FIGURE 5** | Mean overyielding (OY) predicted by the Kriging metamodels for N<sub>ears</sub>, L<sub>perc</sub> and F<sub>tot</sub> as a function of GAI<sub>c</sub>-diff and L<sub>max</sub><sup>B</sup>-diff. Predictions were made with all other parameters set to their mean values.</figcaption></figure>

<figure><img src="image_8.png" /><figcaption>**FIGURE 6** | Density of each input factor for the initial LHS (gray) and for the 50 runs with the highest F<sub>tot</sub> values (red). The medians for each group are represented by a dashed line. To avoid redundancy in the combinations of parameters, results are represented only with positive values of L<sub>max</sub><sup>B</sup>-diff.</figcaption></figure>

less susceptible to tiller regression. However, there is a negative correlation between L<sub>max</sub><sup>B</sup>-diff and GAI<sub>c</sub>-diff (Figure 7), which means that it is preferable for the performance of the stand that the genotype with the largest leaves ceases tillering early in response to the neighboring competition. Furthermore, mixtures with a differential of leaf insertion angle are more likely to have high F<sub>tot</sub> values if the genotype with the most erected leaves is the one with the shortest leaves. Finally, the mean values of H<sub>MS</sub> and PAR<sub>t</sub> in the stand should be rather low to allow for high F<sub>tot</sub> values. This means that canopies of short plants with low tiller regression tend to be the most performant ones.

The addition of new runs of WALTER to the initial LHS design has confirmed the first results. In particular, for the new simulations as for the ones in the initial design, the highest values of F<sub>tot</sub> are obtained for low values of PAR<sub>t</sub>-ref (low tiller regression) and for mixtures with a high diversity in leaf dimensions and in susceptibility to tiller regression (Figures 8, 9). When considering the initial design, the five simulations leading to the highest F<sub>tot</sub> values had variable values of H<sub>MS</sub>-ref and H<sub>MS</sub>-diff (Figure 9). However, the five best-performing simulations within the 50 new runs added by EGO all had low H<sub>MS</sub>-ref values, meaning that high performing
<figure><img src="image_9.png" /><figcaption>**FIGURE 7** | Principal components analysis (PCA) of the 10 input parameters for the 50 simulations of the initial LHS with the highest F<sub>tot</sub>. Projection of the 10 input parameters on the first two axes (A) and on the third and fourth axes of the PCA (B). Projection of the 50 points on the first two axes of the PCA (C) and on the third and fourth axes of the PCA (D). The points are colored according to the corresponding value of F<sub>tot</sub>. To avoid redundancy in the combinations of parameters, results are represented only with positive values of L<sub>max</sub><sup>B-diff</sup>.</figcaption></figure>

canopies are composed of short plants, as suggested by the
results presented in **Figure 6**. Moreover, for the new runs of
WALTer, the mixtures offering the highest F<sub>tot</sub> values presented
an important diversity in plant height, and the shorter genotype
in the mixture was always the one with the larger leaves.
This result is also in accordance with the results presented in
**Figure 6**.

DISCUSSION

FSPM as a Tool to Study Heterogeneous
Canopies

This study illustrates the interest of using FSPMs to explore
the impact of architectural traits on the performance of
heterogeneous canopies (Evers et al., 2019; Gaudio et al., 2019;
<figure><img src="image_10.png" /><figcaption>**FIGURE 8** | Input parameter values for the 50 points added to the design via EGO. The points are colored according to the corresponding value of F<sub>tot</sub>. The additional points were selected in the whole range of variation of the 10 input parameters.</figcaption></figure>

Muller and Martre, 2019; Louarn and Song, 2020; Stomph et al., 2020). Contrary to previous modeling work on wheat cultivar mixtures (Vidal et al., 2018; Barillot et al., 2019), this study focuses on multiple traits and on the entire development cycle of the plants. Thus, it provides a broader insight into the interactions taking place in wheat cultivar mixtures. Nevertheless, as for any modeling exercise (Passioura, 1996), our study presents some approximations that must be considered when interpreting the results. Indeed, even though WALTER allows for a realistic simulation of the tillering dynamics in response to competition for light (Lecarpentier et al., 2019; Blanc et al., 2021), some of the underlying hypotheses of the model should be mentioned. First of all, light is the only environmental factor considered in the model: water and nutrients (nitrogen, phosphorus, potassium, sulfur...) are not explicitly simulated by the model, even though these factors are known to impact the development of wheat. In particular, numerous studies have highlighted the important role of nitrogen in the regulation of tillering (Sparkes et al., 2006; Assuero and Tognetti, 2010; Dornbusch et al., 2011; Alzueta et al., 2012) and taking this factor into account would impact the results. Simulations are thus made with the hypothesis that resources other than light are not limiting. Moreover, plasticity was considered only for tillering, as it is the most plastic trait in response to competition for light in wheat (Lecarpentier, 2017). However, considering plasticity of other traits, and especially

in leaf dimensions, could lead to different results, although it would come at a higher computational cost. Furthermore, the results are also affected by the parameterization of the model. In particular, light interception was computed using only diffuse radiations, because the simulation of more realistic sky conditions (i.e., diffuse and direct radiations) would be significantly more time-consuming. Sky conditions could have an impact on our results and additional information could be provided by the comparison of our simulations with simulations using only direct radiations (Barillot et al., 2019). Nonetheless, our preliminary studies (unpublished data) showed that there were no significant differences between the results obtained with these two conditions when simulating a simple binary mixture. Finally, even though mixtures cultivated by farmers usually combine more than two cultivars, our study focuses on the simulation of binary mixtures, as these simple mixtures allow an efficient exploration of the relations between plant architecture and stand performance.

In contrast to these limitations, the strengths inherent in the modeling approach allowed us to conduct a study that would have been impossible to achieve experimentally. Contrary to experimental studies on wheat mixtures (Borg et al., 2018; Montazeaud et al., 2020), the effect of each architectural trait on the development of the stand could be isolated from any other factor. In addition, it was possible to vary the different traits independently, without including the existing
<figure><img src="image_11.png" /><figcaption>**FIGURE 9** | Parameter values for the five simulations with the highest $F_{tot}$ in the initial design (left) and in the 50 simulations added by EGO (right). The x-axis is an indication of the ranking of the simulations based on the $F_{tot}$ values. Parameter values were scaled and centered so they could all be represented on the same axis. To avoid redundancy in the combinations of parameters, results are represented only with positive values of $L_{max\_diff}^B$.</figcaption></figure>

correlations between architectural traits, such as the one between leaf dimensions and leaf inclination (Ma et al., 2020), even though these biological constraints should be considered when interpreting the results. The sensitivity analysis combined with a metamodeling approach allowed us to take advantage of the complexity of an FSPM while limiting the computational time required. Thanks to the adaptive design, the uncertainty associated with the metamodeling approach was accounted for. Indeed, the addition of new runs of WALTer in areas of the parameter space that led to high values of $F_{tot}$ improved the quality of the metamodel approximation in these areas of interest.

## Leaf Dimensions and Tillering Are Major Drivers of Mixture Performance

Our study highlights the impact of the leaf dimensions and the tillering capability on the performance of the simulated mixtures. Indeed, for all the outputs considered, those two traits were the most influential ones among the four architectural traits explored. Leaf dimensions and tillering capability have an important impact on the plant's leaf area index (LAI) whereas leaf inclination and height don't have a direct impact on the LAI but rather on the spatial distribution of the leaf area and on its orientation. Thus, our results suggest that traits impacting the LAI of plants have more influence on the performance of the stand than traits impacting the arrangement of the

leaves. This result is consistent with the well-established effect of LAI on competition for light interception and productivity (Pugnaire and Valladares, 2007). Furthermore, it is in agreement with several studies using FSPMs to simulate trees (Da Silva et al., 2014; Streit et al., 2016; Perez et al., 2018), or maize intercropped with soybean (Munz et al., 2018), in which the traits related to leaf area were identified as the ones with the most influence on light interception. Finally, it is also consistent with the results of the experimental study of Montazeaud et al. (2020), in which the tillering capability was shown to have an important impact on the performance of binary wheat mixtures. The relatively smaller impact of height and leaf inclination on stand performance is congruent with the findings of Munz et al. (2018) on intercropped maize and with the study of Barillot et al. (2019) on leaf inclination in wheat cultivar mixtures. However, even though their impact on mixture performance was less important than that of tillering and leaf dimensions, plant height and leaf inclination still had an important impact on some of the performance indicators. In particular, leaf inclination had an important impact on the light interception efficiency of the simulated canopies ($L_{perc}$) and on its OY. As for the differential of plant height between the cultivars in mixture, it was one of the most important factors influencing $F_{tot}$ and the OY of $N_{ears}$.

Contrary to studies on pure stands or on single plants (Sarlikioti et al., 2011; Da Silva et al., 2014; Streit et al., 2016; Perez et al., 2018), we were able to distinguish the effect of the trait
differential between cultivars from the effect of the mean value of the trait in the stand, for each trait. Thus, for the leaf dimensions, the mean value in the canopy greatly impacts the light interception efficiency ($L_{perc}$) and the number of ears produced ($N_{ears}$), but it is mainly the differential between cultivars that influences the performance of the stand, with a major impact on all outputs except $L_{perc}$. As for the tillering capability, the differential between cultivars and the mean value in the stand both had an important impact on mixture performance. This is in line with the experimental study conducted by Montazeaud et al. (2020), in which mixture performance was impacted by both average trait values and trait differences. In particular, they showed an important effect of both the average tillering capability and its differential on the OY of grain yield.

Interestingly, modeling studies could guide efforts to address the lack of experimental data on intraspecific mixtures. More precisely, our results suggest that it would be interesting to design experiments mixing varieties with differences in leaf dimensions and tillering capabilities. Although it is difficult to know values of these traits in wheat varieties (due to the lack of phenotypic data and the important investment required for the measurements), this type of experiment could, in turn, help with the biological validation of the models.

## Toward the Identification of Assembly Rules

Our results show that the performance of mixtures is very variable depending on the values of the explored architectural traits. In particular, the best performances were achieved by mixing cultivars with different architectures, which is in agreement with numerous studies linking the diversity of functional traits in plant communities to their productivity (Reiss and Drinkwater, 2018). However, some of the worst performances were also achieved by mixing varieties differing in their aerial architecture and the OY indicators were <1 for a significant proportion of the simulated mixtures. This suggests that diversity is not a sufficient criterion to design efficient mixtures and that random-trait assembly will not necessarily lead to better performances than pure stands, as previously stated by Louarn et al. (2020) for grassland mixtures. Overall, these results highlight the importance of establishing assembly rules to design efficient mixtures.

When focusing on the trait configurations leading to the best performances, our results converge toward a single mixture ideotype for the total fitness of the plot ($F_{tot}$). Indeed, even if some parameters can have variable values in well performing stands, due to their low impact on $F_{tot}$, there is no multimodality in the distribution of the parameters with high impact on $F_{tot}$ for the best performing simulations. In particular, regarding the mean values of traits in the stand, the best-performing canopies are composed of short plants with low levels of tiller regression. Interestingly, the fact that stands with short plants perform better than stands with taller plants is in line with the assumption of a trade-off between vegetative development and reproduction. As a matter of fact, significant yield gains have been achieved by breeding for wheat cultivars with reduced height, as a shorter

stature limits lodging and improves the harvest index (Gale et al., 1985). In WALTER, the negative correlation between plant height and stand performance is probably due to the formalism of tiller regression. Indeed, tiller regression depends on the amount of light intercepted per unit area, and taller tillers have a larger surface of internodes and sheaths, which makes them more likely to regress.

Regarding the differential of trait values between the two varieties in the mixture, the best-performing stands are the ones with an important diversity in both leaf dimensions and tiller regression. Our results tend to encourage the association of a very competitive genotype with a genotype with poor competitive ability, as the genotype with the largest leaves must also be the one less sensitive to tiller regression in order for the stand to perform well. This conclusion must however be qualified as the cultivar with the largest leaves should also be the shorter one and the one with the earliest cessation of tillering in response to the neighboring competition. Thus, in the best-performing plots, the most productive genotype at maturity was usually not the one that produced the most tillers. Regardless, the high $F_{tot}$ values offered by this trait combination is explained by the very important productivity gain of the cultivar with the largest leaves, which outweighs the loss of productivity of the other cultivar in the mixture. However, it is important to specify that no plant mortality is considered in WALTER and that all individuals produce at least one ear, which is not always the case in experimental settings. In the case of the best-performing configurations that were simulated with WALTER, the cultivar with the smallest leaves had a very low productivity that might have been even lower in field conditions because of potential plant death or individuals with null productivity. Under such experimental conditions, the gain of the cultivar with the largest leaves might not have compensated for the loss in productivity of the other genotype, thus calling for caution when mixing cultivars with significant differences in their competitive ability. Moreover, the results of Montazeaud et al. (2020) showed a negative effect of diversity in tillering capability on the performance of binary mixtures of spring durum wheat sown in rows, which is in contradiction with our results. However, even if the results of our simulations are not meant to be directly translated into assembly rules, they improve our understanding of the interactions taking place in mixtures.

To further improve our understanding, it would be interesting to consider phenological traits, such as the number of days to heading or to maturity, in addition to the architectural ones. Indeed, diversity in phenology could lead to temporal complementarity in addition to the spatial complementarity offered by differences in traits of the aerial architecture. Several studies have emphasized the impact of phenological traits on the performance of wheat mixtures (Borg et al., 2018; Montazeaud et al., 2020). However, these traits were not yet considered in our study because it would require the integration of a specific module to WALTER to ensure that modifications in phenology and their implication on plant development are reliably simulated. Future work could also include the exploration of the effects of some agricultural practices. We simulated balanced mixtures (with 50% of each
genotype) sown at a density of 200 plants/m². However, it is expected that both sowing density and the proportions of the two cultivars in mixture should have an impact on the results (Grace and Tilman, 1990). Moreover, the sowing density has been identified, as expected, as one of the most influential inputs in WALTer (Lecarpentier et al., 2019; Blanc et al., 2021), supporting the interest in exploring variations of this parameter. FSPMs are particularly adapted to explore these agricultural practices (Evers et al., 2019; Gaudio et al., 2019). WALTer could readily be used to consider the effects of variations in both sowing density and cultivar proportion as well as their interactions with architectural traits, albeit with an increased computational cost compared to the present study.

Finally, the trait combinations identified in this study are the ones maximizing $F_{tot}$, which is a proxy for the number of grains produced by the stand. However, cultivar mixtures can provide a wide variety of ecosystem services besides yield (Barot et al., 2017). Thus, it would be interesting to consider multiple objectives as exemplified by the study of Louarn et al. (2020) in which the stability of grassland mixtures over time was also considered in addition to yield. Considering the stability of the composition of a mixture over time could also be of interest for wheat, as this crop is often resown from year to year. The mixtures that we identified to maximize $F_{tot}$ present significant differences in performance between their two cultivars, which would lead to a drastic evolution of their composition and possibly to a decrease in stand performance over time. Further work with WALTer would be of high interest to identify trait configurations leading to both high performance and mixture stability over time.

## REFERENCES

Abbate, P. E., Andrade, F. H., and Culot, J. P. (1995). The effects of radiation and nitrogen on number of grains in wheat. *J. Agric. Sci.* 124, 351–360. doi: 10.1017/S0021859600073317

Alinea.astk [package python]. (2019). Available online at: https://github.com/ openalea-incubator/astk.git (accessed September 7, 2021)

Alzueta, I., Abeledo, L. G., Mignone, C. M., and Miralles, D. J. (2012). Differences between wheat and barley in leaf and tillering coordination under contrasting nitrogen and sulfur conditions. *Eur. J. Agron.* 41, 92–102. doi: 10.1016/j.eja.2012.04.002

Assuero, S. G., and Tognetti, J. A. (2010). Tillering regulation by endogenous and environmental factors and its agricultural management. *Am. J. Plant Sci. Biotechnol.* 4, 35–48. http://www.globalsciencebooks.info/Online/GSBOnline/ images/2010/AmJPSB_4(SI1)/AmJPSB_4(SI1)35-48o.pdf

Barillot, R., Chambon, C., Fournier, C., Combes, D., Pradal, C., and Andrieu, B. (2019). Investigation of complex canopies with a functional-structural plant model as exemplified by leaf inclination effect on the functioning of pure and mixed stands of wheat during grain filling. *Ann. Bot.*, 123, 727–742. doi: 10.1093/aob/mcy208

Barillot, R., Escobar-Gutiérrez, A. J., Fournier, C., Huynh, P., and Combes, D. (2014). Assessing the effects of architectural variations on light partitioning within virtual wheat-pea mixtures. *Ann. Bot.* 114(4), 725–737. doi: 10.1093/aob/mcu099

Barot, S., Allard, V., Cantarel, A., Enjalbert, J., Gauffreteau, A., Goldringer, I., and Porcher, E. (2017). Designing mixtures of varieties for multifunctional

## DATA AVAILABILITY STATEMENT

The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

## AUTHOR CONTRIBUTIONS

EB, JE, and PB conceived the study, designed the simulation experiments and established the statistical framework. Developments to the WALTer model were made by CF, CL, CP, and EB. EB wrote the R-scripts for the analyses and the figures with input from PB and CL. EB wrote the manuscript. All authors reviewed and approved the final manuscript.

## FUNDING

This work is supported by public grants overseen by the French National Research Agency (ANR) as part of the Investissement d'Avenir program, through the IDI 2017 project funded by the IDEX Paris-Saclay, ANR-11-IDEX-0003-02 and the MoBiDiv project, ANR-20-PCPA-0006.

## ACKNOWLEDGMENTS

The authors thank Franck Gauthier and Bilâl Derouich for their valuable help in the development of WALTer.

## SUPPLEMENTARY MATERIAL

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fpls.2021.734056/full#supplementary-material

agriculture with the help of ecology. A review. *Agron. Sustain. Dev.* 37:13. doi: 10.1007/s13593-017-0418-x

Blanc, E., Enjalbert, J., and Barbillon, P. (2021). Automatic calibration of a functional-structural wheat model using an adaptive design and a metamodelling approach. *bioRxiv*. doi: 10.1101/2021.07.29.454328

Borg, J., Kiær, L. P., Lecarpentier, C., Goldringer, I., Gauffreteau, A., Saint-Jean, S., and Enjalbert, J. (2018). Unfolding the potential of wheat cultivar mixtures: a meta-analysis perspective and identification of knowledge gaps. *Field Crops Res.* 221, 298–313. doi: 10.1016/j.fcr.2017.09.006

Chelle, M., Andrieu, B., and Bouatouch, K. (1998). Nested radiosity for plant canopies. *Vis. Comput.* 14, 109–125. doi: 10.1007/s003710050127

Chen, T. W., Nguyen, T. M. N., Kahlen, K., and Stützel, H. (2014). Quantification of the effects of architectural traits on dry mass production and light interception of tomato canopy under different temperature regimes using a dynamic functional-structural plant model. *J. Exp. Bot.* 65, 6399–6410. doi: 10.1093/jxb/eru356

Currin, C., Mitchell, T., Morris, M., and Ylvisaker, D. (1991). Bayesian prediction of deterministic functions, with applications to the design and analysis of computer experiments. *J. Am. Stat. Assoc.* 86, 953–963. doi: 10.1080/01621459.1991.10475138

Da Silva, D., Han, L., Faivre, R., and Costes, E. (2014). Influence of the variation of geometrical and topological traits on light interception efficiency of apple trees: sensitivity analysis and metamodelling for ideotype definition. *Ann. Bot.* 114, 739–752. doi: 10.1093/aob/mcu034
Demotes-Mainard, S., and Jeuffroy, M. H. (2001). Incorporating radiation and nitrogen nutrition into a model of kernel number in wheat. *Crop Sci.* 41, 415-423. doi: 10.2135/cropsci2001.412415x

den Dulk, J. A. (1989). *The Interpretation of Remote Sensing: a Feasibility Study.* (Doctoral dissertation, Den Dulk).

Donald, C. T. (1968). The breeding of crop ideotypes. *Euphytica*, 17, 385-403. doi: 10.1007/BF00056241

Dornbusch, T., Baccar, R., Watt, J., Hillier, J., Bertheloot, J., Fournier, C., et al. (2011). Plasticity of winter wheat modulated by sowing date, plant population density and nitrogen fertilisation: dimensions and size of leaf blades, sheaths and internodes in relation to their position on a stem. *Field Crops Res.* 121, 116-124. doi: 10.1016/j.fcr.2010.12.004

Dornbusch, T., Watt, J., Baccar, R., Fournier, C., and Andrieu, B. (2009). "Towards a quantitative evaluation of cereal lamina shape using an empirical shape model," in *2009 Third International Symposium on Plant Growth Modeling, Simulation, Visualization and Applications*, 229-236. China, Beijing:IEEE. doi: 10.1109/PMA.2009.29

Dupuy, D., Helbert, C., and Franco, J. (2015). DiceDesign and DiceEval: two R packages for design and analysis of computer experiments. *J. Stat. Softw.* 65, 1-38. doi: 10.18637/jss.v065.i11

Evers, J. B., Van Der Werf, W., Stomph, T. J., Bastiaans, L., and Anten, N. P. (2019). Understanding and optimizing species mixtures using functional-structural plant modelling. *J. Exp. Bot.* 70, 2381-2388. doi: 10.1093/jxb/ery288

Faraji, J. (2011). Wheat cultivar blends: a step forward to sustainable agriculture. *Afr. J. Agric. Res.* 6, 6780-6789. doi: 10.5897/AJARX11.047

Fischer, R. A. (1985). Number of kernels in wheat crops and the influence of solar radiation and temperature. *J. Agric. Sci.* 105, 447-461. doi: 10.1017/S0021859600056495

Food and Agriculture Organization of the United Nations (2018). FAOSTAT statistical database. Available online at: http://www.fao.org/faostat/en/#data/QC (accessed August 25, 2020).

Fournier, C., Andrieu, B., Ljutovac, S., and Saint-Jean, S. (2003). "ADEL-Wheat: a 3D architectural model of wheat development," in *Plant Growth Modeling and Applications*, 54-63. Berlin: Springer Verlag.

Fournier, C., and Pradal, C. (2012). "A plastic, dynamic and reducible 3D geometric model for simulating gramineous leaves," in *2012 IEEE 4th International Symposium on Plant Growth Modeling, Simulation, Visualization and Applications*. China, Shanghai:IEEE, 125-132. doi: 10.1109/PMA.2012.6524823

FranceAgriMer (2019). *Enquête répartition variétale des céréales*. https://visionet.franceagrimer.fr/Pages/accueil.aspx

Gale, M. D., Youssefian, S., and Russell, G. E. (1985). Dwarfing genes in wheat. *Progress in plant breeding*, 1, 1-35. doi: 10.1016/B978-0-407-00780-2.50005-9

Gaudio, N., Escobar-Gutiérrez, A. J., Casadebaig, P., Evers, J. B., Gerard, F., Louarn, G., et al. (2019). Current knowledge and future research opportunities for modeling annual crop mixtures. A review. *Agron. Sustain. Dev.* 39:20. doi: 10.1007/s13593-019-0562-6

Givnish, T. J. (1982). On the adaptive significance of leaf height in forest herbs. *Am. Nat.* 120, 353-381. doi: 10.1086/283995

Godin, C., and Sinoquet, H. (2005). Functional-structural plant modelling. *New Phytol.* 166, 705-708. doi: 10.1111/j.1469-8137.2005.01445.x

Grace, J. B., and Tilman, D. (eds.). (1990). *Perspectives on Plant Competition*. New York, NY: Academic Press. doi: 10.1016/B978-0-12-294452-9.50005-9

Iooss, B., Da Veiga, S., Janon, A., and Pujol, G. (2017). sensitivity: Global Sensitivity Analysis of Model Outputs. R package version 1.15.0. Available online at: https://CRAN.R-project.org/package=sensitivity

Isbell, F. (2015). Agroecology: agroecosystem diversification. *Nat. Plants* 1, 1-2. doi: 10.1038/nplants.2015.41

Johnson, M. E., Moore, L. M., and Ylvisaker, D. (1990). Minimax and maximin distance designs. *J. Stat. Plan. Inference* 26, 131-148. doi: 10.1016/0378-3758(90)90122-B

Jones, D. R., Schonlau, M., and Welch, W. J. (1998). Efficient global optimization of expensive black-box functions. *J. Glob. Optim.* 13, 455-492. doi: 10.1023/A:1008306431147

Lecarpentier, C. (2017). *Plasticité de l'architecture aérienne du blé en réponse à la compétition pour la lumière au sein de cultures pure ou d'associations variétales : caractérisation expérimentale et développement d'un modèle 3D*. Thèse de doctorat, Université Paris-Saclay.

Lecarpentier, C., Barillot, R., Blanc, E., Abichou, M., Goldringer, I., Barbillon, P., et al. (2019). WALTer: a three-dimensional wheat model to study competition for light through the prediction of tillering dynamics. *Ann. Bot.* 123, 961-975. doi: 10.1093/aob/mcy226

Litrico, I., and Viole, C. (2015). Diversity in plant breeding: a new conceptual framework. *Trends Plant Sci.* 20, 604-613. doi: 10.1016/j.tplants.2015.07.007

Long, S. P., ZHU, X. G., Naidu, S. L., and Ort, D. R. (2006). Can improvement in photosynthesis increase crop yields? *Plant Cell Environ* 29, 315-330. doi: 10.1111/j.1365-3040.2005.01493.x

Louarn, G., Barillot, R., Combes, D., and Escobar-Gutiérrez, A. (2020). Towards intercrop ideotypes: non-random trait assembly can promote overyielding and stability of species proportion in simulated legume-based mixtures. *Ann. Bot.* 126, 671-685. doi: 10.1093/aob/mcaa014

Louarn, G., and Song, Y. (2020). Two decades of functional-structural plant modelling: now addressing fundamental questions in systems biology and predictive ecology. *Ann. Bot.* 126, 501-509. doi: 10.1093/aob/mcaa143

Ma, J., Tu, Y., Zhu, J., Luo, W., Liu, H., Li, C., et al. (2020). Flag leaf size and posture of bread wheat: genetic dissection, QTL validation and their relationships with yield-related traits. *Theor. Appl. Genet.* 133, 297-315. doi: 10.1007/s00122-019-03458-2

Maézieux, E. (2012). Designing cropping systems from nature. *Agron. Sustain. Dev.* 32, 15-29. doi: 10.1007/s13593-011-0027-z

Marrel, A., Iooss, B., Laurent, B., and Roustant, O. (2009). Calculations of sobol indices for the gaussian process metamodel. *Reliab. Eng. Syst. Saf.* 94, 742-751. doi: 10.1016/j.ress.2008.07.008

McKay, M. D., Beckman, R. J., and Conover, W. J. (1979). Comparison of three methods for selecting values of input variables in the analysis of output from a computer code. *Technometrics* 21, 239-245. doi: 10.1080/00401706.1979.10489755

Midmore, D. J., Cartwright, P. M., and Fischer, R. A. (1984). Wheat in tropical environments. II. Crop growth and grain yield. *Field Crops Res.* 8, 207-227. doi: 10.1016/0378-4290(84)90064-9

Montazead, G., Viole, C., Roumet, P., Rocher, A., Ecarnot, M., Compan, F., et al. (2020). Multifaceted functional diversity for multifaceted crop yield: towards ecological assembly rules for varietal mixtures. *J. Appl. Ecol.* doi: 10.1111/1365-2664.13735

Moon, P., and Spencer, D. (1942). Illumination from a non-uniform sky. *Trans. Illum. Eng. Soc.* 37:707712.

Muller, B., and Martre, P. (2019). Plant and crop simulation models: powerful tools to link physiology, genetics, and phenomics. *J. Exp. Bot.* 70, 2339-2344. doi: 10.1093/jxb/erz175

Munz, S., Graeff-Hönninger, S., and Henke, M. (2018). "Functional-structural plant model for testing the effect of maize architecture on hourly light distribution in strip-intercropping systems," in *2018 6th International Symposium on Plant Growth Modeling, Simulation, Visualization and Applications (PMA)*, 115-119. China, Hefei:IEEE. doi: 10.1109/PMA.2018.8611609

Nininemets, Ü. (2010). A review of light interception in plant stands from leaf to canopy in different plant functional types and in species with varying shade tolerance. *Ecol. Res.* 25, 693-714. doi: 10.1007/s11284-010-0712-4

Nix, H. A. (1976). *Climate and Crop Productivity in Australia. Climate and Rice*. Los Baños: IRRI, 495-507.

Passioura, J. B. (1996). Simulation models: science, snake oil, education, or engineering?. *Agron. J.* 88, 690-694. doi: 10.2134/agronj1996.00021962008800050002x

Perez, R. P., Dauzat, J., Pallas, B., Lamour, J., Verley, P., Caliman, J. P., et al. (2018). Designing oil palm architectural ideotypes for optimal light interception and carbon assimilation through a sensitivity analysis of leaf traits. *Ann. Bot.* 121, 909-926. doi: 10.1093/aob/mcx161

Perez, R. P., Pallas, B., Le Moguédec, G., Rey, H., Griffon, S., and Caliman, J. P. (2016). Integrating mixed-effect models into an architectural plant model to simulate inter-and intra-progeny variability: a case study on oil palm (Elaeis guineensis Jacq.). *J. Exp. Bot.* 67(15), 4507-4521. doi: 10.1093/jxb/erw203

Pradal, C., Dufour-Kowalski, S., Boudon, F., Fournier, C., and Godin, C. (2008). OpenAlea: a visual programming and component-based software platform for plant modelling. *Funct. Plant Biol.* 35, 751-760. doi: 10.1071/FP08084
Pradal, C., Fournier, C., Valduriez, P., and Cohen-Boulakia, S. (2015). OpenAlea: scientific workflows combining data analysis and simulation. In: Proceedings of the 27th International Conference on Scientific and Statistical Database Management. Gupta Amarnath, Rathbun Susan (ed.). New York: ACM, 6 p. ISBN 978-1-4503-3709-0. doi: 10.1145/2791347.2791365

Pugnaire, F., and Valladares, F. (eds.). (2007). Functional plant ecology. Boca Raton: CRC press. doi: 10.1201/9781420007626

R Core Team (2018). *R: A Language and Environment for Statistical Computing*. Vienna, Austria: R Foundation for Statistical Computing. Available online at: https://www.R-project.org/.

Reiss, E. R., and Drinkwater, L. E. (2018). Cultivar mixtures: a meta-analysis of the effect of intraspecific diversity on crop yield. *Ecol. Appl.* 28, 62–77. doi: 10.1002/eap.1629

Reynolds, M., Foulkes, J., Furbank, R., Griffiths, S., King, J., Murchie, E., et al. (2012). Achieving yield gains in wheat. *Plant Cell Environ.* 35, 1799–1823. doi: 10.1111/j.1365-3040.2012.02588.x

Roustant, O., Ginsbourger, D., and Deville, Y. (2012). DiceKriging, DiceOptim: two R packages for the analysis of computer experiments by kriging-based metamodeling and optimization. *J. Stat. Softw.* 51, 1–55. doi: 10.18637/jss.v051.i01

Sacks, J., Welch, W. J., Mitchell, T. J., and Wynn, H. P. (1989). Design and analysis of computer experiments. *Stat. Sci.* 4, 409–423. doi: 10.1214/ss/1177012413

Sarlikioti, V., de Visser, P. H., Buck-Sorlin, G. H., and Marcelis, L. F. M. (2011). How plant architecture affects light absorption and photosynthesis in tomato: towards an ideotype for plant architecture using a functional-structural plant model. *Ann. Bot.* 108, 1065–1073. doi: 10.1093/aob/mcr221

Sobol, I. M. (1993). Sensitivity estimates for nonlinear mathematical models. *Math. Modell. Comput. Exp.* 1, 407–414.

Sparkes, D. L., Holme, S. J., and Gaju, O. (2006). Does light quality initiate tiller death in wheat?. *Eur. J. Agron.* 24, 212–217. doi: 10.1016/j.eja.2005.08.003

Stomph, T., Dordas, C., Baranger, A., de Rijk, J., Dong, B., Evers, J., et al. (2020). Designing intercrops for high yield, yield stability and efficient use of resources: are there principles? *Adv. Agron.* 160, 1–50. doi: 10.1016/bs.agron.2019.10.002

Streit, K., Henke, M., Bayol, B., Cournède, P. H., Sievänen, R., and Kurth, W. (2016). “Impact of geometrical traits on light interception in conifers:

analysis using an FSPM for Scots pine,” in 2016 *IEEE International Conference on Functional-Structural Plant Growth Modeling, Simulation, Visualization and Applications (FSPMA)*, 194–203. China, Qingdao:IEEE. doi: 10.1109/FSPMA.2016.7818307

Valladares, F., and Niinemets, Ü. (2008). Shade tolerance, a key plant feature of complex nature and consequences. *Annu. Rev. Ecol. Evol. Syst.* 39:173506. doi: 10.1146/annurev.ecolsys.39.110707.173506

Vidal, T., Gigot, C., de Vallavieille-Pope, C., Huber, L., and Saint-Jean, S. (2018). Contrasting plant height can improve the control of rain-borne diseases in wheat cultivar mixture: modelling splash dispersal in 3-D canopies. *Ann. Bot.* 121, 1299–1308. doi: 10.1093/aob/mcy024

Vos, J., Evers, J. B., Buck-Sorlin, G. H., Andrieu, B., Chelle, M., and De Visser, P. H. (2010). Functional-structural plant modelling: a new versatile tool in crop science. *J. Exp. Bot.* 61, 2101–2115. doi: 10.1093/jxb/erp345

Zhang, B., and De Angelis, D. L. (2020). An overview of agent-based models in plant biology and ecology. *Ann. Bot.* 126, 539–557. doi: 10.1093/aob/mcaa043

**Conflict of Interest:** The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

**Publisher's Note:** All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Copyright © 2021 Blanc, Barbillon, Fournier, Lecarpentier, Pradal and Enjalbert. *This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.*