QTL Analysis of Five Morpho-Physiological Traits in Bread Wheat Using Two Mapping Populations Derived from Common Parents

Traits such as plant height (PH), juvenile growth habit (GH), heading date (HD), and tiller number are important for both increasing yield potential and improving crop adaptation to climate change. In the present study, these traits were investigated by using the same bi-parental population at early (F2 and F2-derived F3 families) and late (F6 and F7, recombinant inbred lines, RILs) generations to detect quantitative trait loci (QTLs) and search for candidate genes. A total of 176 and 178 lines were genotyped by the wheat Illumina 25K Infinium SNP array. The two genetic maps spanned 2486.97 cM and 3732.84 cM in length, for the F2 and RILs, respectively. QTLs explaining the highest phenotypic variation were found on chromosomes 2B, 2D, 5A, and 7D for HD and GH, whereas those for PH were found on chromosomes 4B and 4D. Several QTL detected in the early generations (i.e., PH and tiller number) were not detected in the late generations as they were due to dominance effects. Some of the identified QTLs co-mapped to well-known adaptive genes (i.e., Ppd-1, Vrn-1, and Rht-1). Other putative candidate genes were identified for each trait, of which PINE1 and PIF4 may be considered new for GH and TTN in wheat. The use of a large F2 mapping population combined with NGS-based genotyping techniques could improve map resolution and allow closer QTL tagging.


Introduction
Bread wheat (Triticum aestivum L.) is grown on more than 200 million hectares of land worldwide, and in 2020, the global production reached about 760 million ton [1]. Despite this, by 2050, wheat production may need to be increased by at least 50% relative to current levels [2,3]. Although crop yields continue to increase globally, climate change represents a tremendous challenge for achieving this objective. In fact, recent research suggests that some major crop yields have already stagnated or even been reduced by the impact of climate change [4]. Phenology genes also regulate the physiological development of wheat [5], and some morpho-physiological traits have been identified as effective in breeding drought-adaptive varieties [6]. Traits such as flowering time, plant height, tiller number, growth habit, flag leaf angle, and spike characteristics, which are important for increasing crop yield potential, are also functional in determining the adaptation of wheat to climate change [7]. Genetic improvement has contributed significantly to modify wheat architecture by reducing plant height through selection of specific alleles at Rht-1 loci [8] as well as by adjusting flowering time. It is well known that Ppd genes influence the plant height, by regulating the number of leaves emitted during the vegetative phase of wheat development, the number of tillers, and fewer spikelets per spike [9]. Giunta et al. [10] demonstrated a major role for the Vrn genes on the control of tillering, which was independent of any photoperiod or temperature effects. The wheat ideotype generated during the Green Revolution was also characterized by an erect juvenile growth habit, a reduced number of tillers, and a reduced flag leaf angle [11]. If, on the one hand, this made it possible to increase the investment of plants per unit of area, reducing the intraspecific competition, on the other hand, it made the crop more vulnerable to weeds, reducing interspecific competition. Unfortunately, the trade-off between productivity in weed-free situations and competitive ability is a main obstacle to the release of competitive cultivars [7]. Genetic linkage mapping in segregating populations is one of the traditional approaches to gain insights into the genetic control of key characteristics using quantitative trait locus (QTL) analysis and is a suitable method underlying the molecular marker-assisted selection (MAS) in wheat [12]. The commonly used mapping populations in plants include F 2 , backcross (BC), double haploid (DH), and recombinant inbred line (RIL) populations [13,14]. Among these, the F 2 population, using codominant markers, provides the most abundant genetic information, including effects such as epistasis and dominance [15]. Ferreira et al. [16], in a simulation study, showed that more accurate maps were obtained with F 2 codominant and RILs than with BC, DH, and F 2 -dominant populations.
The RIL populations have several advantages for use in QTL mapping, which have been described by several authors [27][28][29][30]. Multiple selfing processes can increase the number of recombination events [31], which results in a finer mapping of QTLs, but the most important aspect is that once RILs are established they can be repeatedly used for investigating the QTLs of various phenotypes under different environments [10]. However, at least six generations are required to obtain the RILs and they are not suitable to estimate dominance effects of mapped QTLs due to the absence of heterozygous genotypes [32]. Nonetheless, RILs have been extensively used in wheat to identify QTLs for yield [33,34] and other important agronomic traits such as plant height [35], flowering time [36][37][38], and tiller number [39][40][41][42][43]. Otherwise, studies conducted for mapping the juvenile growth habit have been very limited [44].
The knowledge previously acquired by using bi-parental populations and the current availability of high-throughput and cost-effective single nucleotide polymorphism (SNP) genotyping platforms have opened the way for understanding the genetic mechanisms of these important agronomic traits [45].
Stable QTLs are a prerequisite for their use in a marker-assisted breeding program. To validate the QTLs, it is necessary to confirm them in other mapping populations, or in different generations from the same crossing, or using the same RIL population evaluated in multiple locations and in multiple years [46]. QTLs with large effects should be detected consistently across generations, but the increased precision of the RILs should allow the detection of QTLs with smaller effects [16]. In addition, additional opportunities for recombination should improve the genetic resolution of linked QTLs and distinguish linked effects from pleiotropic ones in some instances [47].
In this study, we performed QTL analysis of five morpho-physiological traits, all potentially associated with the yield potential and adaptability in bread wheat, using the same bi-parental population at early (F 2 and F 2 -derived F 3 ) and late (RIL F 6 and F 7 ) stages of inbreeding. We aimed to improve the understanding of the genetic basis (including additive and dominance effect) by comparing the QTLs detected by both populations, which were grown at the same location and evaluated in two years each. In addition, based on a priori knowledge of genes with biological function regulating the investigated traits in wheat and in other species, we also searched for putative candidate genes in chromosome regions tightly linked to QTLs for those traits by exploiting the reference bread wheat genome.

