Estimating the Genetic Parameters of Flowering-Time Related traits in a Miscanthus Sinensis population tested with a staggered-start design

There is a growing interest in the cultivation of miscanthus on marginal land, but biomass yields are much lower there than on good farming land. Therefore, understanding what causes such instability is of primary interest for breeding later-owering but stable miscanthus genotypes. Our objectives were to estimate the genetic parameters -genetic variance and genetic heritability- and genetic correlations for owering-time-related traits in a biparental Miscanthus sinensis diploid population, and to divide the year effect into both age and climate effects using a staggered-start design. The population was established with single plants organized in a staggered-start design and consisted of two genotype groups established twice, in 2014 and 2015, with a total of 159 genotypes and 82 common genotypes between the two groups. All plants were extensively phenotyped for different panicle and anther emergence traits in 2018 and 2019. All traits were delayed by about 20 days in 2019 compared to 2018, which was explained by climatic conditions that occurred before the oral transition, mainly by a 3°C decrease in maximal and minimal temperatures. When dividing the year interaction effect, the genotype × climate interaction was much higher than the genotype × age interaction. The climate effect not only caused a delay in the owering time but also involved differential genotype behavior through climate × genotype interactions, which increased the corresponding genotype × climate interaction variance compared to the genotype × age interaction variance: the climate effect decreased the genetic parameters for all owering-time related traits, up to 20 % for broad-sense heritability. Interestingly, all traits responded similarly to the climate effect, excepting the interval between the start of panicle emergence and that of anther appearance, for which the correlation coecients were lower due to signicant climate interactions, compared to genotype × age interactions. Therefore, M. sinensis breeding for owering-time related traits must be conducted under contrasted climates in order to select more stable genotypes.


Introduction
Biomass is expected to play a major role in the energy transition, which involves switching to a safe and sustainable lowcarbon economy, in order to address the growing challenges due to depleting fossil resources and climate change [1]. In their review, Gabrielle et al. [1] reported that a small proportion of European renewable energy is derived from dedicated bioenergy crops. Most of them are perennials grown to generate electricity and heating, with the most frequent species being miscanthus, willow, reed canary grass, and poplar. In addition, some biomass crops are exploited in response to environmental pollution [2].
Miscanthus is a perennial C4 rhizomatous grass originating from eastern and southern Asia [3]. Notably, Miscanthus × giganteus was introduced in the 1930s in Europe from Japan by the Danish nurseryman, Aksel Olsen [4], and is known as a natural interspeci c hybrid between Miscanthus sacchari orus and Miscanthus sinensis [5]. Miscanthus cropping relies on M. × giganteus clones and has thus been used since 1983 [6]. It was rapidly identi ed as providing promising lignocellulosic biomass due to its high biomass yield per area [3], low nitrogen needs [7] and capability to recycle its nitrogen [8]. However, M. × giganteus is sterile and presents genetic uniformity. This crop currently relies on a single clone, which may cause risks that have been recognized in case of disease, pest, or climate pressure [9]. As early as 1997, Greef et al. [10] reported that genetic diversity was very low within a pool of 31 European M. × giganteus accessions [11]. concluded that 27 American and 6 European accessions were derived from a single clone by vegetative propagation. In contrast, M. sinensis exhibits huge genetic variability, with major genetic groups originating from Asia [3,12,13], which offers numerous opportunities to enlarge the varietal offer.
To produce biomass, crops require the longest possible vegetative growth period because the shift from vegetative to reproductive growth prevents or stops biomass accumulation. This shift is of particular interest on marginal land where biomass yields are very often much lower than on good agricultural land, and where Wagner et al. [14] found that economic sustainability is limited by the biomass yield. In miscanthus, several previous studies using different miscanthus germplasm found that late-owering was associated with higher biomass production. Zub et al. [15] compared M. sinensis, M. sacchari orus, and M. × giganteus and showed that the latest-owering genotypes produced more biomass than the earliest-owering ones. Clifton-Brown and Lewandowski [16] observed that some non-hybrid genotypes of M. sinensis had shorter growing seasons and lower plant heights due to their earlier owering, and their biomass yields were generally lower than those of later-owering hybrids. Regarding the owering diversity of miscanthus, Jensen et al. [17] compared 244 genotypes of two miscanthus species (M. sinensis, M. sacchari orus) and 8 inter-speci c hybrids including M. × giganteus. They found that M. sinensis clones were the earliest ones to ower, but showed the greatest diversity in terms of owering onset, which highlights the importance of this species for breeding new varieties. It also shows the need for new lateowering M. sinensis to increase their biomass production.
The genetic parameters, genetic variance and heritability of traits related to owering time in miscanthus are essential to breed such new M. sinensis varieties. Slavov et al. [18] observed the date on which the rst ag leaf emerged in a population of 138 M. sinensis plants and found that this trait was highly heritable (broad-sense heritability = 0.83). Gifford et al. [19] observed the day of the year of the heading date and 50% anthesis for two successive years in a population of 221 progenies from a cross between two M. sinensis, and found that these traits were highly heritable, with a broad-sense heritability of 0.80 and 0.67 for the heading date and 0.77 and 0.70 for the 50% anthesis. Dong et al. [20] estimated the broad-sense heritability of two diploid F1 M. sinensis populations for the date on which owering rst began (with opened orets, elongated stamens, and pollen shedding) and found high heritability levels (0.76 and 0.85). In all these studies, the owering-related traits were limited to a single trait except in the work of the last authors who studied several traits, but only one related to owering, i.e. the date on which owering rst began.
So far, the climate effect on genetic parameters has not been speci cally addressed in miscanthus. In the present study, the objectives were hence to estimate the genetic parameters and genetic and phenotypic correlations for traits related to owering, with a speci c focus on the year effect. We studied a biparental diploid M. sinensis population using a staggeredstart design [21] established twice, in 2014 and 2015, and extensively phenotyped in 2018 and 2019 for owering-time related traits. The power of such a design in which the population was established twice makes it possible to divide the "growing year" effect into "age" effect and "climate" effects [22,23]. Therefore, we paid attention to what caused the owering-time related traits to be delayed in 2019 compared to 2018. What caused the genetic parameters to change? We predicted that the climate effect would be different from the age effect, which sheds new light on the study of genetic parameters and correlations in M. sinensis.

