Major Contribution of Flowering Time and Vegetative Growth to Plant Production in Common Bean As Deduced from a Comparative Genetic Mapping

Determinacy growth habit and accelerated flowering traits were selected during or after domestication in common bean. Both processes affect several presumed adaptive traits such as the rate of plant production. There is a close association between flowering initiation and vegetative growth; however, interactions among these two crucial developmental processes and their genetic bases remain unexplored. In this study, with the aim to establish the genetic relationships between these complex processes, a multi-environment quantitative trait locus (QTL) mapping approach was performed in two recombinant inbred line populations derived from inter-gene pool crosses between determinate and indeterminate genotypes. Additive and epistatic QTLs were found to regulate flowering time, vegetative growth, and rate of plant production. Moreover, the pleiotropic patterns of the identified QTLs evidenced that regions controlling time to flowering traits, directly or indirectly, are also involved in the regulation of plant production traits. Further QTL analysis highlighted one QTL, on the lower arm of the linkage group Pv01, harboring the Phvul.001G189200 gene, homologous to the Arabidopsis thaliana TERMINAL FLOWER1 (TFL1) gene, which explained up to 32% of phenotypic variation for time to flowering, 66% for vegetative growth, and 19% for rate of plant production. This finding was consistent with previous results, which have also suggested Phvul.001G189200 (PvTFL1y) as a candidate gene for determinacy locus. The information here reported can also be applied in breeding programs seeking to optimize key agronomic traits, such as time to flowering, plant height and an improved reproductive biomass, pods, and seed size, as well as yield.


