Identification of Quantitative Trait Loci for Component Traits of Flowering Capacity Across Temperature in Petunia

For ornamental annual bedding plants, flowering performance is critical. Flowering performance includes the length of the flowering period, the longevity of individual flowers, and the number of flowers produced during the flowering period, or flowering capacity. Flowering capacity is a function of several component traits, including the number of branches producing flowers, the number of inflorescences per flowering branch, and the number of flower buds per inflorescence. We employed an F7 Petunia axillaris × P. exserta recombinant inbred line population to identify QTL for flowering capacity component traits. The population was phenotyped at 14, 17, and 20° over two years. Fifteen robust QTL (rQTL; QTL detected in two or more temperatures/years) were identified across six of the seven Petunia chromosomes (Chr) for total flower bud number (FlBud), branch number (Branch), flowering branch number (FlBranch), and primary shoot flower bud number (FlBudPS). The largest effect QTL explained up to 28.8, 34.9, 36, and 23.1% of the phenotypic variation for FlBub, FlBudPS, Branch, and FlBranch, respectively. rQTL for FlBud and FlBranch co-localized on Chr 1, and rQTL for FlBud, FlBudPS, and FlBranch co-localized on Chr 4. These regions in particular should be useful for identifying genes controlling flowering capacity of this important ornamental plant.

The number of flowers formed per inflorescence is also an important component of flowering capacity. The petunia inflorescence is a cyme producing from one, in the extra petals mutant, to many flowers (Souer et al. 1998). An F 2 Petunia integrifolia · P. axillaris population exhibited a bimodal distribution for flower number on the primary shoot (FlBudPS), and a QTL explaining 43% of the variation for this trait was identified in chromosome 6 ). An F 7 recombinant inbred line (RIL) population derived from that same F 2 population was phenotyped at 14, 17 and 20°. QTL for FlBudPS where identified in similar but non-overlapping regions of chromosome 6 explaining 32, 20, and 14% of observed variation at 14, 17, and 20°, respectively (Guo et al. 2017b).
The recent availability of Petunia spp. genome (Bombarely et al. 2016) and transcriptome (Guo et al. 2015) sequences greatly facilitates genetic mapping and gene discovery for traits of interest in the genus. In this study, an F 7 interspecific Petunia axillaris · P. exserta RIL population was utilized to characterize phenotypic variation and identify potential genetic interactions between total flower bud number and four flowering capacity component traits under a range of temperatures. We previously determined that P. axillaris produces more flower buds at first flowering than P. exserta (Warner 2010). This population was previously genotyped to develop single-nucleotide polymorphism (SNP) markers (Guo et al. 2017b). Identification of quantitative trait loci (QTL) for the component traits of flowering capacity can facilitate the development of marker-assisted breeding strategies to improve breeding efficiency for improved and novel cultivars and aid identification of candidate genes controlling these traits. The objective of this study was to identify QTL associated with flowering capacity component traits using the interspecific P. axillaris · P. exserta F 7 RIL population.

MATERIALS AND METHODS
Seeds of 171 F 7 P. axillaris (PI 667515) · P. exserta (OPGC943) RILs and the two parents were sown on 05 Nov 2014 and again on 20 Nov. 2015 in 288-cell plug trays filled with 50% vermiculite and 50% soil-less media (70% peat moss, 21% perlite, 9% vermiculite [v/v]; Suremix, Michigan Grower Products Inc., Galesburg, MI, USA). These RILs were previously genotyped using a genotyping-by-sequencing approach (Guo et al. 2017b). Seed trays were covered with clear dome lids and kept in a growth chamber at 23°and 50% relative humidity under a 10-h photoperiod (provided by fluorescent lamps) for germination. Dome lids were removed when 75% of the seeds had germinated within a tray. Seedlings were thinned to one plant per cell as needed. When seedlings had developed two true leaves, the air temperature was lowered to 20°.
Twenty-one days after seeds were sown, the trays were moved to the Plant Science Greenhouses at Michigan State University (East Lansing, MI) under ambient light. On 02 Dec. 2014 and 15 Dec. 2015, nine plants per RIL and parent were transplanted into 10-cm diameter round pots (height: 8.5 cm; 450 mL volume) with the soilless media mix described above and moved into treatments.
Three temperature treatments, each consisting of three replications of one plant each per RIL and parent, were arranged in a randomized complete block design within each temperature. Treatment air temperatures were constant 14, 17, or 20°under a 16-h photoperiod. Actual average weekly air temperatures are presented in Fig. S1. All plants received supplemental lighting (95 6 15 mmol m -2 s -1 of photosynthetically active radiation from 0600-2200 HR) provided by high-pressure sodium lamps. Initially, plants were grown pot-tight and were subsequently spaced to 20 cm between pot centers in each row and column 14, 21, and 27 days after initiation of treatments (DAT) at 20, 17 and 14°, respectively. Plants were overhead irrigated as needed with deionized water containing a water-soluble fertilizer (125 ppm N, 30 ppm P, 145 ppm K; MSU Orchid RO Water Special 13 N-3P-15K; GreenCare Fertilizers, Inc., Kankakee, IL).