Materials And Methods
Trial conditions and experimental design The eld trial was established in 2014 at the INRAE unit in Estrées-Mons (49° 72' N, 3° 00' E) in northern France, relying on the frame of the Biomass For the Future (BFF) project. It was repeated twice in 2014 and in 2015 using a staggered-start design [21]. The eld soil type was characterized as deep loam (Ortic luvisol, FAO, classification). The trial received no nitrogen input and weeds were regularly manually removed. Meteorological data were recorded daily at a short distance from the trial (about 1 km) and downloaded from the CLIMATIC database (INRAE AGROCLIM, 2020).

Population and experimental design
The population consisted of the F1 progeny of 159 individuals generated from an intraspeci c cross between two diploid M. sinensis varieties, "Silberspinne" and "Malepartus". Plantlets were initially produced from the seeds and grown in a greenhouse in order to provide clonal replicates of each genotype using in vitro propagation method [24]. The population was established with single plants organized in a staggered-start design [23]. An initial group of genotypes (coded G1), that contained a total of 111 genotypes and included the two parents, was established in 2014 using an initial incomplete block design [25]. A second group of 130 genotypes (coded G2) included the two parents as well and was established in 2015 using a second incomplete block design. 82 genotypes were common to both groups and established twice in 2014 and 2015 to allow the partition of the year effect into age and climate effects. The staggered-start design comprised a total of were grown in each block while the G1 remaining genotypes were randomly spread into four blocks on average and the G2 genotypes into at least three blocks. The two groups of genotypes, which were adjacent in the eld, were separated using 75 plants (15 x 5) and surrounded by 180 border plants, to guarantee similar competition effects. Individuals were randomly planted using the Excel package « BICA » wh,ich was based on VBA programming language and developed for the randomization of incomplete blocks (JC Bastien, personal communication). A density of 1 plant/m² was applied at regular intervals of 1 m × 1 m between the plants.