Plant Material and Field Trials
Two mapping populations derived from common parents, Lankaodali and Rebelde, were developed at the Department of Agriculture and Forest Sciences, University of Tuscia, Viterbo, Italy, in collaboration with the CREA Research Centre for Cereal and Industrial Crops (CREA-CI), Foggia, Italy. The female parent, Lankaodali, is a Chinese bread wheat cultivar with very large kernels [48], low tillering capacity, and a medium flowering time. The male parent, Rebelde (Bologna/CH-710//Bologna), is an Italian hard winter wheat cultivar with excellent agronomic characteristics, medium-late flowering time, and high tillering capacity. The two cultivars were crossed for generating first an F 2 population; then, by advancing through single seed descent, a F 6 RIL population was generated.
A total of 176 F 2 individual plants and F 2 -derived F 3 families of 15 plants (hereinafter referred to as F 3 ) together with the parents were grown in a randomized block design with two replications during the two consecutive growing season 2014-2015 and 2015-2016 at CREA-CI (41 • 27 40.2" N 15 • 30 04.5" E). F 2 and F 3 seeds were individually hand-planted in single rows, 2 m length, and 0.5 m between rows.
The RIL population comprising 178 lines and the two parents were planted during the two consecutive growing seasons 2017-2018 and 2018-2019 (hereinafter referred to as F 6 and F 7 , respectively) at the same experimental site. To simulate a standard open field sowing, a standard sowing density corresponding to approximately 350 seeds per square meter was carried out in single rows, 2 m long and 0.20 spaced, in a randomized block design with two replications. Whole plants of the F 2 , F 3 , and the two RIL generations (F 6 and F 7 ) were manually harvested at the end of June 2015, 2016, 2018 and 2019, respectively. Crop management was performed using standard cultivation methods. Growing season precipitation and temperature data were collected from the meteorological station of the CREA-CI.

Phenotypic Evaluation
Traits were measured both on a row and individual plant basis as follows: Heading date (HD) was recorded as days from the 1st of April in all growing seasons when 50% of the plants in a row were at GS59 [49]. Plant height (PH, cm) of individual F 2 and F 3 plants and ten randomly selected plants from RILs were measured at maturity. The juvenile growth habit trait (GH) was estimated by measuring the tiller angle between the last developed tillers and the ground level with a protractor at the maximum tillering stage (GS25 to GS29) by following UPOV [44] guidelines: 1, erect; 3, semi-erect; 5, intermediate; 7, semi-prostrate; and 9, prostrate. Tillering was monitored throughout the study period. Ten plants randomly selected on each single row were marked to facilitate identification at harvest. This was necessary for the F 6 and F 7 RIL generation with dense planting. At maturity, the number of tillers bearing spikes for each plant was counted, together with the number of fertile spikes. The tillering data were used to calculate the maximum number of tillers (TTN) and the number of fertile tillers (FTN). All phenotypic data were processed in the original form without any transformation. Principal component analysis (PCA) was performed using the mean values of all phenotypic data with the function "prcomp" implemented in the factoextra v. 1.0.7 library under R environment v. 4.0.2. Pearson correlation coefficients (r) were calculated among traits using the "cor" function and the correlogram was constructed and visualized using the Corrplot package v. 0.84 implemented in R v. 4.0.2 [50]. Coefficient of variation was calculated as follows: (1) where "Std" represents the standard deviation and "M" represents the mean for a trait. Broad-sense heritability (H 2 ) was estimated as follows: where σg is genotypic component of variance, σgy is the variance explained by interaction between genotypes and year, σe is the variance of residual effects, τ is the number of the repetitions, and y is the number of the year. The Shapiro-Wilk (w) statistics was used to test the null hypothesis that the phenotypic data were normally distributed [51].

Linkage Maps Construction
Plant materials from F 2 , F 6 generations, and both parents were genotyped with the wheat Illumina 25K Infinium single nucleotide polymorphism (SNP) array that was developed at the Trait Genetics (available online: http://www.traitgenetics.com, accessed on 24 June 2020). Both datasets were filtered by removing missing data (call rate > 10%) and rare alleles (MAF < 1%). Markers showing significant (p < 0.01) segregation distortion were also discarded.
The two filtered datasets were used to construct the genetic maps (hereinafter referred to as F 2 and RIL map) using the package ASMap v. 0.4.1 [52]. Linkage groups were constructed using the MSTmap algorithm [53] based on the p-value threshold set to 1 × 10 −10 and 1 × 10 −12 for F 2 and RIL populations, respectively. The Kosambi function was used to calculate the distance between markers and convert the recombination frequency into map distance (cM). The recombination fraction between all pairs of markers was estimated by the "est.rf" function to verify if each marker was placed in the right chromosome and in the correct order as compared with the Wen's reference map [54]. According to the Wen's reference map, some chromosomes were flipped through the "flip.order" function, whereas different linkage groups belonging to the same chromosome were merged using "mergeCross" function. The colinearity between the F 2 and F 6 linkage maps was evaluated by the Spearman rank correlation coefficient calculated by the R function "cor.test". A comparison between the obtained maps was performed by the R package SOFIA v. 1.0 [55].

QTL Mapping Approaches
The QTLs were detected for each trait in both generations by performing genome-wide composite interval mapping (gCIM), using the R software QTL.gCIMapping.GUI v.2.0 [56]. The first step was the binning of the redundant markers, as these markers have the same segregation in the population and bring no additional information to QTL analysis. The selected model was fixed and a walk speed for genome-wide scanning set to 1 cM. The LOD score threshold was set to 3.0 for both F 2 :F 3 and RILs. A maximum likelihood (ML) function was used, the random seeds were set to 11,001, and the completing CIM in one neighborhood was set to "FALSE" for the F 2 :F 3 population only. Moreover, the effects of QTLs, which include additive and dominance effects, were estimated by empirical Bayesian methods. The QTLs were named based on the recommended rules for gene symbolization in wheat (available online at https://wheat.pw.usda.gov/ggpages/wgc/98/Intro.htm, accessed on 15 January 2021).

Identification of Putative Candidate Genes
Putative candidate genes were identified within the QTL confidence intervals based on the annotation of the Chinese Spring genome, according to the IWGSC RefSeq v1.0 [57] (https://wheat-urgi.versailles.inra.fr/Seq-Repository/Annotation, accessed on 1 February 2021). The acronym of genes was assigned using the UNIPROT database (https://www.uniprot.org/, accessed on 1 February 2021). The sequence of the PINE1 (LOC_Os12g42250) gene isolated in rice was retrieved by Gómez-Ariza et al. [58], whereas PIF4 (LOC_Os12g41650) sequence was searched on the rice genome database available at https://shigen.nig.ac.jp/rice/oryzabase/ (accessed on 10 February 2021). The rice nucleotide sequences were used for a TBLASTN search against the genomes of the durum cultivar Svevo (http://www.gramene.org/, accessed on 10 February 2021). Then, the physical map position of the best sequence hits was compared to that of the genes identified in the QTL intervals, to investigate the correspondence.
Finally, all genes found in the QTL intervals were associated with a metabolic pathway by trait using the software MAPMAN [59].

Phenotypic Evaluation and Correlation Analysis between Early and Late Generations
Mean values, ranges, coefficient of variation (CV), and broad-sense heritability (H 2 ) of the investigated traits were calculated for each of the four generations (Table 1). Significant differences between the mean values of the parents were observed in all generations for HD, GH, TTN, and FTN, except for PH in the late generations. The mean values of HD were higher in F 2 and F 3 than in F 6 and F 7 , while the opposite was observed for PH. F 2 had the lowest mean value for GH, F 2 and F 6 showed the highest values in TTN and FTN. For the five traits, the CV ranged from 10.13 (HD) to 71.94% (GH) in the early generations and from 14.74 (HD) to 37.93% (FTN) in the late generations. The H 2 values ranged from 0.17 to 0.76 and from 0.37 to 0.83 in the early and in the late generations, respectively. The H 2 values higher than 0.69 were observed for HD and PH in early and RIL generations, while GH had H 2 higher in the RILs (0.82) than in the F 2 /F 3 (0.30). The heritability of the TTN and FTN was low in both generations with the lowest value in the early generations (0.17 and 0.20 for TTN and FTN, respectively).
The Shapiro-Wilk test confirmed the normality of PH, TTN, and FTN data in F 2 and F 7 populations, with a right-skewed distribution for the last two traits. In addition, normality was also observed for the trait FTN in the F 3 and for the traits, GH and PH, in the F 6 generations. A bimodal distribution for HD was observed in early generations (F 2 and F 3 ) (Supplementary Figure S1).
The PCA analysis, based on the five morpho-physiological traits, showed a large variability within the four populations ( Figure 1). The first two dimensions explained 70.1%, 69.9%, 66.5%, and 63.6% of the total variance in the F 2 , F 3 , F 6 , and F 7 generations, respectively.

Construction of F2 and F6 Linkage Maps and Their Colinearity
After filtering, 3684 and 3958 SNP variants were retained for F2 and F6, respectively, and used to build the two genetic maps. Thirty-nine and thirty-eight linkage groups were identified in F2 and F6 populations, respectively. The smaller linkage groups were joined to form the twenty-one chromosomes of T. aestivum L., based on the information available in the consensus map of Wen et al. [54]. The two genetic maps spanned 2486.97 and 3732.84 cM in length, for the F2 and RILs, respectively, with an average marker density of 1.37 and 1 marker/cM (Supplementary Table S3).
The chromosomes length ranged from 5.01 (chr. 3D) to 231.56 cM (chr. 1B) with an The first component discriminated PH from the other traits in all populations, except for the HD in F 2 . The second dimension discriminated HD and GH from PH, FTN, and TTN in all generations except for F 7 .

Construction of F 2 and F 6 Linkage Maps and Their Colinearity
After filtering, 3684 and 3958 SNP variants were retained for F 2 and F 6 , respectively, and used to build the two genetic maps. Thirty-nine and thirty-eight linkage groups were identified in F 2 and F 6 populations, respectively. The smaller linkage groups were joined to form the twenty-one chromosomes of T. aestivum L., based on the information available in the consensus map of Wen et al. [54]. The two genetic maps spanned 2486.97 and 3732.84 cM in length, for the F 2 and RILs, respectively, with an average marker density of 1.37 and 1 marker/cM (Supplementary Table S3).

ARF
TraesCS7A01G160800 GRF TraesCS7A01G163400 FRL1 TraesCS7A01G165600 a LOD, logarithm of the odds; b R 2 , phenotypic variance explained (%) for each QTL; c Add, estimation of the proportion of genetic variance due to additive effects; d Dom, estimation of the proportion of genetic variance due to dominance effects; e Transcript ID refers to the Chinese Spring genome, according to the IWGSC RefSeq v1.0 [57], (https://wheat-urgi.versailles.inra.fr/Seq-Repository/Annotation, accessed on 10 February 2021).

Candidate Genes and Their Metabolic Pathway
Several genes modulating photoperiodic flowering pathway were annotated within QTLs (Table 2 and Supplementary File S1). Two transcriptor factors (TF) involved in flowering control, the zinc-finger premature internode elongation 1 (PINE1, TraesCS4B01G342300) and the phytocrome-interacting bHLH factor (PIF4, TraesCS5A02G04960), were identified within QGh-4B and QTtn-5A, respectively. In addition, TFs that play a role in the regulation of flowering time and tillering, such as DELLA protein gai (GAI), flowering-promoting Figure 3. A graphical representation of the 2B and 4B linkage groups derived from the F 2 and F 6 maps and their comparison with the 90K consensus map. The distance between markers is in centimorgans (cM). The markers flanking the QTLs are colored in red (QHd-2B.1 and QPh-4B.1) and blue (QTtn-4B.1) in F 2 , while in green for the QTLs (QHd-2B.2, QPh-4B.2, and QTtn-4B.2) in F 6 . Solid lines connect the genetic loci in common between the F 2 and F 6 genetic maps and the 90K consensus map [60].
For HD, six QTLs were identified on chromosomes 2B, 2D, 4D, 5A, and 7D, with individual QTL explaining 3.70-30.98% of the phenotypic variance. Of these, two QTLs were detected on chromosome 5A (QHd-5A) and 7D (QHd-7D), both in F 6 and F 7 populations, with individual QTLs explaining the highest phenotypic variance in the F 7 (R 2 = 30.98 and 16.68% for QHd-5A and QHd-7D, respectively). Two distinct QTLs (QHd-2B.1 and QHd-2B.2) were found on chromosome 2B in the F 2 and RILs. However, the markers flanking the two QTLs were adjacent each other, with a distance of 2.24 cM according to the 90K consensus map [60] (Figure 3). One QTL (QHd-4D) was identified on chromosome 4D both in F 2 and F 3 , exhibiting 1.82% and 6.10% of the phenotypic variance, respectively. One QTL (QHd-2D) was identified on chromosome 2D in F 2 , with a negative additive effect.
Six QTLs were found for GH in both mapping populations (Table 2). Among all regions, a major QTL was found on chromosome 5A in late generations (F 6 and F 7 ), explaining the highest LOD scores and phenotypic variance (LOD = 7.23 and 10.21, R 2 = 13.1% and 22.6% in F 6 and F 7 , respectively). On chromosome 1B, two adjacent regions were identified explaining a high fraction of phenotypic variance (12.37% and 9.31% for F 6 and F 7 , respectively). Two QTLs were identified on chromosome 2D, of which, QGh-2D.1 was detected both in F 3 and F 6 , while the QGh-2D.2 was mapped only in F 6 . Another QTL (QGh-4B) was found on chromosome 4B in F 3 , with a dominance effect of −0.67.
Supplementary Figure S3 summarizes the metabolic pathways in which the candidate genes are involved. Six different pathways were identified for HD, GH, PH, TTN, i.e., photosynthesis, secondary metabolism, tetrapyrrole, carbohydrate, lipid, and N-metabolism. The functional categories more represented for all traits were classified as fermentations and light reactions for photosynthesis, OPP, and the nitrogen assimilation (Supplementary Figure S3).

Discussion
Traits such as plant height, heading date, juvenile growth habit, and tiller number, historically, have been subjected to strong selective natural and artificial pressure, to improve the adaptation of bread wheat to different climatic conditions and to increase the grain yield [66][67][68][69]. However, these same traits are not only important for increasing crop yield potential, but they are also useful in determining the adaptation to climate change [7].
In the present work, the genetic control of five morpho-physiological traits was investigated, using two mapping populations derived from the same parents, to identify associated QTLs and candidate genes.
To achieve our aims, plants were deep phenotyped during four growing seasons. High heritability was recorded for PH and HD in both mapping populations. This was consistent with previous reports in wheat and also in other plant cereal species such as rice and barley, indicating a high response to selection of these traits [70][71][72][73]. By contrast, low heritability was found for GH, TTN, and FTN, as also recently reported by Marone et al. [44] and Bilgrami et al. [74] for juvenile growth habit and tiller number, respectively. Continuous distribution or absence of discrete segregating classes for PH, TTN, and FTN suggested that its inheritance is either determined by a large number of genes with small effects or by a few major genes with substantial environmental effects. The presence of transgressive segregants in all traits investigated suggested that each of the two parental cultivars had desirable and undesirable alleles in various proportions for loci governing these traits, as shown in Table 2.
More than 3000 high quality SNP markers were used to build the two genetic maps, and as expected, most of them were placed on genome A and B, in line with previous results [75][76][77]. Wen et al. [54] showed that the D genome had fewer markers than the A and B genomes in the high-density consensus map in common wheat.
The resulting maps were different in length, being 2486.97 cM in the F 2 and 3732.84 cM in the F 6 . Similar to our results, Price and Tomos [78] obtained a map length of roughly 1280 cM in the F 2 population in rice, whereas a longer one (1680 cM) was observed in the F 6 produced by self-pollination of the same F 2 [79]. The observed differences were probably due to the different levels of heterozygosity, ranging about 50% and 5% in F 2 and F 6 RIL generations, respectively.
A total of 28 and 10 QTLs were found in the F 2 /F 3 and RIL populations, respectively, for all traits (Table 2). Among all QTLs identified in the present study, only one was detected across both generations. However, the comparative QTL analysis of chromosomes 2B and 4B between F 2 and F 6 populations showed that three QTLs for HD, PH, and TTN could be considered to be adjacent and nearly overlapping. The results of the joint analysis indicate F 2 /F 3 progenies have a higher power for detecting QTLs than the RIL populations. This is probably because the QTLs with dominant effects were only detected in the F 2 :F 3 but not in the RIL generations. For quantitative traits with low heritability, such as TTN and FTN, the QTLs were identified only in early generations. This contrasted with what was expected, given that for quantitative traits, the precision of QTL mapping with an F 2 population may be relatively poor and, to solve this problem, a corresponding F 3 progeny was phenotyped, as proposed by Zhang and Xu [80]. However, in our study a possible explanation could be linked to the different sowing method adopted; spaced plantings in early generations allowed a better and more accurate evaluation of the tillering capacity of the wheat lines [81]. Furthermore, a higher marker density was observed in F 2 than in F 6 , in the QTL regions. The low density of markers could affect the QTL analysis in the RIL populations [12,82]. Indeed, Yi et al. [83], in a previous study, which included QTLs consistency across generations derived from the same cross, suggested that results largely depend on many factors including generations, environment, and genetic backgrounds.
Pearson rank between the assessed traits revealed that HD was correlated with the different morpho-physiological traits such as GH, PH, and TTN, in agreement with Bilgrami et al. [74], Rabbi and Hisam [77], and Mecha et al. [84]. In our study, several QTLs for different traits co-localized on the same chromosome, suggesting that they were not distributed evenly in the wheat genome, but they tended to cluster in particular chromosome regions ( Table 2). According to their overlapping intervals, the regions on chromosome 5A, 2D, and 4D contained QTLs that were coincident for HD and GH, HD and PH, as well as HD and TTN, respectively. The two overlapping QTLs for TTN and FTN, found on chromosomes 1B and 7A, such as the common QTL for PH and TTN on chromosome 4B were expected. Other studies have confirmed that several QTLs/genes associated with canopyarchitecture-related traits shared regions associated with yield-grain-related traits [85][86][87][88].
According to the physical map information obtained in this study and those of the major genes known in the literature to be responsible for the main morpho-phenolic wheat traits (i.e., Vrn, Ppd, and Rht), several co-localizations were found.
In particular, two overlapped QTLs for HD (QHd-5A) and GH (QGh-5A) contained the well-known gene Vrn-A1 (TraesCS5A01G391700). Similarly, QGh-2D.2 was found close to the photoperiod gene (Ppd-D1). This means that adaptation genes such as Vrn and Ppd influenced the expression of juvenile growth habit trait either because of linkage, or because they exert a pleiotropic effect on this trait, as previously suggested by Marone et al. [44]. As expected, we observed a co-localization among the QTLs QHd-2B.1 and QHd-2B.2 and the photoperiod gene Ppd-B1 [89]. Similarly, the QHd-2D was only 4.5 Mb away from the Ppd-D1 (TraesCS2D02G079600) gene. In addition, QHd-4D was found to be 7.2 Mb from the Rht-D1 gene (TraesCS4D02G040400), suggesting that genes controlling plant height might also affect flowering time. In general, all the QTLs identified for HD had been extensively reported by the literature [38,61,[90][91][92], confirming the strength of the use of two mapping populations in our work.
Similar results were also shown for the PH trait, which is one of the most studied [65]. Two QTLs (QPh-4B.1 and QPh-4B.2) included the well-known dwarfing Rht-B1 gene (TraesCS4B02G043100) and the region overlapped also with the one recently identified for plant height by Wang et al. [64], and with another QTL associated with maturity date [90]. In addition, the Rht-D1 gene on chromosome 4D co-localized with the region identified by the QPh-4D in both early generations, sharing the common SNP marker RAC875_rep_c105718_304 with Zanke et al. [65].
Regarding TTN, all QTLs were identified in the early generations (with the exception of QTtn-4B.2 in F 6 ) and distributed on almost all chromosomes (Table 2), confirming previous studies conducted in different genetic backgrounds and agronomic contexts [74,93]. QTtn-4D was detected very close to Rht-D1 [65], similar to QHd-4D. The QFtn-4B detected on chromosome 4B was also interesting since it covered a genomic region close (6 Mb) to the Rht-B1 gene. In the past, Rht-B1b and Rht-D1b were associated with increased productive tillers, which also contributed to increased yield [94][95][96]. Indeed, in our study we found the QTL for TTN in the same regions where the QTL for grain yield (QTtn-1A.2) and NDVI at tillering (QTtn-7B) were identified [97,98]. In this study a significant association between Rht genes and tillering traits was confirmed. The QTtn-1A.1 identified by the SNP marker IWA6644 has been confirmed by a previous work by Bilgrami et al. [74].
It was observed that favorable alleles for each QTL were contributed by both parents, which was indicated by the additive effect of the QTLs. In addition, QTLs for TTN and FTN in the F 2 population showed both additive and dominant effects, suggesting the importance of both modes of gene action in conditioning tiller number and the greater difficulty in selection in late generations.
In addition to well-known genes involved in flowering time regulation (Vrn, Ppd, and Rht), QTL mapping detected candidate genes implicated in plant architectural traits including the regulation of the GA pathway, perception of light, photosynthesis, transport of photosynthates, and source-sink relationship. The role of these genes has been widely dissected in species such as rice and Arabidopsis, while they are still little investigated in wheat.
Gibberellin is a hormone that plays a central role in plant growth, organ development, and stress responses. In general, the existence of a regulatory network coordinating flowering and GA-dependent growth is well known. In our work we discovered genes coding enzymes of gibberellin biosynthesis (Ga3ox, GA20ox) and its receptors (DELLAinteracting protein GID1 and SD1) associated with GH and TTN [99,100]. A barley GA 20-oxidase gene (Hv20ox2) has been proposed as a candidate for sdw1/denso, likely ortholog to the rice sd1/Os20ox2 gene [101][102][103]. In barley, plants possessing the semi-dwarfing sdw1/denso gene have been characterized by prostrate GH, confirming our results [101]. In addition, further analysis revealed that the barley sdw/denso gene, has a pleiotropic effect on several agronomic traits such as plant height, heading, and flowering time [104].
Many zinc finger transcription factors have been found on chromosomes (1B, 2D, 4B, 4D, and 7A) such as CCCH-domain zinc finger proteins, which are known to be involved in plant growth habit and in leaf angle regulation [105]. Our findings confirmed the role of these transcription factors in determining the bread wheat GH, as reported in durum wheat [44].
Among its various functions, GA promotes the phenomenon of stem and leaf elongation [99]. In the flowering plants, altering internode length permits plants to vertically shift the distribution of their leaves and buds attached at nodes, changing their growth habit. Benefits of internode elongation at flowering include inflorescence escape from soilborne pests, outcompeting neighbors, and potential for diversified pollination and wider seed dispersal. In our work, the gene PINE1 was a strong candidate associated with GH. It was recently identified in rice as a transcription factor involved in coordinating internode elongation and photoperiodic signals. It regulates the reduction of the sensitivity of the stem to gibberellin, with a consequent promotion of flowering [58].
In addition, we found a TF belonging to the NAC family associated with all traits. The NAC genes coding for TF have been shown to influence plant height through regulating the key genes in the GA pathway, as well as flowering time [106]. It has been observed that the OsNAC2 promotes shoot branching; indeed, its overexpression resulted in more tillers [107]. Transgenic plants that constitutively expressed OsNAC2 had shorter internodes, shorter spikelets, and were more insensitive to gibberellic acid. In addition, the ERF genes that we identified in QTLs related to TTN, trigger the upregulation of the internode elongation, through an increase in bioactive GA levels [108].
Plants alter their morphology to avoid the risk shading, through changes associated with the light spectrum and intensity. The responses to shade avoidance include variation in the hypocotyl and stem and activation of multiple hormones such as auxin, ethylene, and GA. The transcription factor PIF4, which we found within QTLs for TTN, is activated by the photoreceptors signal via separate pathways in leaf responses in A. thaliana. PIF4 is a target for GA signaling through interaction with DELLA proteins, but also GA20ox and GA3ox expression, with a mechanism regulation currently unknown. It has been demonstrated that by destabilizing DELLAs, GA signaling promotes PIF4 function [99,109].
In addition, during the early growth phase and under abiotic stress conditions, PIF4 promotes auxin biosynthesis genes such as YUCCA, enhancing tissue auxin level that regulates small auxin up-RNA (SAUR) genes responsible for hypocotyl elongation when grown under light conditions [110][111][112].
Another important gene is the sucrose-signaling metabolite synthase trehalose 6phosphate synthase (TPS6) located in the interval of the Q.Ttn-3A, which showed the higher LOD (17.91). This gene encodes what is considered to be a potent signaling molecule in plants for embryogenesis and normal postembryonic growth and development. Wahl et al. [113] indicated that the TPS6 gene is required for timely initiation of flowering. In Arabidopsis, Fichtner et al. [114] suggested that the TPS6 gene has two different functions depending on some parts of the plant and some developmental stages, both as a signal and as a negative feedback regulator of sucrose levels.
Finally, we found genes involved in nitrogen metabolism and carbohydrate assimilation (RUBISCO) and metabolism (SUS) associated with all traits. In Arabidopsis, it has been demonstrated that there is a significant interaction among nitrate, carbohydrates, RUBISCO, and the starch accumulations that affect plant growth and development, stress responses, and signal molecules [111,115,116]. The selected genes underlying the QTL regions represent a partial selection of genes that could regulated the traits investigated, therefore, further studies are needed to validate the causative genes.

Conclusions
The results of the present study show that QTL mapping results obtained from populations derived from founding parents at different generations are reliable. The exploitation of the two mapping populations made it possible to identify all the main regions of the morphological characters and even identify new regions of interest especially for low heritability characters such as GH ant tiller number. The QTLs detected in early generation are more than those in late generation, especially for those minor QTLs. Further studies are needed to confirm the chromosomal regions responsible for low heritability traits such as TTN and FTN. For these traits, it will be necessary to utilize spaced plantings to properly assess tillering habits of wheat cultivars because plants grown at high population densities in grower production fields are not likely to represent full tillering potential of a line/genotype. It was particularly interesting to observe the strong interaction between the adaptation genes and the morpho-physiological characteristics of the plants. It is evident, in fact, that the genes, in the past, that were fundamental for increasing productivity, will also be so in the near future as they are also responsible for the characters indirectly associated with the new climate scenario. This requires an in-depth study of the effect of allelic variants of adaptation genes by exploiting modern genomic technologies.
Supplementary Materials: The following are available online at https://www.mdpi.com/article/10 .3390/genes12040604/s1, Table S1: Shapiro-Wilk test results. Normality test on the trait distributions among the generations using the Shapiro-Wilk test, Table S2: Phenotypic correlations among the five investigated traits. r = Pearson correlation coefficients. HD, heading date; PH, plant height; GH, growth habit; TTN, total tiller number; FTN, fertile tiller number, Table S3: Genetic map features. Main map characteristics concerning genetics maps constructed on F 2 and F 6 generations, Table S4: Genetic map data summary. Distribution of SNP markers across the twenty-one chromosomes in the F 2 and F 6 genetic maps. Chromosome length (cM), distance (cM), number of markers in common between F 2 and F 6 , and coefficient of collinearity (Rho) were reported, Figure S1: Phenotypic distribution for the five traits across the F 2 (grey bars), F 3 (orange bars), F 6 (red bars), and F 7 (green bars). The x-axis represents phenotypic values, whereas the y-axis corresponds to density. Phenotypic values of parental lines Lankaodali (L) and Rebelde (R) are also shown. HD, heading date; PH, plant height; GH, growth habit; TTN, total tiller number; FTN, fertile tiller number, Figure S2: Pearson's correlation matrix for the five traits measured across the F 2 (A), F 3 (B), F 6 (C), and F 7 (D) populations. The square blank refers to the Pearson's correlation test yielded a non-significant p-value (p > 0.05). HD, heading date; PH, plant height; GH, growth habit; TTN, total tiller number; FTN, fertile tiller number, Figure S3: MapMan analysis. The blue squares indicated the major pathways of metabolic processes in which candidate genes were involved. HD, heading date; PH, plant height; GH, growth habit; TTN, total tiller number; FTN, fertile tiller number, File S1: List of candidate genes within QTLs annotated using the Svevo durum wheat genomic sequence. HD, heading date; PH, plant height; GH, growth habit; TTN, total tiller number; FTN, fertile tiller number.
Author Contributions: P.D.V. and D.L. conceived the work; P.D.V. coordinated and supervised the activities; D.L., S.P. and F.S. provided seeds of the genotypes; F.F. and I.P. phenotyped the mapping populations; S.E. and P.V. performed statistical, QTL analyses, and data visualization; S.E. and F.T. identified candidate genes; P.D.V., F.T., P.V., S.E. and N.P. contributed to interpretation and presentation of data; P.V., P.D.V., F.T. and S.E. wrote the manuscript. All authors critically revised the manuscript. All authors have read and agreed to the published version of the manuscript.