INTRODUCTION
Common bean (Phaseolus vulgaris L.) is the most important food legume for direct human consumption. It is considered to be a rich source of proteins, micronutrients, and calories for human daily needs (Broughton et al., 2003). It is grown over a wide range of latitudes; however, there is an adaptation of each cultivar to a relatively narrow range of latitudes. Wild bean forms present indeterminate growth habits and require day-lengths of <12 h to initiate flowering. Together, earlier flowering and determinate growth habit genotypes under day-lengths longer than 12 h allowed for the adaptation to higher latitudes (Gepts and Debouck, 1991). Domesticated and wild common bean display notable differences in growth habit types. Thus, there is a considerable variability in domesticated cultivars, which show pronounced differences in growth form, i.e., determinate vs. indeterminate, and growth habit, i.e., I, II, III, and IV (Debouck and Hidalgo, 1986;Debouck, 1991). Determinate common bean cultivars generally flower and mature early, and the transition of the terminal shoot meristem from vegetative to reproductive state results in a terminal inflorescence in the axil of the older leaf primordia. By contrast, in indeterminate cultivars, the terminal shoot meristem continuously produces modular units until senescence, each one consisting of a leaf and an inflorescence. Thus, the plant will have a terminal shoot meristem that remains in a vegetative state throughout the production of vegetative and reproductive structures (Ojehomon and Morgan, 1969;Tanaka and Fujita, 1979). It has been documented that stem termination has great effects on plant height, flowering and maturity period, amount of branching, length of internodes on the main stem, and node production, which conditions how many flowers and leaves, and therefore pods and seeds, are produced. Thus, understanding the genetic control of vegetative growth and flowering time in common bean will enable genetic manipulation of major components of yield.
Previous studies demonstrated that the FIN locus mainly regulates stem growth habit in common bean and that the indeterminate growth habit is dominant over the determinate one. This gene was mapped on the Linkage Group (LG) Pv01 (Norton, 1915;Koinange et al., 1996). Despite the fact that FIN is a monogenic locus, it is possible to find a wide range of stem termination types among common bean cultivars, which may be regulated by a second unnamed locus mapped on Pv07 (Kolkman and Kelly, 2003). In addition, control of twining has been attributed to the TOR gene, distinct from the FIN locus (Norton, 1915) although either FIN has a pleiotropic effect on twining or that TOR is tightly linked to FIN. Additionally, other loci have been reported to control flowering time and other flowering-related traits in common bean (Norton, 1915;Wallace et al., 1993;Jung et al., 1996;Bassett, 1997;McClean et al., 2002;Tar'an et al., 2002;Kolkman and Kelly, 2003;Blair et al., 2006;Checa and Blair, 2008;Chavarro and Blair, 2010;Pérez-Vega et al., 2010). The majority of these loci have been detected as Quantitative Trait Loci (QTL), although they have been treated as Mendelian factors in some studies. Probably, different works have detected the same loci; however, the lack of common markers among different mapping studies makes difficult to determine whether or not they are the same locus.
During the past two decades, the model species Arabidopsis thaliana has been mainly used to study the process of phase transition from vegetative to reproductive growth at developmental, environmental, genetic, and molecular levels (Weigel, 1995;Yanofsky, 1995;Bradley et al., 1997;Ma, 1998;Pidkowich et al., 1999). In this species, floral meristem identity is determined by two different pathways. The heterodimer FLOWERING LOCUS T (FT)/FLOWERING LOCUS D (FD) is proposed to promote flowering initiation through activating the APETALA1 (AP1) expression (Abe et al., 2005;Wigge et al., 2005). Furthermore, another key regulator of flowering initiation is TERMINAL FLOWER1 (TFL1), a floral repressor and a regulator of inflorescence meristem development, which acts by repressing the expression of AP1 and LEAFY (LFY) floral identity genes (Bradley et al., 1997;Ohshima et al., 1997;Nilsson et al., 1998;Boss et al., 2004). FT and TFL1 have closely related sequences, although key amino acids have been found, which have been proposed as responsible for making that these two proteins perform opposite functions. Within the legumes, the value of a comparative approach to candidate gene identification has led to the characterization of the molecular identity of TFL1 co-orthologs such as FIN in common bean (designated as PvTFL1y; Repinski et al., 2012); Dt1 in soybean (GmTFL1; Tian et al., 2010); and DETERMINATE (DET; PsTFL1a) and LATE FLOWERING (LF; PsTFL1c) in pea (Foucher et al., 2003;Weller and Ortega, 2015). Recent findings have clearly shown that in several legume species, determinate inflorescence architecture is conferred by mutation of specific TFL1 genes (Benlloch et al., 2015). The determinate growth habit caused by mutations within specific TFL1-homologs in other grain legumes indicates that the determinate function is conserved in these species (Kong et al., 2010;Tian et al., 2010;Kwak et al., 2012;Repinski et al., 2012;Dhanasekar and Reddy, 2014).
In this study, the aim was to identify the genetic determinants of vegetative growth and its relationship to days to flowering and fruit maturation in two Recombinant Inbred Line (RIL) populations by using a mixed-model based composite mapping method for QTL detection. Both populations were derived from inter-gene (Andean × Mesoamerican) pool crosses between determinate and indeterminate genotypes, and shared an Andean determinate common parent, which allowed for the comparison of the results in reference to a tester line. Parents of each RIL population also differ in rate of plant production (leaf, flower and fruit size, and yield components), allowing for the dissection of the genetic architecture of these traits, as well as the study of the relationship among vegetative growth and these traits. Comparative QTL mapping indicated that, in both RIL mapping populations, the genomic region where the FIN locus is located not only is involved in the regulation of vegetative growth and rate of plant production but also affect flowering time, suggesting a pattern of pleiotropic effects accounting for the genetic bases of these traits.

Plant Material and Trials
Two F 2:8 RIL populations derived by single-seed descent from an F 2 population from crosses between the Mesoamerican (M) and Andean (A) gene pools were used. The MA population was generated from the cross between lines PHA0419 (Mesaomerican) and Beluga (Andean), whereas the AM population was obtained from the Beluga and PHA0399 (Mesoamerican) cross. Beluga is a large white kidney seed which seed weight averages 60 g 100 seeds −1 (range: 52-66 g 100 seeds −1 ), 55 cm in length of main stem with a type I determinate growth habit, and blooms 30 days after planting (DAP). Beluga is resistant to bean common mosaic virus (BCMV), as it bears the autosomal dominant hypersensitive I gene, and possesses the Andean Co-1 gene for resistance to races 65 and 73 of anthracnose (Kelly et al., 1999). Both PHA0419 and PHA0399 are great northern beans which averages 86 g 100 seeds −1 (range: 76-103 g 100 seeds −1 ), and shared a type IV climbing indeterminate growth habit (averages 278 cm of length of main stem) and late flowering (46 DAP as average). PHA-0419 possesses the Mesoamerican Co-4 2 , Co-6, and Co-10 genes that condition resistance to races 23, 39, 102, 448, and 1545 of anthracnose, while PHA-0399 carries Co-4 2 , Co-4 3 , and Co-6 genes that condition resistance to races 23, 55, 102, and 1545 of anthracnose (Santalla et al., 2004;Figure 1

Field Evaluation and Data Collection
Data were collected for time to flowering and fruit maturity evaluated as days to flowering, days to immature pod harvest, and days to physiological pod maturity (Supplementary Table 1). The days to flowering (FT) trait was scored as the number of days from planting date to the opening of the first flower. The days to young-green pod (PGT) trait was recorded when a plant had 50% of the immature pods. The days to physiological pod maturity (PST) trait was scored as the number of days from planting date to the appearance of a mature dry pod on a primary branch.
Measurements were collected for growth ability reordered as length of main stem, number of primary stem branches, and internode length, in five random plants from each RIL. Length of main stem (LMS) was measured in cm from the base of the plant to the uppermost leaflet on the longest branch. The number of primary stem branches (NPB) was recorded as the number of the stems winding around the support strings at the midpoint of the plant length. Internode length (LI) was recorded at the fifth internode on the main stem in cm.
To quantify crop growth and productivity, flower, pod and seed size, number of seeds per pod, number of pods per plant, and seed yield were determined in ten random flowers and fruits on a plot basis. Bracteole length (BL) measured in mm along the midrib of the lamina. Bracteole width (BWI) recorded in mm between the widest lobes of the lamina perpendicular to the lamina mid-rib. Maximum leaflet width (LWI) was measured in cm at the largest point perpendicular to the midrib. Leaflet length (LL) was recorded in mm from the lamina tip to the point of petiole intersection along the midrib. Pod length (PL) was recorded in mm as the exterior distance from the peduncle connection point to the apex excluding the beak. Pod width (PWI) determined in mm as the distance at right angles to the sutures at the level of the second seed from the apex. Pod thickness (PT) was recorded in mm as the distance between sides at the level of the second and third seed from the apex. Seed width (SWI) was the longest distance across the hilum, in mm. Seed thickness (ST) was measured in mm from hilum to opposite side. Seed length (SL) was measured in mm parallel to the hilum. Seed weight (SW) was determined on 100 dry seeds per plot. The number of seeds per pod (NSP), number of pods per plant (NPP), and seed yield (SY) were recorded and expressed in kilograms per hectare at moisture content of 140 g kg −1 .

Experimental Design and Statistical Data Analysis
Comprehensive statistical analysis (mean value, standard error, and range of variation) and normality test (Kolmogorov-Smirnov test) were carried out for each quantitative trait and environment. LMS, NPB, LI, BL, BWI, SW, NSP, NPP, and SY traits failed to meet normality assumptions and the Box-Cox transformation was used to improve normality, while LOG transformation was used for LL and LWI traits. Significant variation in the expression of traits through the environment conditions was analyzed using PROC MIXED (SAS Institute Inc. 9.04, Cary, NC, USA). The estimates of variance components were obtained by the REML method with Proc MIXED in SAS9.04 and used to calculate the broad-sense heritability on a progeny-mean basis (h 2 = σ 2 λ /[(σ 2 t /e) + σ 2 λ + (σ 2 e /re)] where: σ 2 λ = genetic variance of the trait; σ 2 t = variance due to environmental factors; σ 2 e = error variance; r = number of replications and e = number of environments). In order to increase the precision of the entry mean basis heritability estimate, it was used the harmonic means for the number of replications and environments, where each experimental line was tested (Holland et al., 2003). Approximate standard errors of heritability estimates were obtained by the delta method (Holland, 2006). Phenotypic Pearson correlation coefficients for the different traits were calculated using PROC CORR through the environment conditions in SAS9.04.

DNA Isolation and Marker Analysis
Young leaves from individual plants of both RIL mapping populations (60 and 179 lines for AM and MA, respectively) were collected for genomic DNA isolation as described by Chen and Ronald (1999). Total DNA was stored in sterile water and examined by electrophoresis in 1% agarose gels in 1X SB buffer (10 mM sodium boric acid). The concentration and quality of extracted DNA was determined by reading at 230, 260, and 280 nm using Nanodrop Thermo ScientificTM (NanoDrop 2000). DNA was diluted in sterile water to get a working dilution of 5-10 ng/µL, which was used for PCR analysis.

Linkage Map and QTL Analysis
JoinMap 4.0 software (van Ooijen, 2006) was used to construct the genetic linkage maps for both MA and AM mapping populations. A minimum logarithm of odds ratio (LOD) of 6.0 was considered to establish significant linkage. Locus order Frontiers in Plant Science | www.frontiersin.org within the LOD grouping was generated for each LG using the regression mapping algorithm with the following JoinMap parameters: Rec = 0.3, LOD = 2.0, and Jump = 5. The Kosambi map function (Kosambi, 1944) was used to calculate the genetic distance between markers.
LGs were designated according to Pedrosa-Harand et al. (2008). JoinMap 4.0 (van Ooijen, 2006) was also used to generate pairwise recombination frequencies and LOD scores for the selected sets of representative loci for each LG, which were then combined into a single group node in the navigation tree. The regression mapping algorithm was used and the LG lengths for the consensus map of all the representative markers were calculated.
The physical position of genetic markers was obtained by sequence similarity analysis using BLASTN (Altschul et al., 1997) against the common bean genome (Phytozome v.11: Pv1.0; Schmutz et al., 2014) in the Phytozome database (http://www.phytozome.net/). The correlations between physical distance and genetic map in each LG were calculated by Spearman's rank correlation coefficients.
QTLNetwork 2.0 software (Yang et al., 2008) was used for multi-environment QTL analyses. In order to identify putative single-locus QTLs and their environment interactions (QTLs × Environment, QE), a mixed-model based composite interval mapping method (MCIM) was carried out for one-dimensional genome scan. In addition, with the aim to detect epistatic QTLs (E-QTL) and their environment interaction effects (E-QTLs × Environment, E-QE), a two-dimensional genome scan was performed. A QTL was declared significant when the F-value was higher than the F-value threshold determined by a 1000 permutation test at 95% confidence level. Markov Chain Monte Carlo method was used to estimate the effects of QTLs and environment interactions (Wang et al., 1994). Both testing and filtration window sizes were set at 10 cM, with a walk speed of 1 cM. Candidate interval selection, putative QTL detection, and QTL effect was estimated with an experimental-wise significance level of 0.05. MapChart 2.2 software (Voorrips, 2002) was used to draw the genetic map and the detected QTLs. QTL regions were positioned onto the consensus map. QTL designations were made using abbreviations for the quantitative trait, and followed by LG number at which the QTL was mapped.

Vegetative Growth and Time to Flowering Variation
Genetic variation for vegetative growth as length of main stem and rate of plant production has been studied, together with its relationship to days to flowering and maturity in both inter-gene pool RIL populations. The populations were grown in different years and locations (field vs. greenhouse). Both populations segregated for different levels of growth ability: indeterminate vs. determinate growth habits. Classification of the 60 and 179 lines of AM and MA RIL populations for growth habit identified 37 and 63 lines as homozygous determinate type I, and 23 and 116 lines as homozygous indeterminate type IV, respectively. The observed growth habit segregation fitted to a 1:1 ratio (χ 2 = 3.27, P ≤ 0.05) for the AM population, indicating that a single gene determined the trait. However, growth habit distribution appeared distorted in the MA population (χ 2 = 15.69, P ≥ 0.05). On the basis of segregation analysis results, the gene for growth habit (FIN) was mapped along with the molecular markers. Supplementary Figure 1 shows phenotypic distribution of the RIL populations based on line means. The large range of variation and transgressive segregations observed for most traits in both RIL populations suggested a complex control of these traits, with positive alleles shared between the two parents of the RIL populations. Transgressive segregation in both directions was observed for days to flowering and maturity in both populations. For LMS, bimodality was observed in the AM population, with a clear separation of phenotypic classes that would indicate monogenic inheritance, although this separation was not evident in the MA distribution. Lines shorter or taller than the height of the parents were found for LMS trait in both populations. Likewise, for NPB, the number of branches produced by many of the RILs was higher than the parental lines in both populations, which indicated a positive transgressive segregation for this trait, although some skewing was observed in MA to low values. For LI, the histograms showed a similar pattern across both populations, indicating positive transgressive segregation. Hence, in general, the phenotypic segregations for rate of production traits in these two RIL populations exhibited normal distribution, and transgressive segregation in both directions, a typical phenomenon of a quantitative trait, regulated by several genes and influenced by the environment.
Mean values, standard errors, and ranges of variation for the quantitative traits in each RIL mapping population for each environmental condition have been summarized in Supplementary Tables 2, 3. Mesoamerican PHA0419 and PHA0399 lines were late in flowering and taller, with larger rates of plant production compared to Beluga. In both RIL populations and for all evaluated traits, it was found differences among environments in mean values and ranges of variation, although environment × line interactions were not significant for most of the evaluated traits in both RIL populations. Significant differences among parents were detected except in some environments for PST, LL, LWI, PL, and PWI traits in the MA RIL population (Supplementary Table 2) and for BL, PWI, and ST traits in the AM RIL population (Supplementary Table  3). Despite of that, it was observed significant differences among RILs for all quantitative traits except for LL in MA population under F108 and F109 environmental conditions (Supplementary Table 2).
High broad-sense heritability estimates (≥ 0.50; Supplementary Table 4) were detected for most of the quantitative traits across the four environments except for NPB, LI and BWI in both populations, and PST, BL, and PT in MA population. Higher heritability estimates and correlations were observed in AM population than in MA population. LMS and FT were significantly correlated in both AM (r = 0.50, P ≤ 0.001) and MA (r = 0.38, P ≤ 0.001) populations. Their inter-relationship suggests that some genomic regions influence both traits. Finally, a negative and significant correlation was observed between FT and NPB (r = -0.43, P ≤ 0.001) and positive and significant correlation between FT and LI (r = 0.33, P ≤ 0.001) in AM population, while no significant correlation values were observed in MA population. These correlations indicated that, on average, genotypes with a low number of primary branches, and high length of internodes showed a later flowering date and a longer main stem. As expected from the ontogenetic pattern of vegetative growth, LMS-values were positive and significant correlated with plant parts size (bracteole, leaf, pod, and seed) in both populations except for BWI in AM RIL, suggesting a pleiotropic effect. Positive and significant correlations were observed between FT, PGT, and PST; these last traits were positively correlated with LI in both RIL populations. No significant or low correlation values were observed between FT and all other rate of plant production traits measured. Significant and positive correlations were found between PST and rate of plant production traits. In contrast, the correlations with PGT were negative. SY presented a significant and low positive correlation value with FT in the MA and AM RILs (r = 0.16, r = 0.15, P ≤ 0.01). The correlations between SY and PGT, and PST were not significant except for PST in the MA RIL (r = 0.37, P ≤ 0.001).

Marker Segregation Analysis and Consensus Genetic Map Construction
Six hundred and thirty-four markers were screened for DNA polymorphisms in the parents of MA RIL, which rendered a 36% polymorphism rate. Sixty-two (27%) out of the 228 polymorphic markers evaluated in the MA RIL population exhibited segregation distorsion and thus, they could not be mapped. Finally, the genetic map of the MA RIL population (Figure not shown) consisted of 166 loci (158 SSRs, 1 SCAR, 6 SNPs and the FIN locus) distributed in 11 LGs. Out of 166 markers, 56 and 110 were dominant and codominant, respectively.
LGs were named as reported in Pedrosa-Harand et al. (2008) using for the assignment of LG number and orientation 55 common SSR markers previously mapped (Freyre et al., 1998;Blair et al., 2008Blair et al., , 2010Blair et al., , 2011Cichy et al., 2009;García et al., 2011). The map covered a genetic distance of 1188.9 cM, with an average of 7.2 cM among markers, which ranged from 5.3 cM (Pv02) to 16.0 cM (Pv05). The average genetic distance per LG was 108.1 cM, ranging from 70.6 cM (Pv10) to 134.1 cM (Pv07). A detailed description of this map is provided in Supplementary Table 5.
Likewise, 364 markers were screened for DNA polymorphisms in the parents of AM RIL, which rendered a 39% polymorphism rate. A total of 245 polymorphic markers were evaluated in the AM RIL population. Sixty-five (26%) out of these 245 markers, could not be mapped as they were not linked to other markers. Thus, the AM genetic map was constructed with a total of 180 loci (179 SSRs and FIN locus), from which 73 were dominant and 107 codominant. These loci were distributed among 11 LGs (Figure not shown) that covered a genetic distance of 1175.5 cM. The density of markers ranged from 3.4 cM (Pv01) to 12.6 cM (Pv04), with an average of 6.9 cM per LG. The longest LG was Pv08 (144.4 cM), whereas Pv11 was the shortest LG (48.2 cM), with an average genetic distance of 106.9 cM per LG. A complete description of the AM map is shown in Supplementary Table 6.
The construction of a consensus map (Figure 2) was performed by connecting the MA and AM mapping data. To integrate both MA and AM maps in a single consensus map, 103 common markers were used as anchor points. As a result, marker segregation data were assembled for a total of 202 marker loci (196 SSRs, 5 SNPs, and the FIN locus) into 11 LGs. The total length of the consensus genetic map was 1156.2 cM and had a marker average density of one marker per 6.1 cM ( Table 1). The marker order of the integrated map was largely collinear with the two individual maps, although a few local inversions and marker rearrangements over short intervals were observed. Most of the markers from both RIL populations showed a good linear relationship between their position on the genetic map and on the physical map of the common bean genome (Supplementary Table 7).

Multiple Environment Single-Locus QTL Analysis
QTL analysis based on MCIM mapping using QTLNetwork 2.0 was undertaken to identify single-locus QTLs across all environments. The positions of QTLs and their confidence intervals for the traits on the consensus map are shown in Figure 2 and Tables 2, 3. Thus, multi-environment QTL analyses allowed for the detection of 43 and 40 single-locus QTLs with significant additive main effects and/or QE effects for the evaluated quantitative traits in MA and AM RILs, respectively.
The distribution of the flowering and maturity time QTLs in the MA RIL population ( Table 2) varied from 3 on Pv01 (FT-1 MA , PGT-1 MA , PST-1 MA ), to 2 on each Pv02 (PGT-2 MA , PST-2 MA ), and Pv09 (PGT-9 MA , PST-9 MA ). All these QTLs had a positive additive value, indicating that alleles from the PHA0419 parent increase flowering and maturity times, whose main additive effects accounting for 3.84 (PGT-9 MA ) to 18.96% (FT-1 MA ) of the phenotypic variance for these traits. Five out of the seven QTLs identified showed QE interactions effects, ranged from 0.64 (for PGT-2 MA in G108) to 2.79% (for PGT-1 MA in F108). In addition, for the AM RIL population (Table 3), six single-locus QTLs were detected, three on Pv01 (FT-1 AM , PGT-1 AM , PST-1 AM ), two on Pv05 (FT-5 AM , PST-5 AM ), and one on Pv03 (PGT-3 AM ). The additive effects of these QTLs explained a phenotypic variance up to 30.74% (FT-1 AM ). Three of them exhibited significant QE interaction effects, which varied from 1.43 (for PST-1 AM in F109) to 18.89% (for PST-1 AM in G108). Most of the QTLs had a negative additive value, except for the QTL PGT-3 AM , which indicated that alleles from the PHA0399 parent mainly enhance flowering and maturity times. Combining both QTL mapping results, two genomic regions stood out, BMD045-FIN and FIN-BMC224 in MA and AM RIL populations, respectively, which contained QTLs controlling FT, PGT, and PST traits that were correlated in both mapping populations (Supplementary Table 4).

Epistatic QTL Interactions
Epistatic and environment interactions among QTLs were detected by means of a two-dimensional genome scan using QTLNetwork 2.0. Thus, a total of 33 and 49 significant E-QTLs involved in 17 and 25 epistatic interactions were identified for the MA and AM RIL populations, respectively (Tables 4, 5). Most of these interactions were due to loci without detectable QTL additive main-effects, and only five and three E-QTLs had both individual additive and epistatic effects in the MA and AM RIL populations, respectively. No significant epistastic interactions were detected for PST, LMS, BL, and PT traits in both RIL populations; whereas E-QTLs were not identified for NPB, ST, NSP, NPP, and SY, and for PWI and SWI traits in the MA and AM RIL population, respectively (Tables 4, 5).
Regarding vegetative growth traits, only one epistatic interaction was identified in the MA RIL population (Table 4), which did not show E-QE effects and explained the 2.34% of the phenotypic variance for the LI trait; whereas four epistatic interactions were detected in the AM RIL population ( Table 5) whose effect ranged from 1.16 (eLI-2 AM × eLI-3 AM ) to 4.73% (eNPB-2 AM × eNPB-10 AM ). These four epistatic interactions displayed E-QE effects that explained up to 5.32% (eNPB-7 AM × eNPB-9 AM in F108) of the phenotypic variance.

DISCUSSION
Many traits of agronomic or biological importance undergo dynamic phenotypic changes during vegetative growth. In fact, the temporal control of flowering initiation determines the time invested in vegetative growth and, consequently, the vegetative resources available during reproduction. In addition, the rate of leaf, flower, and fruit production is a major component of vegetative growth. There is a close association between flowering initiation and vegetative growth; however, how these processes are coordinated during plant development remains unexplored. In this study, the variation for vegetative growth, the rate of leaf, flower, and fruit production, alongside their relationships with flowering and fruit time were investigated in two RIL populations from inter-gene pool crosses of an indeterminate Mesoamerican race Durango and a determinate Andean cultivar race Nueva Granada with the aim to unravel the genetic dynamics underlying vegetative growth and time to flowering.

The Genetic Architecture of Vegetative Growth and Flowering Time in Common Bean
A prerequisite for QTL mapping is the assessment of the quantitative trait in multiple environments. In this study, agronomic evaluation was assessed in four different environments, under open field and greenhouse conditions. The analysis of the two RIL populations under these environments indicates predominantly quantitative inheritance rather than qualitative genes controlling vegetative growth and time to flowering. However, not only additive main effects are responsible for the phenotypic variation observed in our RIL mapping populations, but also epistatic interactions play an important role on the genetic control of flowering time and rate of plant production traits.
In both MA and AM RIL populations, normal distributions were found for most traits, although it was detected that the length of main stem could be controlled by a single gene in the AM RIL population since a bimodal distribution was found in the determinate type I Andean × indeterminate type IV Mesoamerican cross. This bimodal trend is adjusted with the discrete classes and the proportion expected for an autosomal major gene. Simple and complex genetic model has previously been proposed for plant height in common bean. Thus, Kornegay et al. (1992) observed that plant height is a trait of simple inheritance and high heritability. However, Frazier et al. (1958) stated that to reach again the erect trait in a plant of typically determinate growth habit, in addition to the FIN locus, it is needed the action of at least three recessive genes, or a set of minor action genes. Likewise, Davis and Frazier (1966) predicted several genes for internode length, as well as Checa et al. (2006) in indeterminate/indeterminate crosses of Andean and Mesoamerican beans; these authors found that the inheritance of plant height and internode length was mostly additive with only a few genes involved in the expression, and that these genes were most likely modified by interaction with minor genes and with the environment. Checa and Blair (2008) observed a quantitative inheritance rather than qualitative genes controlling plant height in an indeterminate type IV Mesoamerican × indeterminate type II Andean cross, although they also detected a relatively major   gene or single locus with pleiotropic effects on plant height and internode length. Transgressive segregation in both directions was observed for most of the evaluated traits in both RIL populations, except for the tendency of most lines to present averages for LMS and LI traits in the MA RIL population, and for SY and NPP traits in both RIL populations closer to the determinate Andean genotype, which might be related to segregation distortion. Thus, positive and negative transgressive segregations observed suggest that parental lines bear alleles that contribute to vegetative growth and time to flowering variation at several loci. High and significant heritability estimates, as well as positive correlations between LMS with FT and rate of plant production were found. FT was associated with initiation to immature or green pod and physiological pod or dry pod. In spite of the fact that FT and PST correlated with rate of plant production traits, it was not found correlation with PGT. Tar'an et al. (2002) also revealed positive correlations between components of plant height with FT and rate of seed production. Likewise, Koinange et al. (1996) reported that the gene FIN controlling determinacy has pleiotropic effects on FT, PST, and rate of plant production. In this work, it is shown that there is a cause-effect relation among LMS, FT, and rate of plant production traits, suggesting that physically linked or pleiotropic genes might be involved in the regulation of these traits (Aastveit and Aastveit, 1993).
Most of the single-locus QTLs detected in this work overlap with QTLs identified in some quantitative analyses of flowering time and vegetative growth carried out in common bean Checa and Blair, 2008;Wright and Kelly, 2011). In general, these single-locus QTLs were responsible for the majority of the genetic variation for rate of plant production and time to flowering traits within common bean populations. However, it should be noted the presence of a considerable number of epistatic interactions for vegetative growth, plant production and flowering time traits that have not been reported so far. Taken together, these epistatic interactions reinforce the hypothesis that epistasis is involved in the genetic control of agronomical traits as well as epistatic interactions are more frequent in inter-gene pool crosses than in intra-gene pool crosses of common bean (Borel et al., 2016). Thus, for example, significant genetic interactions were found for the genomic region where the FIN locus is located (Pv01) with other genomic regions on Pv02, Pv06, and Pv11 in the MA RIL population, as well as on Pv04 in the AM RIL population. Epistatic interactions for flowering time traits were also detected in other genomic regions, whose effects explained up to 5.73% (eFT-7 AM and eFT-9 AM on Pv07 and Pv09, respectively) of the phenotypic variance.
Seed yield is mainly determined by factors such as number of seeds per pod and seed weight. In this study, the skewness to the lower values shown by the two related-productivity traits, SY and NPP, is in agreement with previous studies (Singh and Urrea, 1995;Johnson and Gepts, 1999;Bruzi et al., 2007), where the biological constraints of the inter-gene pool crosses, Andean and Mesoamerican, hamper reaching the maximum possible yields. Said limitation might result from the loss of favorable epistatic combinations or low probability of recovering superior genotype combination, among other reasons (Johnson and Gepts, 2002;Moreto et al., 2012). According to this hypothesis, the high number of epistatic interactions detected for all yield components is remarkable. Thus, epistatic interaction effects accounted for more than 10% of the phenotypic variance for PL (12.54%) and SL (13.7%) traits in the MA RIL population, as well as for ST (17.29%), SW (15.4%), and NSP (14.63%) traits in the AM RIL population. Hence, results of this research showed the importance of the epistatic effects in the genetic regulation of yield component traits. Thereby, both main and epistatic interaction effects should be considered for a successful application of marker assisted selection (MAS) programs in order to increase yield in common bean.

Rate of Plant Production and Time to Flowering Genetic Links
In order to determine the genetic basis of the rate of plant production during the vegetative growth and time to flowering, QTL mapping was carried out with the average traits estimated in different conditions from flower to seed components. In both  Frontiers in Plant Science | www.frontiersin.org MA and AM RIL populations, it was detected a QTL located on Pv01 at the FIN locus, which showed large additive relative effects on time to flowering traits (up to 19% for FT and 32% for PST of the phenotypic variation in the MA and AM RIL populations, respectively). Comparative QTL analyses showed that this genomic region on Pv01 was also involved in the regulation of vegetative growth traits. Thus, QTLs for LMS, NPB, and LI traits were detected, which explained large additive effects (50-66% of the phenotypic variance for LMS). Within this genomic region, additional QTLs controlling pod size and productivity components (NSP and NPP) which explained more than 8% of the phenotypic variance. Furthermore, in both RIL populations, the FIN locus also displayed large effects on seed yield (19 and 15% of the phenotypic variance for SY trait in the MA and AM RIL populations, respectively). In the MA RIL population, the FIN locus also affected bracteole, leaf, pod, and seed size with a small to moderate additive effect (5 QTLs up to 3% of the phenotypic variance); whereas in AM RIL population, it affected pod size (9% of the phenotypic variance for PWI).
In addition to the FIN locus, other genomic regions were also involved in the regulation of rate of plant production and time to flowering traits, although with a minor effect. Thus, for example, QTLs for PGT and PST in the MA population and for PT in the AM population were detected on Pv02 (PVESTBR006-GATS91), which colocalized with a QTL for seed weight previously mapped by Blair et al. (2006). Taken together, comparative QTL analysis results indicated that vegetative growth has a large effect on time to flowering and the rate of plant production traits, explaining the pleiotropic effects observed for FT and LMS traits. It is mainly due to the phenotypic effect of the recessive allele at the FIN locus, which is present in the determinate type I genotype (Beluga) and controls the meristem switch from a vegetative to a reproductive state. The fin allele reduces the plant growing period, causing a reduction in length of the main stem and number of branches, and an increase in internode length, as well as small bracteoles and leaves, large pods, and seeds that give rise to a lower yield (Norton, 1915;Ojehomon and Morgan, 1969;Koinange et al., 1996), resulting in common bean cultivars with a reduced flowering period since they mature more rapidly (Cober and Tanner, 1995;Koinange et al., 1996).

Genetic and Molecular Mechanisms
Underlying the Link between Time to Flowering and Vegetative Growth Currently, there are genes underlying flowering time QTLs which have been isolated in other crops; said crops have a growth pattern similar to common bean through successive series of modular units. In a common bean plant with a determinate type I growth habit, after floral initiation, the terminal shoot meristem produces a terminal inflorescence and ceases its vegetative growth. However, in a common bean plant exhibiting an indeterminate type IV growth habit, the terminal shoot meristem produces stem nodes, each one composed by one compound leaf and an inflorescence in its axil; thus, vegetative, and reproductive structures are continuously produced until maturity and senescence. This growth pattern through successive series of modular units is similar to that of tomato (Sage and Webster, 1987;Schmitz and Theres, 1999). In tomato, SELF-PRUNING (SP) and FALSIFLORA (FA) control meristem identity. SP gene suppresses the transition of vegetative to reproductive state, keeping a plant indeterminate (Pnueli et al., 1998). FA is responsible for floral meristem identity and promotes flowering (Molinero-Rosales et al., 1999). SP is the homolog of the Antirrhinum majus CENTRORADIALIS (CEN) and A. thaliana TERMINAL FLOWER1 (TFL1) genes (Pnueli et al., 2001). FA is an ortholog of A. thaliana LFY and A. majus FLORICAULA (FLO). LFY in A. thaliana activates directly AP1, causing flowering (Komeda, 2004;Saddic et al., 2006). Foucher et al. (2003) found that a pea homolog of the A. thaliana TFL1, PsTFL1a, corresponds to the determinacy locus (DETERMINATE; DET), and that another pea TFL1 homolog, PsTFL1c (LATE FLOWERING; LF), acts as a repressor of flowering. Likewise, Mir et al. (2014) revealed the orthologous nature of CcTFL1 gene for determinacy in pigeonpea; whereas in soybean, Tian et al. (2010) showed evidence that Dt1 (GmTFL1) is a homolog of the Arabidopsis TFL1 gene, which has a high-level of conservation with the common bean PvTFL1y gene. This gene has been proposed as a candidate gene for the FIN locus since mutations at the PvTFL1y locus were found to cosegregate with the determinate growth habit phenotype (Repinski et al., 2012). In this study, the position of PvTFL1y (Phvul.001G189200) was found within the QTLs positioned in the marker intervals of FIN-BMC224 and BMD045-FIN on Pv01. This finding is consistent with evidence at Andean PvTFL1y haplotype of the determinacy locus (Repinski et al., 2012). However, more information is needed to know whether different PvTFL1y haplotypes derived independently in each gene pool or whether the determinacy locus arose in a single gene pool, as happens in rice. In this species, the determinacy locus arose in a unique gene pool (indica or japonica) and was later transferred to the other pools (Sweeney et al., 2007). Furthermore, the correlated effects of the PvTFL1y locus on other plant traits, such as length of main stem or internode length, and productivity have been demonstrated in this work.

FUTURE REMARKS
The information herein reported could be used not only to establish different breeding strategies combining loci from the different gene pools of common bean, but also to look for associations of genetic variation in determinacy candidate genes in other legume crops with varieties bred for determinate growth habit, such as P. coccineus (runner bean), P. lunatus (lima bean), and Cajanus cajan (pigeonpea) (Waldia and Singh, 1987;van Rheenen et al., 1994;Huyghe, 1998). Furthermore, exploring the interaction and linkage of loci for vegetative growth, flowering time, and the rate of plant production may allow for the expansion of common bean to the geographic locations in which novel adaptation traits can be evaluated.

AUTHOR CONTRIBUTIONS
AG carried out molecular marker and QTL analysis and drafted the manuscript. FY supported mapping methodologies and contributed to a critical review of the manuscript. SS performed the phenotypic data collection. SB contributed to molecular marker analysis. AD contributed to the review of the manuscript. RL collaborated in experimental design and critical review of the manuscript. MS planned the research work, assisted in analysis and interpretation of the data, and edited the manuscript. All authors have read and approved the final version of the manuscript.