Phenotyping of traits related to owering
The phenotypic observation of the individuals within the population for owering-related traits was conducted in 2018 and 2019. The following traits were recorded on each plant and targeted stems belonging to the higher part of the plant canopy (above 80%): the start of panicle emergence (SPE), which corresponds to the date on which the rst panicle was visible at least 1 cm above the plant stem; anther appearance (AA), which corresponds to the date on which anthers appeared and began to shed pollen -this observation was carried out on the rst panicle that emerged in the plant; the mid-point panicle emergence (MPE), which corresponds to the date on which 50% of the panicles had emerged in the plant, and nally, the complete panicle emergence (CPE), which corresponds to the date on which over 80% of the panicles of the entire plant had emerged. In addition, one variable was calculated: the interval between heading and owering (IHF) that was de ned as the difference between the heading and owering dates observed on the rst panicle that emerged and owered. To reduce bias and increase accuracy in the data To guarantee data accuracy, phenotypic observations of the traits were made three times a week, starting with the rst plant that showed panicle emergence. The rst panicle to appear on the plant was marked and used for further determination of the date of anther appearance to facilitate determination and avoid sampling bias. In addition, other traits were observed in the population: the emergence date (ED), the senescence date (SD) that corresponded to the date on which the rst stem cohort had become completely senescent, and the date of the rst winter frost (FD).
To further calculate weather indicators, extensive observations were conducted on Malepartus, one of the two population parents, in order to determine the date of the oral transition. This clone was established in 2014 in a speci c trial not far from the population trial, at a plant density of 2 plants/m² (see Leroy et al., submitted to BioEnergy research). The sampling was carried out on three plots of about 75 m² each. The oral transition corresponded to the differentiation of the terminal vegetative apex into a panicle and was observed by dissecting 30 apices per day on average from May 24 th onwards. It was determined as the date of the rst oral transition appearance, which was expected to correspond to the oral transition of the rst stem that showed a panicle above the canopy.
Based on previous results obtained in maize, leaf appearance rate and duration responded to temperature in approximately the same way as for primordia and were constant across temperature regimes, at least for the appearance of the rst 12 leaves [26]. According to Hesketh and Warrington [27], the use of a degree-day model must predict the appearance rate of new organs, such as leaf primordia, tassel, and ears, as a function of air temperature or degree days. These results relating to constant degree-day durations were applied here in miscanthus to predict oral transition (referring to leaf primordia results in maize) and owering (referring to tassel and ear results in maize). We found that the sum of temperatures calculated from emergence to owering for the Malepartus clone with a base temperature between 6 and 8.5°C were almost equal in 2018 and 2019 when the base temperature used was of 8°C. This meant that the base temperature was close to 8°C for Malepartus and could be applied for further calculations. Therefore, growing-degree days (GDD) were estimated using 8°C as a base temperature for Malepartus and its offspring. Using the oral transition date determined from the apices observations in Malepartus, the ratio between the GDD for the emergence-to-oral transition period and the GDD for the emergence-to-owering period was calculated for Malepartus in 2018 and reached 0.17. This ratio was then applied to the population in order to estimate oral transition dates based on owering date observations. Determination of plant growth periods Five periods of the plant cycle were taken into account to calculate the weather indicators ( Table 1). The rst period corresponded to the pre-emergence phase (coded Period 1) and the calculations started at senescence, which occurs in the autumn, and precedes shoot emergence. The second period (Period 2) started with shoot emergence and ended with oral transition. The third period (Period 3) corresponded to the interval between oral transition and anthesis. The fourth period (Period 4) was from anthesis to senescence. Lastly, the fth period (Period 5) covered the duration from senescence to the date of the rst winter frost, the plant cycle being stopped by the rst frosts. Weather indicator description In order to explain the climate effect, the following meteorological indicators were determined: air minimum temperature (MinT), air maximum temperature (MaxT), air mean temperature (MeanT), soil minimum temperature (MinTs), precipitation (Prec), cumulated precipitation (Prec_S), mean humidity (MeanH), vapour pressure de cit (VPD), maximum VPD (VPD_M), Penman evapotranspiration (ETPP), photosynthetically active radiation (PAR), cumulated growing degree-days (CGDD). The summary of the meteorological indicators per growing period over the two years is shown in table 2.