Data collection
The number of nodes on the primary shoot were counted 0 and 14 DAT. Day 0 started on 06-08 Dec. 2014 and 19-21 Dec. 2015, depending on treatment. Development rate (DRate) was calculated as the increase in node number per unit time and expressed in nodes d -1 . The following data were determined for each plant when the first flower opened on the main stem: days to anthesis ( Internode length (Internode) was calculated as the average distance between nodes (cm).

Data analysis
Data were analyzed using Statistical Analysis Software v9.4 (SAS Institute, Cary, NC). Broad-sense heritability (H 2 ) was calculated for all evaluated traits as described by Fehr (1987)  where s 2 gy is the variance among the genotype by year, s 2 gt is the variance among the genotype by temperature, s 2 gty is the variance among genotype by temperature and year, s 2 e is the residual, y is the number of years in the study, t is the number of temperature treatments, and r is the number of replicates. Broad-sense heritability was calculated at individual air temperature treatments using the above equation, however the variance of the environmental effect was calculated as s 2 e ¼ s 2 gy y þ s 2 e ry and terms as described above.

Linkage map construction
Genotyping of the population was described by Guo et al. (2017b). Of the 171 RILs phenotyped, 158 had genotypic data available and were utilized for linkage map generation and QTL mapping. A total of 6,291 single nucleotide polymorphisms (SNPs) were converted into 368 bins based on recombination breakpoints (Xu 2013). A genetic linkage map was generated using JoinMap 4.0 (Van Ooijen 2006) with the bin markers. Bins with a similarity value of 1.00 were removed from the calculations. The bin markers were placed into individual linkage groups using the LOD (logarithm of the odds) thresholds from 2.0 to 10.0 and linkage groups were determined using LOD thresholds of 4.0 to 6.0. Marker order and map distance were calculated using the regression module with the Kosambi mapping function (Kosambi 1943). The mapping parameters were set to a recombination frequency of 0.30, a LOD score of 3.00, and a goodness-of-fit jump threshold of 5. The linkage groups were oriented and assigned chromosome (Chr) numbers according to a previous study (Bossolini et al. 2011).

QTL mapping
The 158 RILs and a total of 356 bin markers were used for QTL mapping. QTL analysis was performed using the composite interval mapping (CIM) Model 6 algorithm in QTL Cartographer v2.5 software (Wang et al. 2012). The forward-backward regression method was used with five markers to control for genetic background, as described in the software manual (Wang et al. 2012). The control parameters were set to a window size of 10.0 cM and the marker probability threshold was defined at 0.05. A walk speed of 1.0 cM and a LOD threshold of 3.6 (Lander and Kruglyak 1995) were used to identify significant QTL. LOD values for each QTL were calculated from the likelihood-ratio (LR) statistics. The proportion of total phenotypic variation explained n■ axillaris, PE = average for P. exserta, t-value = t-test comparing RIL means to parental line means, Transgression= t-tests comparing the highest RIL mean to the higher parental mean (Upper) and the lowest RIL mean to the lower parental mean (Lower) for each trait. x Ã and ÃÃ indicate significance at P , 0.05 and 0.001, respectively. w na = test was not performed because LLeng and LWid values for P. axillaris were lower than for any RIL at this temperature.
(%VE) by each QTL was estimated using R 2 values. QTL were visualized using MapChart v2.32 software (Voorrips 2002) using a subset of markers to facilitate visualization. Markers were filtered for visualization with the criteria that markers must be a minimum of 1 cM apart. QTL for the same trait with overlapping confidence intervals that were detected in two or more temperatures or years were considered the same QTL and were denoted as robust QTL (rQTL). The peak position of the rQTL potentially falls between the positions of the previous QTL. QTL names were determined by denoting "q" for QTL, followed by the trait abbreviation, the chromosome where the QTL was detected, and the order within the chromosome.

Data availability
Supplemental tables and figures, as well as phenotypic and genotypic data used to conduct statistical analyses, are available at GSA Figshare. Table S1 lists trait correlation coefficients within each temperature treatment. Figure S1 shows actual average weekly greenhouse temperatures over the experimental periods. Figures  S2-S4 show population distributions for measured traits at 14°, 17°, and 20°, respectively. Figure S5 shows the full P. axillaris · P. exserta F 7 RIL genetic linkage map.

RESULTS
The population exhibited transgressive segregation for all evaluated traits at least one temperature (  (Table 1). Development rate was negatively correlated with DTA and positively correlated with Nodes even though DTA was also positively correlated with Nodes at all air temperatures (Table 3). Development rate was positively correlated with Branch at 14°however it was negatively correlated with FlBud (Table S1).
Petunia exserta exhibited the earliest flowering time of the two parents at all air temperatures ( Table 1). Six of the RILs (AE11, AE20, AE230, AE301, AE315, and AE81) flowered earlier than either parent at all air temperatures in both years. However, there were 67 more lines that flowered earlier than either parent in 2015-16 (data not shown). Average DTA for the population was 68, 49, and 43 d at 14, 17, and 20°, respectively. DTA was positively correlated with FlBud and FlDiam at all air temperatures, however it was only positively correlated with Branch at 17 and 20° (Tables 3 and S1). Additionally, DTA was positively correlated with FlBudPS at all air temperatures but negatively correlated with FlBudLS at 17 and 20°. FlBud was positively correlated with FlBudPS and FlDiam at all air temperatures and positively correlated with Branch at 17 and 20°while negatively correlated at 14°(Tables 3 and S1). Mean FlBud was 36, 28, and 23 at 14, 17, and 20°respectively, which represents a 36% decrease in flower number from 14 to 20° (Table 1). Petunia axillaris had higher FlBud at all temperatures, higher Branch at 20°and higher FlBudPS at 17°Compared to P. exserta.

Broad-sense heritability estimates
Broad-sense heritability was relatively high for all measured traits (Table 4). Similar heritability estimates were observed across the different air temperatures for all traits excluding DRate, which was 46% and 44% lower at 14°Compared to 17 or 20°, respectively. With the exception of DRate, FlBudLS, FlBranch, and FlBudPS, all traits had high heritability (.0.7) across the air temperature treatments.

Linkage map
A total of 356 out of 368 bins were mapped to the seven Petunia Chrs (Fig. S5). The linkage map contained an average of 51 bins per Chr (Table 5). Chr 5 had the fewest markers with 23, while Chr 3 had the most with 92 bin s. The linkage map spanned a total genetic distance of n■

QTL analysis
Cumulatively, 15 QTL were detected for FlBud on Chrs 1-4 ( Table 6) and six of these were rQTL (Figure 1). The rQTL qFB1.1 was detected in five of the six environments across the two years. Two rQTL, qFB4.1 and qFB4.2, on Chr 4 were detected in four of the six environments and explained up to 27.2 and 28.8% of the phenotypic variation, respectively. The additive effects for the FlBud QTL ranged from 1.15 to 4.54. P. exserta contributed the beneficial alleles for two rQTL but P. axillaris contributed the beneficial allele for the remaining QTL, including four rQTL. For the FlBud component traits FlBudPS and FlBudLS, 14 and seven QTL were detected, respectively (Table 6). For FlBudPS, QTL were detected on all Chr except on Chr 5. The QTL for FlBudLS were detected on Chr 2, 3, 4, and 6. There were four rQTL for FlBudPS but no rQTL was detected for FlBudLS. The rQTL qFBP4.1 and qFBP4.4 for FlBudPS co-localized to the same regions on Chr 4 as the rQTL qFB4.1 and qFB4.2 for FlBud, respectively, (Table 6; Figure 1). Additionally, two rQTL for FlBudPS explained more than 25% of the phenotypic variation in at least one environment, whereas none of the QTL for FlBudLS explained more than 10%. The additive effects ranged from 0.15 to 0.48 and 0.47 to 0.78 for FlBudPS and FlBudLS, respectively. For FlBudPS P. exserta contributed the beneficial alleles for three QTL, however, P. axillaris contributed the beneficial alleles for the remaining QTL including the four rQTL. Additionally, P. axillaris contributed the beneficial alleles for five of the seven QTL for FlBudLS.
A total of 17 QTL each were detected for Branch and FlBranch, with 13 of these QTL co-localizing on Chr 1, 3, 4, 5, and 6 ( Table 6). There were two rQTL for Branch, and both were detected on Chr 3 (Table 6; Figure 1). The rQTL qBR3.3 explained from 7.5 to 36% of the phenotypic variation, depending on temperature and year. This rQTL also had the greatest additive effect on Branch. Petunia axillaris contributed the beneficial alleles for six Branch QTL while P. exserta contributed the beneficial allele for the remaining QTL, including the two rQTL. Three rQTL were detected for FlBranch on Chr 1, 4, and 5, respectively (Table 6; Figure 1). None of the QTL for FlBranch explained more than 25% of the phenotypic variation, but six explained 10-20%. Additionally, the QTL for FlBranch have additive effects ranging from 0.29 to 0.77 and the beneficial alleles were equally contributed by P. axillaris and P. exserta.
A total of 15 rQTL were detected for four traits on Chr 1-6 (Table 6; Figure 1). Six of these rQTL were detected on Chr 4, with only one rQTL detected each on Chr 5 and 6. Three rQTL detected for FlBud, FlBudPS, and FlBranch co-localized to a 5 cM region, whereas the three rQTL detected for FlBud, and FlBudPS co-localized to a region of approximately 3 cM on Chr 4. Additionally, two rQTL detected for FlBud and FlBranch co-localized to a 1 cM region on Chr 1.

DISCUSSION
Flower number is an important trait that influences the aesthetic value of ornamental plants. Desirable flower characteristics include increased flower number and repeat or continuous blooming. However, quantitative analysis and candidate gene identification for these traits have not been comprehensively studied in ornamental crops. In this study, QTL for flowering capacity component traits of an F 7 P. axillaris · P. exserta RIL population were identified following phenotypic evaluation across multiple temperature environments and years. The QTL results are presented on a genetic linkage map, although the P. axillaris genome sequence is available (Bombarely et al. 2016). This is due to the level of fragmentation of the P. axillaris genome and the employment of bin markers for mapping, in which bins are often comprised of multiple SNPs that map to more than one genomic scaffold. The physical location of every SNP in each bin marker utilized for this study was previously reported (Guo et al. 2017b). The total linkage map distance reported here (270.1 cM) is shorter than would be expected. However, previous linkage maps generated for Petunia have often resulted in short linkage groups due to a low frequency of recombination (Strommer et al. 2002;Galliot et al. 2006;Vallejo et al. 2015;Guo et al. 2017b).
Fifteen QTL were detected for FlBud with both parents contributing favorable alleles (Table 6). The flowering capacity of a plant is a product of multiple traits that influence total flower number, including the number of branches, the number of inflorescences per branch, and the number of flowers per inflorescence. Dissecting the genetic control of these traits is challenging because a single genotype may exhibit a wide range of phenotypic variation in differing environments. The complex interaction between genotype and environment is compounded because multiple genes could be in linkage within the genetic region associated with the trait (Darvasi and Pisante-Shalom 2002).
n■ n■ In a two-year field evaluation of this same P. axillaris · P. exserta RIL population (Cao et al. 2018), and a P. integrifolia · P. axillaris RIL population (Cao et al. 2019), QTL for flower count (number of open flowers per plant, collected weekly for seven weeks) were identified on Chr 1, 2 and 4 in the AE population and 1 and 2 in the IA population. QTL for FlBud identified in the current study on Chr 1, 2 and 4 did not co-localize with the flower count QTL. However, for both traits, positive additivity was provided by P. exserta for the QTL on Chr 1 and by P. axillaris on Chr 4, while for the two FlBud QTL on Chr 2 each parent contributed positive additivity at one locus (Table 6).
x LOD values calculated from likelihood-ratio statistics.
w LOD threshold determined at 0.05 probability based on 1,000 permutations. v Additive effect of QTL, positive values indicate beneficial alleles from P. axillaris u Percentage of variation explained by QTL estimated using R 2 statistics. flower bud number-related QTL, P. axillaris contributed more beneficial alleles, whereas P. exserta contributed more beneficial alleles for branching-related QTL. However, both parents contributed favorable alleles for all traits. These results support the utility of incorporating wild species into breeding programs to introgress alleles that may have been lost during breeding to improve flower component traits, although linkage drag on other important traits is of concern.
Four branching QTL, including one rQTL each for Branch and FlBranch, were detected on the same chromosome as the previously identified QTL for branch number in a P. integrifolia · P. axillaris F 2 population . Additionally, four QTL including one rQTL for FlBudPS was detected on the same chromosome as the previously identified QTL for flower buds on the main stem in the F 2 population. The QTL for total number of flower bud on the primary stem on Chr 6 (FBP6.1) and total branch number on Chr 1 (BR1.1) explained 43 and 26% of the variation, respectively, in the F 2 population . Conversely, in this study, the QTL for these traits that were detected on the same Chr explained only 6-11% of the variation (Table 6). Additionally, in the RIL population, two major QTL for each trait FlBud and FlBudPS were detected on Chr 4 and one major QTL on Chr 3 and one on Chr 4 for Branch. Whereas the QTL (FB1.1) for flower number on Chr 1 in the P. integrifolia · P. axillaris population was not a major QTL and only explained 13% of the variation, but it had a large additive effect (17.78 flowers) from P. axillaris . While the largest additive effect for any total flower number QTL in the P. axillaris · P. exserta population was also inherited from P. axillaris, the effect was much lower at 4.54.
Additionally, in the RIL population, P. exserta contributed the beneficial alleles for the QTL on Chr 1 for FlBud, which indicates that both parents can provide beneficial alleles for this trait.
Across all temperatures and within each temperature, FlBud was consistently most highly positively correlated with FlBranch and FlBudPS (Table 1). Additionally, FlBranch and FlBudPS were highly positively correlated at each temperature, and a large effect rQTL for FlBranch (qFBN4.1) co-localized with a large effect rQTL for FlBudPS (qFBP4.1), suggesting potential for a common mechanism regulating vegetative and inflorescence branching. Some genes impacting both branch number and flower number per inflorescence have been identified. For example, the tomato BLIND gene encodes a MYB transcription factor that controls lateral meristem initiation, with blind mutants exhibiting reduced numbers of lateral shoots and flowers per inflorescence (Schmitz et al. 2002).
Several plant hormones have been implicated in regulating branching, including auxins, cytokinins, and strigolactones (Shimizu-Sato et al. 2009;Drummond et al. 2009). Auxins maintain shoot apical dominance and repress axillary outgrowth by downregulating cytokinin biosynthesis (Eklof et al. 1997;Nordstrom et al. 2004). In contrast, cytokinins promote axillary bud outgrowth even in the presence of auxin at certain developmental stages (Müeller and Leyser 2011). One rQTL for FlBud co-localized to the same region as the rQTL for FlBranch and FlBudPS on Chr 4 (Figure 1). In rice, a QTL for spikelets per panicle and primary branch number co-localized (Balkunde et al. 2013). One of the four candidate genes within the QTL region was a putative expressed nitrilase, which converts indole-3-acetonitrile (IAN) to the auxin indole-3-acetic acid (IAA) through hydrolysis (Kobayashi et al. 1993).
Strigolactones are carotenoid-derived plant hormones that have been identified as inhibitors of axillary bud outgrowth and shoot formation Gomez-Roldan et al. 2008;Kretzschmar et al. 2012). In chrysanthemum, phenotypic variation for shoot branching was associated with allelic variation in genes in the strigolactone pathway (Klie et al. 2016). The nearest bin marker to the branching QTL qBR4.2 and qFBN4.3 (AE_bin_210_117), which explained ca. 20% of the variation for Branch and FlBranch, contains five SNPs located on the P. axillaris genome scaffold containing PhDAD2 (Peaxi162Scf000081; 2.9 Mb; (Bombarely et al. 2016)), which encodes an a/b hydrolase involved in strigolactone perception. Orthologs of the branching-and strigolactone pathway-related genes MORE AXILLARY BRANCHING (MAX), CAROTENOID CLEAVAGE DIOXYGENASE (CCD), and TCP have been identified in petunia, including PhMAX2B, PhCCD7, PhCCD8 and PhTCP1-3 (Drummond et al. 2009;Drummond et al. 2012;Drummond et al. 2015). No branching QTL identified in this study localized to PhCCD7 or PhCCD8. However, PhMAX2B is located on scaffold Peaxi162Scf00384 of the P. axillaris genome. This 1.35 Mb scaffold contains a marker (AE_bin_301_62_14_156_2_2) flanking the rQTL qFBN.5.1 for flower branch number. Additionally, PhTCP1 is located on genome scaffold Peaxi162Scf00086, which contains a marker (AE_bin_89_94_49_1), located ca. 300 kb from PhTCP1, flanking QTL for both branch number (qBR1.2) and flowering branch number (qFBN1.2). Understanding the potential role of these genes, and identifying additional genes of interest in these QTL regions, will help develop a more thorough understanding of the quantitative mechanism for branching regulation in petunia, and the contribution of branching to flowering capacity.
The current study of flower production and its component traits provides new insight into its complex genetic control. Co-localization of rQTL for flower number and flowering capacity component traits on Chr 1 and 4 provide attractive targets for future studies to fine map these candidate regions to identify genes controlling flower capacity component traits and molecular markers for improving flower production in petunia through marker-assisted breeding.