Statistical analyses
Two different analyses for unbalanced datasets were performed, according to the "age" effect and the "climate" effect models ( Fig. 1). 159 genotypes were included in each of these analyses, with 82 in common between the G1 and G2 genotype groups.
An initial round of analyses of Linear Mixed Models consisted of a full statistical model which took into account block and spatial effects using the "breedR" package implemented in the R software [28]. It accounted for the environment heterogeneity within each trial among the blocks and diagnosed spatially distributed patterns of remaining residual variations using variograms. All variance components of the progeny were estimated by using the remlf90 function of the package, based on restricted maximum likelihood estimations. The equation of Model 1 was as follows: where Yijkl represents the phenotypic value measured on plant k of genotype i at age j in block l; µ is the overall mean; αi is the random effect of genotype i; βj is the xed effect of age j; (αβ)ij is the random interaction between genotype i and age j; δl is the effect of block l and εijkl is the random residual for plant k of genotype i at age j in block l. An autoregressive spatial component was also considered in the model based on x and y coordinates in each plot in order to partition the residual εijkl into a spatially dependent parameter, θikl, for plant k of genotype i in block l, and an independent remaining residual [29].
This full model was compared with two sub-models for each trait: a sub-model with no decomposition of the residual term into spatially dependent and independent effects and another sub-model with no block effect but with the decomposition of the residual term into spatially dependent and independent effects. The performance of the three models was compared on their Akaike information criterion, for which a lower value corresponded to a better performance. For all traits, the last submodel, which dropped the block effect but took into account the spatial effect, was found to have a slightly better performance than the full model or than the sub-model which dropped the spatial effect but took into account the block effect. Accordingly, it was retained afterwards, which simpli ed the models.
The "age" effect was assessed using the following Linear Mixed Model, tted for each studied trait and considering each year in which the progeny grew. It was referred to as "age effect model": where Y ijk represents the phenotypic value measured on the plant k of genotype i at age j, μ is the overall mean; α i is the random effect of the genotype i; β j is the xed effect of age j; (αβ) ij is the random interaction between genotype i and age j; θ ik is the random spatial effect for plant k of genotype i; and ε ijk is the random residual error for plant k of genotype i at age j. This modeling was carried out by using two subsets of the data (rectangles colored in blue in Fig. 1). The year 2018 was considered for G1 and G2 and the age effect was taken into account by 4-year old plants with G1 and 3-year old with G2.
Plants of different ages grew in the same climate during a same year (2018), which makes it possible to take into account the age effect in the same climate conditions. The age effect was also analyzed in 2019 in the same manner when the genotypes were one year older.
The performance of the full model including the interaction term was compared with the corresponding sub-model with no interaction by using the Akaike information criterion. Based on the AIC, the analysis showed that a model that included interaction effect was always more performant than the sub-model for all traits.
Then, the "climate" effect was assessed using the following Linear Mixed Model, for each age of the progeny. It was referred to as "climate effect model": where each term is similar to model (1), except that the age effect β j is replaced by γ l , the climate effect of year l and the interaction between genotype i and age j by the interaction between genotype i and climate l . It can be noted that the spatial effect term was a different estimate in comparison to the estimate from model 1. In this modeling, only 4-year old plants were considered and corresponded to G1 genotypes which grew in 2018 and G2 genotypes in 2019 (orange rectangle in Fig. 1). Plants of the same age grew in two different climate conditions, related to each year considered. The performance of the full model including the interaction term was tested in the same manner as previously and the results showed a higher performance of the model that included interaction effect for all traits.
According to these two models, genotypic broad-sense heritability was estimated as follows: For the age effect, where is the variance attributed to the genotype, is the variance of the genotype × age interaction, is the spatial variance and is the residual variance.
For the climate effect, Here, each term is similar to the previous formula, except for which is the variance of the genotype × climate interaction.
According to the two previous models, progeny-mean heritability [30,31] was also estimated for the age and climate effect models, respectively as: where J are the two ages considered, L are the two climates (of each year) and K is the mean number of repetitions per genotype.
Due to the unbalanced feature of the design, as all genotypes were not common between G1 and G2 (i.e. not established twice in 2014 and 2015), most previous models were applied twice as a second round of analyses was restricted to the 82 common genotypes between G1 and G2. It allowed a more precise comparison of the climate and age interaction effects, thus avoiding biases due to non-common genotypes between the comparisons.
Genetic and phenotypic correlations were then assessed for each trait, considering a given age or a given year. Genetic correlations were estimated as described by Howe et al [32]. To calculate phenotypic correlations, Pearson correlation coe cients were computed using R package "stats" and visualized using R package "corrplot" [33].
A Principal Component Analysis (PCA) was performed on the weather indicators determined in Table 2 using the R "FactoMineR" package. This PCA was followed by a clustering that was applied to the coordinates of the individuals in the two rst components of the PCA. These individuals corresponded to combinations between the periods and the two years. Within a period, the number of clusters was de ned when the two years of a given period were separated. Then, the climatic features of the separated combinations were pointed out.

Results
A delay of about 20 days was observed in 2019 compared to 2018 for all owering-time related traits When comparing the two years of observation, the owering time was signi cantly delayed -by around 20 days on average -in 2019 in comparison to 2018 (Table 3). In 2018, the median DOY values for the start of spike emergence (SPE), anther appearance (AA), mid-point of spike emergence (MPE), and completion of spike emergence (CPE) were 212, 218, 218, and 222, respectively. In 2019, the median DOY values for SPE, AA, MPE, and CPE were all delayed to 233, 238, 239, and 243, respectively. After discarding extreme odd values (corresponding to a 95% interval from 2.5-97.5% DOY quantiles), the range of the population was rather stable in both years for each owering trait. It was around 30 days excepting for CPE observed in the rst group of established genotypes (G1). The climate effect induced lower genetic parameters than the age effect The genetic and environmental variance traits were calculated for the owering-time related traits using the two models that were distinguished by the age and climate effects, and by the corresponding interaction effects with the genotype (models 1 and 2). These traits included SPE, AA, and MPE, CPE, and IHF (Fig. 2).
Considering the climate effect model, an initial dataset was used, which corresponded to G1 tested in 2018 and G2 in 2019 ( Fig. 2a): all these plants presented the same plant age of 4 years old, but grew in different years, i.e. different climatic conditions, which made it possible to test the climate effect on the variance components. A second dataset was also used as it corresponded to the previous one which was restricted to the common genotypes between G1 and G2 (Fig. 2b). The results showed that genetic variance most in uenced each owering trait (Fig. 2a, b). The variance corresponding to the interaction between genotype and climate was substantial, while the residual variance was lower than the two previous ones. Therefore, the broad-sense heritability estimations for all owering-time related traits were moderate and ranged from 0.55 to 0.62, whereas the progeny-mean heritability estimations were higher than expected, and ranged from 0.77 to 0.82. The second dataset, which was restricted to the common genotypes between G1 and G2, provided similar trends. Therefore, estimates of heritability were relatively similar for balanced designs using these common genotypes compared to the unbalanced dataset (i.e. the rst dataset): broad-sense heritability estimations for all owering-time related traits were moderate and ranged from 0.53 to 0.61 while the progeny-mean heritability estimations were high and varied from 0.79 to 0.85.
Regarding the age effect model, two datasets were observed: G1 (established in 2014) and G2 (established in 2015) that were observed either in 2018 (dataset corresponding to Fig. 2c) or in 2019 (dataset corresponding to Fig. 2e). G1 and G2 corresponded to the plant groups at different ages but grew during the same year. In 2018, when the plant age was 4 for G1 and 3 for G2, all owering-time related traits showed large genetic variances but displayed small variances for the interaction of the genotype with age, as well as for the spatial effect (Fig. 2c). Due to these small environmental component variances, the broad-sense heritability estimations for all owering-time related traits were higher than those observed with the climate effect model and ranged from 0.73 to 0.77, while the progeny-mean heritability estimations were very high and varied from 0.9 to 0.92. As for the climate effect model, the same trend was observed when the analysis was restricted to the common genotypes as the same estimates of genetic parameters were obtained (Fig. 2d).
In 2019, older plant ages were analyzed and corresponded to plant age 5 for G1 and plant age 4 for G2 (Fig. 2e). In comparison to 2018, there was a noticeable shift in the total variance and genetic variance towards lower values while the genetic variance remained the largest contribution to the total variance. Larger interaction effects between genotype and age were also noticeable in 2019, in particular for SPE and AA. Consequently, broad-sense heritability was lower in 2019 than in 2018 and varied from 0.61 to 0.69, and so were progeny-mean heritability estimations that ranged from 0.81 to 0.83 for all owering-time related traits.
When comparing Fig. 2a with 2c in order to better understand the difference between the climate and age effects on owering-time related traits, the results showed that the total variance of each owering trait was almost similar under the effects of climate and plant age (observed in 2018). For all owering-time related traits, the interaction variance between genotype and climate (with model 2) was much higher than the interaction variance between genotype and plant age (with model 1): the climate increased these interaction variances from three to sevenfold according to the trait. The interactions between genotype and climate (on the x-axis, Fig. 3) were indeed higher than interactions between genotype and age (on the y-axis). In contrast, the spatial and residual variances were rather stable between the climate effect model and the age effect modelled in 2018. To avoid biases in these comparisons, the datasets were restricted to the common genotypes between G1 and G2 (Fig. 2b, d) and provided very similar results to those obtained from the whole datasets by comparing Fig. 2a to Fig. 2c. When comparing the climate effect model to the age effect modelled in 2019 (Fig. 2a, e), the trend was the same, although the total variance of the age effect modelled in 2019 was lower and the increase due to the climate was of a lesser extent: nevertheless, the interaction variance between the genotype and the climate for all owering-time related traits was still higher than the interaction variance between the genotype and the age for all owering-time related traits modelled in 2019, but increased from two to fourfold.
High positive genetic correlations where observed among most owering-time related traits when the plant age or climate effects were modelled Genotypic and phenotypic correlations among the owering-time related traits were estimated using climate and age effects models (reported above and below the diagonal of Fig. 4, respectively). The results showed that the phenotypic correlation coe cients for all owering-time related traits in all data sets were lower than the corresponding genotypic correlation coe cients, but all values coincided. Therefore, the present study focused on comparing correlations at the genetic level. The same traits as before were studied (SPE, AA, MPE, CPE, and IHF).
In 2018, the estimates found under the age effect modelling process showed strong and positive genetic correlations among SPE, AA, MPE, and CPE (Fig. 4). For instance, CPE showed the highest correlation (0.99) with MPE among these higher correlated owering-time related traits, and the smallest correlation coe cient (0.93) was observed for AA with SPE, MPE, and CPE. IHF showed a positive and signi cant correlation (0.45) with AA and a positive but not signi cant correlation with other owering-time related traits (0.08 to 0.18).
Similarly, in 2019, the age effect model of the population corresponded to older plants, the genetic correlation estimates remained strong among SPE, AA, MPE, and CPE. IHF showed a shift toward negative and weak values with SPE, MPE, and CPE. A positive (but not signi cant) correlation with AA can also be noted. In addition, the genetic correlation estimates relating to the climate model at plant age 4 also performed a very similar result as the age effect model of the population.
Comparing the genetic (above the diagonal of Fig. 4) and phenotypic (below the diagonal of Fig. 4) correlations of all traits, the result provided a very similar information, where the genetic correlation coe cient coincided with the phenotypic correlation coe cient, and the genetic correlation coe cient a little larger than the phenotypic correlation coe cient.
Climatic conditions that occurred before the oral transition contributed to the differences between the two years Due to the previously highlighted high climate effect, a Principal Component Analysis (PCA) was conducted to identify the weather indicators that contribute to the climate effect for the owering-time related traits observed in the F1 population. The rst two principal components (PCs) accounted for a large proportion of the total variation (93.0 %), with PC1 explaining 84.2 % and PC2 explaining 8.8 %. (Fig. 5).
The PC1 rst component showed that 8 weather indicators had high positive loadings (> 0.9), on the right-hand position of Fig. 5a: indicators related to temperatures (MinT, MeanT, MaxT for the air and MinTs for the soil), vapour pressure de cit (VPD and maximal VPD_M), Penman evapotranspiration (ETPP), and photosynthetically active radiation (PAR). Moreover, these indicators were highly correlated with each other. PC1 also indicated that MeanH had a high negative loading (<-0.9), on the left-hand position of Fig. 5a. PC2 was explained by a single variable that was related to the sum of precipitation calculated for each period and each year (Prec_S).
The clustering of the coordinates on the rst two PCA components of the individuals (i.e. periods according to year) clearly divided them into seven distinct groups (Fig. 5b). Most groups included both years for the periods, excepting two periods: Period2 which corresponds to the shoot emergence to oral transition phase (determined using the Malepartus clone) and Period4 which corresponds to the owering to senescence phase. Interestingly, Period2 gathered climatic events which preceded owering: it can be noted that both years received little precipitation while the temperature was greater in 2018 compared to 2019. During period 2, all indicators related to temperature were higher for 2018 compared to 2019. which points to warmer weather conditions in 2018: maximum temperature of 18.3°C against 15.3 ° and minimum temperature of 8.2°C against 5.2°C. The duration of period 2 in 2018 was 35 days, which was 30 days less than the same period in 2019, but it was rather constant in both years when estimated in cumulative growing degree-days (209.2 for 2018 and 203.1 for 2019).

Discussion
In the present study, the owering related traits were extensively phenotyped within an F1 diploid population of M. sinensis for two years (2018 and 2019) and included the start of spike emergence (SPE), anther appearance (AA), mid-point of spike emergence (MPE), completion of spike emergence (CPE) as well as the interval between the start of spike emergence and anther appearance (IHF). All owering-time related traits differed signi cantly within the population according to the genotypes and years. Interestingly, the owering time was signi cantly delayed -by around 20 days in 2019 compared to in 2018 (Table 3). A climate effect on the genetic parameters was also observed but traits responded similarly to this effect. Therefore, the discussion focused on three points: 1) the owering time delay may be explained by climatic conditions that occurred before the oral transition, 2) the genotype × climate variances were greater than the genotype × age variances for all owering-time related traits, which implied a decrease in heritability and, 3) although the existence of these previously signi cant effects of the climate, most owering-time related traits of the population behaved similarly to changes in climatic conditions. The delay in the owering time was explained by climatic conditions that occurred before the oral transition Two periods were highlighted based on an extensive analysis of weather indicators according to different periods of the plant cycle: a principal component analysis followed by a clustering of the rst two component coordinates of the periods for both years showed that Period 2 and Period 4 were subjected to different weather conditions in 2018 and 2019. Period 2, which corresponds to the shoot emergence to oral transition period, is interesting because it concerned climatic events that occurred before owering, which may have an impact on the owering time. In particular, cooler temperatures during the period preceding oral transition (a 3°C decrease in maximal and minimal temperatures) contributed to the 20-day delay for the owering time between the two years. In this study, we also observed that all weather indicators for the subsequent period from oral transition to owering were very similar in both years, especially those related to the temperature that were almost identical. Therefore, differences in the temperatures observed during the emergence to oral transition period accounted for the delay in the owering time in 2019. A signi cant delay was also observed between years (including both age and climate effects) in the mean annual owering day in M. sinensis [17]. Growing seasonal precipitation, degree days, and the temperature calculated for the whole growing season partially accounted for this delay. In our study, we highlighted that the period preceding owering, and in particular the period preceding the oral transition stage is critical and has an impact on the owering of miscanthus. Regarding cooler conditions, Clifton-Brown et al. [9] found that, in a multienvironment experiment, M. sinensis owered later in the growing season at lower-temperature experimental sites. Nunn et al. [34] compared the growth of 15 miscanthus germplasm types, including M. sinensis, M. sacchari orus and, M. × giganteus, in six European countries with different climates and found that the growing season was accelerated in the hottest experimental site, with earlier emergence, owering, and senescence than at other experimental sites. Here, we noticed that a three-degree decrease during the plant cycle preceding the oral transition can imply a substantial delay in the owering time by about 20 days, since the cumulative growing degree-days of that period were constant between the two years. Regarding the cumulative growing degree-days, according to Deuter [25], M. sinensis is a plant that will ower upon an accumulation of a certain number of degree days, which was veri ed with our population. But our experiment did not allow for the testing of an additional effect of photoperiod on the owering time in M. sinensis.
This highlights the importance of the emergence to oral transition period in explaining differences in earliness/lateness from one year to another. In particular, the differences in temperature indicators observed during that period contributed to the year effect on the owering time in the M. sinensis population observed in the present study. This involved a substantial delay in the owering time, moreover within the same site, as observed here.
The genotype × climate variances were greater than the genotype × age variances for all owering-time related traits, which implied a decrease in the genetic parameters In the present study, we modelled the age effect (model 1) and the climate effect (model 2) and the corresponding interaction effects on genotypes, using a staggered-start design. In both modelling processes, the spatial effects were adjusted and contributed to obtain more accurate estimates of the traits by decreasing the residual term. In all cases and for all owering-time related traits, the genetic variance contributed to the largest part of the variance, which implied rather high heritability levels (broad-sense heritability above 0.55 and progeny-mean heritability above 0.77). Interestingly, the variance due to the genotype × climate interaction was much greater than that of the variance due to the genotype × age interaction (observed for the 4 − 3 year and the 5 − 4 year ages), thus highlighting the e ciency of the experimental staggered-start design in dividing the year effect into age and climate effects. Therefore, the climate not only implied a delay in the owering time but also involved the differential behaviors of the genotypes, which led to higher genotype × climate interaction variances and consequently lower genetic parameters (genetic variances and heritability) for all owering-time related traits.
Such an environmental in uence on phenotypic traits has also been reported in previous studies. Clark et al. [35] estimated the heritability of 569 M. sinensis genotypes and 9 M. × giganteus for 14 yield component traits at six experimental sites, respectively. They found that the heritability estimates for the same trait were signi cantly different across experimental sites due to different environmental conditions. But in the present study, we also highlighted such environmental effect in a single location. Hazard et al. [36] observed that perennial ryegrass established under outdoor conditions had lower heritability estimates than ryegrass grown in greenhouses, and they suggested that the different temperature conditions during the growth processes may result in heritability changes. Regarding the year effect, Gifford et al. [19] estimated the heritability of 13 phenotypic traits associated with biomass productivity for 221 M. sinensis and found an effect of the year (which includes both age and climate effects). Segura et al. [22] observed that the heritability of most architecture-related traits showed highly signi cant age effects by using a staggered-start design on apple trees, which was interpreted as the result of signi cant changes in climatic conditions. When modelling the age effect on owering-time related traits in our study, signi cant changes were observed between 2018 and 2019, which can be interpreted as a result of signi cant different climatic conditions as well, in particular during the plant cycle prior to the oral transition.
The results above highlighted the e ciency of the staggered-start design in dissecting the year effect into age and climate effects in miscanthus. Consequently, we could observe that the genotype × climate interaction effects were much greater than the genotype × age interaction effects. Therefore, the climate not only implied a delay in the owering time but also implied the differential behavior of the genotypes leading to higher genotype × climate interaction variances, and consequently lower genetic parameters (genetic variances and heritability). Such a result would provide further evidence that climate in uences owering-time related traits, resulting in an up to 20 % decrease in broad-sense heritability according to the owering trait. Therefore, the breeding of such traits has to be conducted under contrasted climatic conditions in order to take into account the interactions of the genotypes with changing climate conditions. Most owering-time related traits of the population behaved similarly to changing climate conditions In order to better understand the effects of climatic conditions and plant age on the owering-time related traits in miscanthus, our study further analyzed the genetic and phenotypic correlations among all owering-time related traits (Fig. 3). The results showed that genetic and phenotypic correlations for owering-time related traits provided very similar information, where genetic correlation coe cients coincided with phenotypic correlation coe cients, the former being slightly larger than the latter. This indicates that the different owering-time related traits were mainly genetically determined. Genetic correlation is a component of phenotypic correlation, and the phenotypic correlation coe cient for two traits with high genetic correlation is largely determined by the genetic correlation coe cient, which is more reliable than phenotypic correlation because genetic correlation can eliminate the in uence of the environment. In this study, there were positive and strong correlations between the owering-time related traits of SPE, AA, MPE, and CPE under the climate effect and age effect models. Interestingly, the differences in the correlation estimates were very small when comparing the climate and age models (comparison of Fig. 4a to Fig. 4c or 4d) except for the variate IHF, which corresponds to the interval between the start of panicle emergence (SPE) and that of anther appearance (AA): the lowest values were indeed obtained when the climate was modelled in comparison to the age model, and this was due to smaller genotype × age interactions than genotype × climate interactions. It seemed that all traits behaved in the same way within the population when subjected to climate condition changes, excepting the interval between the start of panicle emergence and that of anther appearance. This interval tended to be shorter under warmer conditions such those observed in 2018, and showed large differences among genotypes. Such changes according to years, regrouping age and climate effects, have been observed by Gifford et al. [19] who compared the genetic and phenotypic correlations of the two years for biomass yield, plant height, and plant circumference observed on the 221 progeny of Miscanthus sinensis.
Although changes in climatic conditions implied a signi cant delay in the owering time and large variations in the genetic parameters of owering-related traits, most owering-time related traits of the population behaved similarly to varying climate conditions. The single exception was observed for the interval between the start of panicle emergence and that of anther appearance.

Conclusion And Prospects
In conclusion, the climatic conditions preceding the oral transition of M. sinensis led to changes in its owering time, and temperature was one of the main variables that affected the owering time observed on a diploid M. sinensis progeny. In our study, a staggered-start design was e cient in dividing the year effect into plant age and climate effects and made it possible to discover that the genotype × climate variance was much greater than the genotype × age variance for each owering trait. This implied that the differential behavior of genotypes led to genotype × climate interactions, which resulted in altered genetic parameters (genetic variances and heritability). In addition, most of the owering-time related traits observed in our population -start of panicle emergence, anther appearance, mid-point of emergence, and completion of emergence -behaved similarly to changing climatic conditions, excepting the interval between the panicle emergence and anther appearance. All the results above brought strong evidence of the effect of changing climatic conditions on oweringtime related traits in M. sinensis, moreover observed within a single location. This resulted in non-stable earliness/lateness of the genotypes, which may impact the stability of their biomass yield itself. Therefore, breeding for owering-time related traits in M. sinensis must be carried out under contrasted climatic conditions in order to consider the interaction of the genotypes with changes in climatic conditions and select more stable genotypes for the owering time. This is critical for targeting new stable varieties, especially on marginal land where yields are known to be not only lower but highly variable. This study provides new insights for the breeding of more stable genotypes in M. sinensis, which is an important consideration for the cultivation of miscanthus in new areas, in particular on marginal land. Figure 1 The staggered-start design was composed of the establishment of two groups of genotypes: G1 in 2014 and G2 in 2015. It was analyzed using an initial model which accounted for the age effect modeling per year (a) and a second model which accounted for the climate effect modeling at 4 year-old plants (b). In case (a), the year 2018 was considered for G1 and G2, with the age effect that was taken into account by 4-year old plants with G1 and 3-year old with G2. Plants of different ages grew in the same climate during a same year (2018) which makes it possible to take into account the age effect in the same climate conditions. The age effect was also analyzed in 2019 in the same manner when the genotypes were one year older.

Declarations
In case modeling (b), only 4-year old plants were considered: they corresponded to G1 genotypes which grew in 2018 and G2 genotypes in 2019. In this case, plants of the same age grew in two different climate conditions, related to each year considered.

Figure 2
Estimations of genetic and environmental variances and heritability for owering-time related traits considering age effect (model 1) or climate effect (model 2) models. The owering-time related traits included the start of panicle emergence (SPE), anther appearance (AA), mid-point panicle emergence (MPE), complete panicle emergence (CPE), and the interval between anther appearance and the start of panicle emergence (IHF). 5 datasets were modelled: the climate effect model using the whole dataset (on a) or common genotypes between G1 and G2 (b), the age effect model observed in 2018 using the whole dataset (c) or common genotypes between G1 and G2 (d) and nally, the age effect model observed in 2019 using the whole dataset (e). Vgeno stands for genotypic variance, Vg × a for interaction variance between genotype and plant age, Vg × c for interaction variance between genotype and climate, Vspat for spatial variance and Vres for residual variance Figure 3 Illustration of the higher interactions between genotype and climate (cl_2018 and cl_2019 on the x-axis) than between genotype and age (3-yr and 4-yr on the y-axis). BLUPs of the interactions were estimated with the 82 common genotypes. The owering-time related traits (SPE, AA, MPE, CPE, and IHF) were the same as in Fig. 2.