Genomic dissection of plant development and its impact on thousand grain weight in barley through nested association mapping.

Flowering time is a key agronomic trait that plays an important role in crop yield. There is growing interest in dissecting the developmental subphases of flowering to better understand and fine-tune plant development and maximize yield. To do this, we used the wild barley nested association mapping (NAM) population HEB-25, comprising 1420 BC1S3 lines, to map quantitative trait loci (QTLs) controlling five developmental traits, plant height, and thousand grain weight. Genome-wide association studies (GWAS) enabled us to locate a total of 89 QTLs that genetically regulate the seven investigated traits. Several exotic QTL alleles proved to be highly effective and potentially useful in barley breeding. For instance, thousand grain weight was increased by 4.5 g and flowering time was reduced by 9.3 days by substituting Barke elite QTL alleles for exotic QTL alleles at the denso/sdw1 and the Ppd-H1 loci, respectively. We showed that the exotic allele at the semi-dwarf locus denso/sdw1 can be used to increase grain weight since it uncouples the negative correlation between shoot elongation and the ripening phase. Our study demonstrates that nested association mapping of HEB-25 can help unravel the genetic regulation of plant development and yield formation in barley. Moreover, since we detected numerous useful exotic QTL alleles in HEB-25, we conclude that the introgression of these wild barley alleles into the elite barley gene pool may enable developmental phases to be specifically fine-tuned in order to maximize thousand grain weight and, potentially, yield in the long term.

Introduction

In the past few decades the phenotypic characterization and genetic dissection of flowering time has been achieved for numerous model species and crops (Blümel ). Barley (Hordeum vulgare ssp. vulgare L.) has been established as a model species for temperate cereals and serves as a solid base for dissecting the genetic basis of flowering time regulation. In this regard, day length (photoperiod) and sensitivity to cold temperatures (vernalization) have been identified as being two major determinants of flowering time. Vrn-H3 (Yan ) is the key gene controlling flower initiation in barley. It is an orthologue of the Arabidopsis FT (FLOWERING LOCUS T) gene. The FT protein moves from leaves to the shoot apical meristem, as postulated decades ago for the mobile signal ‘florigen’ (Chailakhyan, 1937). Its function as a promoter of flowering is assumed to be preserved across different plant species (Turck ).

The key vernalization genes Vrn-H1 (Yan ) and Vrn-H2 (Yan ) have a major impact on Vrn-H3 expression. Vrn-H2 was determined to be a repressor of Vrn-H3, which prevents flowering under long days before vernalization. Vrn-H1 responds to low temperatures (Oliver ) as a result of cis-regulatory elements in its promotor region (Alonso-Peral ). It is up-regulated after vernalization and, thus, promotes flowering through direct binding to the promoters of Vrn-H2 (repression) and Vrn-H3 (activation) (Deng ). A differentiation between winter and spring barley can be made based on their response to vernalization. The latter lacks the vernalization requirement as a result of a natural deletion of Vrn-H2 (von Zitzewitz ).

Flowering is furthermore promoted by Ppd-H1 under long days (Turner ) and by Ppd-H2 under short days (Kikuchi ). In addition to photoperiod and vernalization, light quality (Nishida ; Pankin ), circadian rhythms (Campoli ; Faure ; Zakhrabekova ; Campoli ; Calixto ), and phytohormones like gibberellic acid (GA) (Jia ; Jia ; Boden ) and cytokinins (Mrízová ) also contribute to the induction of flowering in barley. Many loci that were previously classified as earliness per se genes that act independently from external signals (Laurie ) were recently assigned to those classes. The complexity of floral networks (Blümel ) is caused by the interplay between numerous genes and external signals. It is worth taking a closer look at specific developmental subphases to help unravel these networks.

The life cycle of barley consists of several subphases. The most basic differentiation divides the life of a barley plant into vegetative, reproductive and grain-filling phases (Sreenivasulu and Schnurbusch, 2012). Vegetative plant organs develop during the vegetative phase. The reproductive phase starts with the initiation of spikelets, which develop over time. This phase is terminated with anthesis resulting in the onset of the grain-filling phase. The length of different preanthesis subphases has been shown to be under genetic control (Kernich ; Borras ; Borras-Gelonch ; Borras-Gelonch ; Alqudah ) and impacts yield-related traits (Miralles ; Alqudah and Schnurbusch, 2014, 2015). The postanthesis phase is also assumed to have a major impact on yield by determining the time frame for grain-filling (Evans and Wardlaw, 1976; Egli, 2004), which controls the yield component grain weight (Distelfeld ). Although the timing and duration of developmental phases is important for a plant’s yield potential, most genomic studies dealing with the regulation of plant development and its impact on yield focus solely on flowering time as the only developmental parameter.

Nested association mapping (NAM) has been shown to be a valuable tool in the dissection of the genetic architecture of many traits in maize, sorghum, and barley (Buckler ; Jordan ; Maurer ). A NAM population is the result of wide crosses of highly diverse donor genotypes with a recurrent elite cultivar, followed by several rounds of selfing. It combines the advantages of association mapping (i.e. high allele richness and mapping resolution) and linkage mapping (i.e. high statistical power). It also provides excellent opportunities to evaluate the performance of untapped wild alleles in an elite background as a result of its exclusive mating design (Ogut ).

The aim of this study was to characterize the genetic architecture of barley development. To do this, we used the NAM population HEB-25 (Maurer ) consisting of 1420 highly divergent BC1S3 lines. We also wanted to shed more light on the flowering time pathway in barley by comparing the impact flowering time genes had on different developmental subphases. In addition to time to flowering (HEA), we investigated time to shooting (SHO), duration of the shoot elongation phase (SEL), duration of the ripening phase (RIP) and time to maturity (MAT). In order to gain insight into additional physiological functions of flowering time genes and their impact on yield formation, we also investigated plant height (HEI) and thousand grain weight (TGW, Table 1). Furthermore, we looked at whether there was useful variation present in the wild barley germplasm that could be used to fine-tune specific developmental phases in order to increase yield potential in future barley breeding programmes.

Table 1.

List of evaluated traits

Abbreviation	Trait	Unit	Method of measurement	Years studied
SHO	Shooting	days	Number of days from sowing until first node palpable at least 1cm above the tillering node for 50% of all plants of a plot (BBCH 31; Lancashire et al., 1991)	2011–2014
SEL	Shoot elongation phase	GDD	Time from SHO to HEA	2011–2014
HEA	Flowering	days	Number of days from sowing until first awns visible (BBCH 49; Lancashire et al., 1991) for 50% of all plants of a plot	2011–2014
RIP	Ripening phase	days	Time from HEA to MAT	2012–2014
MAT	Maturity	days	Number of days from sowing until hard dough: grain content solid and fingernail impression held (BBCH 87; Lancashire et al., 1991) for 50% of all plants of a plot	2012–2014
TGW	Thousand grain weight	g	Calculated after harvest by use of MARVIN seed analyser (GTA Sensorik GmbH, Neubrandenburg, Germany) based on a 200 seeds sample of each plot. Before, seeds were cleaned and damaged seeds were sorted out	2011–2013
HEI	Plant height	cm	Recorded at maturity as the distance from ground to tip of the erected ear (without awns), taken as an average across the ears of a plot	2011–2013

Materials and methods

Plant material

The NAM population HEB-25 (Maurer ), consisting of 1420 individual BC1S3 lines in 25 wild-barley-derived subfamilies, was used in this study. HEB-25 was the result of initial crosses between the spring barley cultivar Barke (Hordeum vulgare ssp. vulgare) and 25 highly divergent exotic wild barley accessions (Hordeum vulgare ssp. spontaneum and agriocrithon). F1 plants of the initial crosses were backcrossed with Barke. For detailed information about the population design, see Maurer .

Collecting single nucleotide polymorphism data

Single nucleotide polymorphism (SNP) genotype data were collected at TraitGenetics, in Gatersleben, Germany, for all 1420 individual BC1S3 lines and their corresponding parents. The barley Infinium iSelect 9K chip consisted of 7864 SNPs (Comadran ). The genotype data were processed and stored as indicated in Maurer and 5709 informative SNPs, which met the quality criteria, could be utilized in this study. An identity-by-state approach was used to differentiate between the HEB genotypes. Based on parental genotype information, the exotic allele could be specified in each segregating family, and homozygous exotic genotypes were assigned a value of 2. HEB lines that showed a homozygous Barke genotype were assigned a value of 0. Consequently, heterozygous HEB lines were assigned a value of 1. Numbers can therefore be interpreted as a quantitative variable representing the dose of the wild allele.

HEB-25 field trials

Between 2011 and 2014, four field trials were conducted at the ‘Kühnfeld Experimental Station’ of Martin Luther University Halle-Wittenberg (51°29′46.47″N; 11°59′41.81″E) to gather phenotype data. In 2011, the field trial was conducted with selfed progenies of BC1S3 lines (so-called BC1S3:4), arranged as a single randomized block. The majority of plots (92%) included one or two rows per HEB line whereas the remaining 8% of plots included three to five rows, depending on the number of available BC1S3:4 seeds. Plots in 2011 had a length of 1.50 m, a distance of 0.20 m between rows, and were separated by 0.50 m to reduce competition between plots. In 2012 and 2013, the field trials were conducted with the selfed progenies in BC1S3:5 and BC1S3:6, respectively. Two replications per HEB line, arranged in two randomized complete blocks, were cultivated in 2012 and 2013. The plots consisted of two rows (30 seeds each) with a length of 1.50 m, a distance of 0.20 m between rows, and a spacing of 0.50 m between plots. In 2014, BC1S3:7 seeds were sown in a single randomized block. The plots consisted of two rows (50 seeds each) with a length of 1.50 m, a distance of 0.20 m between rows and a spacing of 0.50 m between plots. Barke was integrated as a check line in all trials. All field trials were sown in spring, between March and April, with fertilization and pest management carried out according to local practice. No additional fertilizer was applied in 2014.

Phenotypic data

Table 1 shows a list of all of the investigated traits and a description of their measurements and the years studied. This information is supplemented insofar as all developmental traits were recorded both as ‘days from sowing’ and ‘growing degree days’ (GDD), which were calculated with a base temperature of 0 °C in accordance with equation (1) in McMaster and Wilhelm (1997). Thus, the mean daily temperatures of all of the days with a mean temperature above 0 °C were cumulated.

Statistical analyses

We performed a one-step phenotypic data analysis for all traits with SAS 9.4 (SAS Institute Inc., Cary, NC, USA), based on a linear mixed model (PROC MIXED) with effects for genotype (i.e. 1420 HEB lines), environment (i.e. 4 years) and interaction of genotype and environment. To estimate variance components, all effects were assumed to be random (PROC VARCOMP). Heritabilities across years were calculated as: ,

where V G, V GY and V R represent the genotypic components, genotype × year, and error variance components, respectively. The terms y and r indicate the number of years and replicates, respectively. Best linear unbiased estimates (BLUEs) were calculated with PROC MIXED for each genotype assuming fixed genotype effects. BLUEs were used to calculate Pearson’s correlation coefficients (r) with PROC CORR.

Genome-wide association study

We applied Model B on trait BLUEs as outlined in detail by Liu . This model was found to be most suited for carrying out genome-wide association study (GWAS) with multiple families (Würschum ) and has already been shown to work properly in HEB-25 (Maurer ). It is based on multiple regression, taking into account a quantitative SNP effect and a qualitative family effect in addition to quantitative cofactors that control both population structure and genetic background (Würschum ). Cofactor selection was carried out on this model and included all SNPs simultaneously by applying PROC GLMSELECT in SAS and minimizing the Schwarz Bayesian Criterion (Schwarz, 1978). PROC GLM was used to perform the genome-wide scan for the presence of marker–trait associations. Cofactors that were linked closer than 1 cM to the SNP under investigation were excluded. The Bonferroni–Holm method (Holm, 1979) was used to adjust marker–trait associations for multiple testing. Significant marker effects were accepted with P BON-HOLM<0.05. The proportion of the phenotypic variance explained by a marker was determined by estimating R 2 after modelling the marker solely in a linear model. Additive effects for each SNP were taken as the regression coefficient of the SNP directly from the GWAS model. Family-specific effects were calculated for all markers based on a simple linear model, including a general family term and the marker effect as nested within the family.

To increase the robustness of the method, the entire procedure was applied 200 times on random subsamples of the full dataset. Each subsample included 80% of the lines, randomly selected per HEB family. We recorded the significant (P BON-HOLM<0.05) markers detected, which is referred to as the detection rate. Markers that were detected in at least 10% of subsamples were accepted as putative QTLs, following Ogut . Significant markers were merged into a single QTL if they were linked by less than 4 cM. Additive effects, P BON-HOLM values and R 2 were averaged across all runs, in which the respective marker was significant. In order to evaluate the explained phenotypic variance, the unbiased estimator R 2 adj (Draper and Smith, 1981) was determined for each subsample by simultaneously modelling all of the significant markers in a linear model. In order to determine the predictive ability R 2 pred, the estimated additive effects of each subsample were used to predict the phenotypic value of the remaining 20% of the lines. We then calculated R 2 pred to be the squared Pearson product–moment correlation between predicted and observed phenotypic values. The means of R 2 adj and R 2 pred, measured over 200 runs, were ultimately recorded as the final values.

We used the BARLEYMAP pipeline (Cantalapiedra ) to identify potential candidate genes to explain the QTLs. BARLEYMAP enables markers to be mapped and gene sequences to be aligned against sequence-enriched genetic (Mascher ) and physical frameworks (International Barley Genome Sequencing Consortium, 2012). This represents a very precise way to find positional coincidence of QTLs with putative candidate genes. The results are presented in Supplementary Table S1 at JXB online. We also compared our wild allele effects to those reported in a barley BC2DH population and a set of wild barley introgression lines (Wang ). Here several known flowering time genes were sequenced in wild and cultivated barley. Since the wild allele effect was also produced directly from a cross of wild barley and cultivated barley, we assume that similar effects indicate the same candidate genes.

Results and discussion

Phenotypes

A broad variation in phenotypes of HEB-25 lines was observed for all traits across and within years. This resulted in high coefficients of variation (Supplementary Table S2). After calculating BLUEs, which were corrected for year effects, we observed large differences of more than 100% between the most extreme HEB lines for each trait except MAT (Table 2 and Supplementary Fig. S1). For instance, a difference of 51.4 days between the earliest and the latest genotypes could be observed for HEA, and likewise TGW varied between 19.4g and 60.2g. The lowest coefficient of variation was obtained for the trait MAT, where the most extreme genotypes nevertheless showed a difference of 33.4 days.

Table 2.

Descriptive statistics for best linear unbiased estimates (BLUEs) and heritabilities across all environments

Traita	Nb	Meanc	SDd	Mine	Maxf	CV%g	h 2 h
SHO	1422	53.3	5.6	38.6	82.6	10.4	0.93
SEL	1422	237.8	42.6	108.9	396.0	17.9	0.75
HEA	1422	67.9	6.3	50.4	101.9	9.2	0.94
RIP	1420	32.7	2.6	19.2	40.5	7.9	0.81
MAT	1420	101.3	4.5	88.5	121.9	4.5	0.91
TGW	1422	46.5	5.0	19.4	60.2	10.8	0.57
HEI	1420	63.9	8.5	41.0	100.5	13.3	0.88

Trait abbreviations are given in Table 1.

Number of observations (genotypes).

Arithmetic mean.

Standard deviation.

Minimum.

Maximum.

Coefficient of variation (%).

Heritability.

Phenotypic correlations

We calculated Pearson’s correlation coefficients (r) in order to gain basic insights into the relationships between the different phases of plant development, HEI, and their influence on TGW. We observed very high correlations (0.88–0.93) between SHO, HEA, and MAT (Table 3). This indicates that early shooting lines also tend to be early for other stages. Another peculiarity is that in the case of RIP correlation coefficients were negative for all other developmental stages. Together with the relatively low coefficient of variation for MAT, this indicates that the time of maturity may be predetermined or limited to some extent by environmental factors such as heat and drought. As both SEL and RIP were calculated to be the difference between two other stages, we can draw conclusions about their main determinant. Following this, HEA had the greatest impact on the duration of SEL and RIP, since its correlation coefficients outperformed those of SHO and MAT, respectively. The duration of SEL greatly impacted HEI (r=0.45), indicating that the occurrence of (semi-)dwarf plants was based more on a shortened period of SEL than on a reduced growth rate. No correlation with TGW could be observed for SEL. TGW was positively correlated with RIP (r=0.37), which may be due to an extension of grain-filling, where starch is being stored in the grains (Distelfeld ). Interestingly, HEI was also positively correlated with TGW (r=0.31).

Table 3.

Pearson’s correlation coefficients (r)

	SEL	HEA	RIP	MAT	TGW	HEI
SHO	0.32	0.92	–0.67	0.88	–0.38	–0.01
SEL		0.66	–0.60	0.57	–0.07	0.45
HEA			–0.79	0.93	–0.35	0.17
RIP				–0.54	0.37	–0.19
MAT					–0.28	0.13
TGW						0.31

Bold values indicate significant correlations at P<0.0001. Trait abbreviations are given in Table 1.

Heritabilities

Heritabilities for all investigated traits were calculated over 3–4 years. We observed heritabilities >0.5 for all traits (Table 2). Heritabilities for the developmental traits SHO, HEA, and MAT were almost identical, regardless of whether they were measured in days or growing degree days (GDD) (Supplementary Table S3). However, days outperformed GDD (0.81 vs 0.68) in the case of RIP, while GDD outperformed days (0.75 vs 0.69) in the case of SEL. Thus, for SEL, GDD rather than days offers a better estimate for the relative time needed to fulfil this stage, since plant growth rate and, hence, plant development is based on the interplay between different temperature-dependent biochemical processes (Atwell ; Atkin and Tjoelker, 2003). This is of particular importance, especially when dealing with data from different years and when taking into account that the onset of SEL differs greatly (from 38.6 to 82.6 days after sowing) between individual lines in the HEB-25 population. In contrast, the use of days instead of GDD resulted in a higher heritability for RIP. This observation may be attributed to the fact that plant development does not benefit from higher temperatures if a certain temperature threshold is met or the plant has reached a critical physiological state. McMaster and Smika (1988) compared temperature-dependent and -independent models to predict growth stages in winter wheat and pointed out that the best model for predicting a developmental stage varied depending on the respective stage. Therefore, we decided to concentrate on GDD for SEL and days for all other developmental traits in our analyses.

GWAS

We conducted GWAS for each trait in order to further analyse the above-mentioned correlations between traits and to elucidate which QTLs are responsible for controling trait variation in HEB-25. We applied a multiple linear regression model, including a population main effect and selected markers as cofactors, to account for genetic background and relatedness. This model was recently shown to perform best in HEB-25 (Maurer ).

We were able to detect numerous associated genomic regions for all of the traits studied using our GWAS method (Fig. 1 and Supplementary Table S4). A total of 89 QTLs could be defined (Table 4 and Supplementary Table S5). Most QTLs were shared by multiple traits. However, we could also detect trait-specific QTLs for all of the seven traits (Supplementary Fig. S2). We obtained broad variation for wild allele QTLs, with an increase or decrease in trait values compared with the Barke control allele. In this regard, we were able to identify 10 QTLs where the wild allele increased TGW by up to 6.6g (Supplementary Table S5). Modelling all significant markers of one trait simultaneously resulted in explained percentages of phenotypic variance (R 2 adj) ranging from 63.6% (TGW) to 82.3% (HEI). The fact that R 2 pred values, resulting from a prediction of phenotypes of an independent sample, were also comparatively high (Table 4) confirms the robustness of the method and indicates that a large fraction of phenotypic variation can be explained by genotypic variation.

Fig. 1.

Comparison of GWAS results across developmental traits, thousand grain weight and plant height. Barley chromosomes are indicated as coloured bars on the inner circle, and centromeres are highlighted as transparent boxes. Grey connector lines represent the genetic position of SNPs on the chromosomes, which is given in centimorgans on the outer circle. Each track represents one trait, and these are (from inside to outside) SHO, SEL, HEA, RIP, MAT, TGW and HEI. Trait abbreviations are given in Table 1. Black boxes indicate the QTL positions. The height of histogram bars above represent the detection rate across 200 repeated random subsamples. The blue and red colours of the bars indicate trait-reducing and trait-increasing effects, respectively, exerted by exotic QTL alleles. Candidate genes of major QTLs are indicated outside the circle.

Table 4.

Number of QTLs and total explained phenotypic variance

Traita	QTLs (n)b	R 2 adj (%)c	R 2 pred (%)d
SHO	49	81.6	59.9
SEL	28	65.9	42.5
HEA	43	79.0	55.5
RIP	30	66.6	37.5
MAT	32	76.9	56.2
TGW	43	63.6	42.8
HEI	37	82.3	45.7
No. of unique QTLs	89

Trait abbreviations are given in Table 1.

Number of QTLs detected for the respective trait.

Mean explained phenotypic variance by GWAS.

Mean ability to predict phenotypes of an independent sample.

Below, we present eight common QTLs (starting with chromosome 1H in ascending order) and discuss their relevance in controlling the five developmental traits SHO, SEL, HEA, RIP, and MAT. Then, in order to draw conclusions about their potential QTL function, we study their effects on HEI and TGW. Since HEB-25 enables high genetic precision in mapping of QTLs (Maurer ), we obtained strong positional coincidence with plausible candidate genes (Supplementary Table S1). We therefore directly refer to these candidate genes in the headings of the subsections. The results of these eight major QTLs are summarized in Table 5 and illustrated in Fig. 2. Finally, we discuss the developmental phase-specific QTLs that were obtained in this study.

Table 5.

Major developmental QTLs and their impact on further traits

			cM intervalc		SHOd	SEL	HEA	RIP	MAT	TGW	HEI	Candidate gene/locus with reference
QTL	Chra	Peak markerb	From	Until
QTL-1H-10	1H	i_SCRI_RS_150786	128.0	133.1	–1.9	–10.2	–2.3		–3.4		–3.2	HvELF3	(Faure et al., 2012; Zakhrabekova et al., 2012)
QTL-2H-4	2H	i_BK_12	22.2	23.8	-7.4	–46.7	–9.3	2.8	–6.8		–7.3	Ppd-H1	(Turner et al., 2005)
QTL-2H-7	2H	i_12_30265	53.3	60.8	–1.9	–39.8	–3.6	1.9	–2.6		–4.7	HvCEN	(Comadran et al., 2012)
QTL-3H-9	3H	i_11_11172	103.8	109.8	–5.7	44.9	–4.3	0.9	–4.0	4.5	12.3	denso/sdw1	(Jia et al., 2009)
QTL-4H-1	4H	i_12_31458	0.6	14.9	2.0	15.4	2.5	–0.7	2.0	–3.5	–2.6
QTL-4H-9	4H	i_SCRI_RS_216897	110.2	114.3	7.8	11.3	8.2	–3.1	5.5	–1.5	4.2	Vrn-H2	(Yan et al., 2004)
QTL-5H-10	5H	i_11_10783	122.4	128.5	7.9		8.5	–3.5	5.2	6.6		Vrn-H1	(Yan et al., 2003)
QTL-7H-3	7H	i_12_30895	29.8	34.3	3.4	36.8	5.7	–2.1	3.5		2.5	Vrn-H3	(Yan et al., 2006)

Negative values are indicated in bold. Blank cells indicate that the respective QTL was not detected for the trait. Trait abbreviations are given in Table 1. For a complete table of all QTLs, see Table S5.

Chromosome on which the QTL was detected.

iSelect name of marker with the highest significance for HEA.

Genetic interval (cM) with lower and upper threshold of QTL, based on the map of Maurer et al. (2015).

Most extreme effect (absolute difference of homozygous wild genotype and homozygous cultivated genotype) of all significant SNPs in respective QTL interval.

Fig. 2.

Spider diagram of major QTL effects across traits. Different traits are represented by corner points of the net. Trait abbreviations are given in Table 1. Effects of wild alleles are indicated by differentially shaped lines for the respective QTL. The blue-shaded and red-shaded backgrounds of the spider net indicate trait-reducing and trait-increasing effects, respectively, exerted by exotic QTL alleles. To enable a comparison of traits within the same scale, values of SEL have been divided by 16.3, which represents the equivalent of one GDD to one day during SEL. Values of HEI have been divided by 2.

QTL-1H-10 (HvELF3)

We observed a QTL that showed significant effects for all seven traits studied except for RIP and TGW. The QTL is located close to the telomere of chromosome arm 1HL (QTL-1H-10). The exotic allele at this locus is associated with a slight acceleration of plant development. The time to reach each stage was shortened by 2–3 days in contrast to the allele of the spring barley cultivar Barke. This QTL region harbours the earliness-inducing eam8/mat-a locus (Franckowiak ) and may correspond to HvELF3 (EARLY FLOWERING3), which is orthologous to the Arabidopsis circadian clock gene ELF3 (Faure ; Zakhrabekova ). HvELF3 was recently shown to influence flowering by regulating GA production in barley (Boden ).

QTL-2H-4 (Ppd-H1)

QTL-2H-4 exerted significant effects on all developmental traits and HEI, thereby explaining up to 34% of the phenotypic variance (R 2). The most significant SNP marker (i_BK_12) is directly located within the Ppd-H1 gene, which is the main determinant of response to long day conditions in barley (Turner ). Most wild barley accessions were shown to carry a highly responsive Ppd-H1 allele, accelerating development under long day conditions (Cockram ). In contrast, Barke, like most spring barley from Northern Europe, carries an allele that exhibits a reduced response to long days. In our study, the wild allele led to a strong acceleration of all developmental phases by up to 9.3 days, except for RIP, which was delayed by up to 3 days. In contrast, HEI was reduced by 7.3cm when the wild Ppd-H1 allele was present. After reaching the reproductive stage, relatively more energy is put into reproductive organs. This leads to reduced vegetative growth. Thus, the wild allele at Ppd-H1 may accelerate floral development at the expense of growth and biomass production.

QTL-2H-7 (HvCEN)

We observed a QTL, located next to the centromere of chromosome 2H, that showed significant effects on all traits except TGW. Comadran identified HvCEN (CENTRORADIALIS), a homologue of Arabidopsis TFL1 (TERMINAL FLOWER1), as the possible gene behind this locus. Recently, Loscos found evidence that HvCEN plays a central role in the induction of flowering in barley. The effects of HvCEN on all developmental phases were similar to the effects of Ppd-H1. However, the scale of the effects was clearly reduced. At this point we should note that the HvCEN effect is presumably not due to genetic linkage to Ppd-H1, since Ppd-H1 was selected as a cofactor in the GWAS model for calculating the HvCEN association.

QTL-3H-9 (denso/sdw1)

QTL-3H-9 was detected for all five developmental traits. This QTL may correspond to the denso/sdw1 locus, causing a semi-dwarf phenotype. HvGA20ox2, coding a GA-20 oxidase enzyme, is a candidate gene for explaining its function (Jia ; Jia ). The denso allele is commonly present in modern European malting barley cultivars like Barke (Jia ). Therefore, it is not surprising that we also observe a strong effect on HEI for this locus, explaining 41% of phenotypic variance (R 2). The presence of the wild allele increased HEI by up to 12.3cm. This QTL simultaneously affected all developmental traits and TGW. The wild allele reduced the time required to reach SHO, HEA and MAT by 5.7, 4.3 and 4.0 days, respectively. At the same time it delayed the time between these stages (SEL and RIP). Furthermore, the wild allele increased TGW by up to 4.5g. As the absolute effects steadily diminished throughout the developmental stages (Table 5), we presume that denso/sdw1 plays a major role in the very early stages of development. All in all, the influence of denso/sdw1 on all investigated traits underlines its complex and important role in plant physiology throughout a plant’s life cycle—a known feature of GA (Taiz and Zeiger, 2007). In addition to the importance of GA in floral regulatory networks (Mutasa-Göttgens and Hedden, 2009), elevated levels of GA also play a role in the delay of senescence, as shown for pea (Wang ). Furthermore, GA levels were shown to be crucial for endosperm differentiation and for the growth of barley grains (Weier ).

QTL-4H-1

QTL-4H-1 is another ubiquitous QTL that has recently been identified as being a major QTL for flowering time in barley (Maurer ). The position of the most significant peak markers differs between 3.5 and 14.9 cM depending on the investigated trait. Nevertheless we assume that both peaks correspond to the same QTL due to a lack of marker density in this region. This QTL has also been found for TGW and HEI, underscoring its crucial role in plant physiology. This QTL delayed all developmental stages by approximately 2 days, except for RIP, which was accelerated by approximately 2 days. Simultaneously, TGW was reduced by 3.5g and HEI was reduced by 2.6cm. To date, no clear candidate gene for this QTL has been found. However, after aligning the sequences of significant markers in this region against sequence-enriched genetic and physical frameworks via BARLEYMAP (Cantalapiedra ) a CCT (CONSTANS, CO-like, and TOC1) domain gene was identified in this region (Genbank accession number AK354746). Many known flowering time regulators have CCT domains (Cockram ). Furthermore, the LOG (LONELY GUY) gene (MLOC_45038.2) could be assigned to this region. LOG encodes a cytokinin-activating enzyme that is required to maintain meristem activity. Its loss of function causes pre-mature termination of the shoot meristem in rice development (Kurakawa ).

QTL-4H-9 (Vrn-H2)

QTL-4H-9, which may harbour Vrn-H2, the main repressor of flowering in vernalization-dependent barley (Yan ), was reliably detected for SHO, HEA, RIP and MAT. The wild allele delayed the respective developmental stages by up to 8 days, whereby RIP was shortened by up to 3 days. However, the high variation of marker effects among significant SNPs indicates that effects vary greatly depending on the segregating families of the respective SNPs (Supplementary Table S4). This locus is naturally deleted in spring barley cultivars like Barke. The result is that no vernalization is required to induce flowering (von Zitzewitz ). This explains why the alleles from wild barley, which are predominantly winter types (Cockram ), delayed flowering in our spring-sown field trials. The fact that this QTL was infrequently detected for the traits SEL, TGW and HEI leads us to assume that the impact of Vrn-H2 on these traits may be diminished or biased by family-specific effects.

QTL-5H-10 (Vrn-H1)

This QTL was shared between SHO (+7.9 d), HEA (+8.5 d), RIP (–3.5 d) and MAT (+5.2 d) and may correspond to the Vrn-H1 locus with HvBM5A (MADS-box 5A) being a candidate gene (von Zitzewitz ). Vrn-H1 is known to be involved in the vernalization pathway of flowering time regulation by responding to the low temperatures required for vernalization (Oliver ). HvPhyC (PHYTOCHROME C) is another gene that affects flowering time (Nishida ; Pankin ). It is closely linked to HvBM5A, which makes it difficult to distinguish between their effects. In addition, a QTL for HEI (QTL-5H-11) has been defined approximately 3 cM distal from the position of HvBM5A and HvPhyC. This further complicates the interpretation of the QTL effects in this genomic region. Again, as already shown for Vrn-H2, the high variation of marker effects points to family-specific effects for this locus. An extraordinarily high variation exists for the trait TGW, which is caused by a single SNP marker (i_11_11090) that is only polymorphic in HEB-F08. The wild allele caused an estimated increase in TGW by 6.6g, whereas all other significant markers in the Vrn-H1 region reduced TGW by 1.5–1.8g. Interestingly, there is low variation for the locus affecting HEI, located 3 cM distal from Vrn-H1. This indicates that this locus might be independent of Vrn-H1.

QTL-7H-3 (Vrn-H3)

The most significant SNP marker at this QTL (i_12_30895) is directly located within Vrn-H3. The Vrn-H3 locus in barley has been shown to correspond to HvFT1, the orthologue of Arabidopsis FT (Yan ; Faure ). This gene plays a central role in the flowering pathways as it is involved in the switch from vegetative to reproductive growth under long day conditions (Turck ). In our population, significant effects of this QTL were observed for all developmental traits and HEI. On average, the wild allele in the Vrn-H3 region delayed all developmental phases by 3–6 days, except RIP, which was shortened by approximately 2 days. However, ample variation in marker effects exists for all developmental traits (Supplementary Table S4). For instance we observed effects for SHO ranging from –3.3 to +3.4 days depending on the families in which the respective marker segregates. This clearly shows the presence of family-specific effects, most likely due to functionally different haplotypes of Vrn-H3 among the HEB-25 donor accessions. However, for the trait HEI (increased by 2.5cm) the effect of the wild allele was relatively constant. Interestingly, gene-specific markers of Vrn-H3 did not impact TGW. Instead, the QTL for TGW was defined at around 4 cM proximal of Vrn-H3 at 38.8 cM (QTL-7H-4) causing a reduction of TGW by up to 2.2g.

Developmental phase-specific QTLs

In addition to the majority of QTLs that were shared across developmental traits, we also found certain QTLs that impacted a specific stage of plant development. These QTLs are a potential source for fine-tuning the phases of plant development and unravelling the physiological pathways.

QTLs affecting early development

QTLs that affect early development are characterized by their influence on SHO and SEL. SEL is strongly related to the phase from awn primordium to tipping, which was shown to be the most decisive developmental phase for spikelet survival (Alqudah and Schnurbusch, 2014) and leaf growth rate (Alqudah and Schnurbusch, 2015) in barley. The length of this phase is therefore assumed to play a key role in determining a plant’s yield potential.

We found three striking QTLs (QTL-1H-6, QTL-3H-14 and QTL-5H-16) that exhibited exclusive effects on SHO and SEL. In the case of QTL-3H-14, SHO was accelerated by approximately 2 days and TGW increased by 1.3g. In the case of QTL-1H-6 and QTL-5H-16, SHO was delayed by 1–2 days without any effect on TGW. These QTLs were not shown to affect HEA since the effect on SEL was contrary to the effect on SHO. This therefore compensated for the effect of SHO in the subsequent developmental phases. Thus, we assume that these QTLs play a role in very early plant development by affecting the time to reach SHO and the duration of SEL.

We also observed QTLs specific to the SEL trait. For instance, QTL-2H-10 extended SEL by 19.7 GDD, while QTL-4H-5, QTL-5H-2 and QTL-5H-7 shortened SEL by 12.2, 12.1 and 18.0 GDD, respectively. Interestingly, some of the above-mentioned QTLs also affected TGW and HEI in addition to affecting SEL. In the case of QTL-5H-2, the shortened period of SEL correlates with reduced HEI. Although QTL-4H-5 and QTL-2H-10 show contrasting effects on SEL, they co-localized with QTLs that decrease TGW. In the case of QTL-2H-10 this effect is likely caused by the six-rowed spike locus Vrs1 (Komatsuda ) originating from the Hordeum vulgare ssp. agriocrithon donor of HEB-F24 in our population. Six-rowed spikes generally show decreased TGW due to a higher grain number. This causes single grains to compete for assimilates. Interestingly, Vrs1 also seems to extend SEL. This corroborates the observation that six-rowed barley has a higher leaf area and leaf dry weight than two-rowed barley, since leaf biomass is mainly produced during SEL (Alqudah and Schnurbusch, 2015).

QTLs affecting late development

QTLs that affect late development are thought to influence RIP and MAT. Following Davies and Gan (2012), whole plant senescence is initiated through a physiological transition at flowering. This leads to a remobilization of nutrients to the developing seeds (Distelfeld ). In the present study, RIP represents the period of whole plant senescence and grain-filling. Thus, it is considered to be highly important in determining yield, particularly as a result of the grain weight component (Distelfeld ). Egli (2004) also emphasized that the duration of this period played an important role in grain crop yields.

The most obvious QTL to specifically affect late development was detected in the centromeric region of chromosome 6H (QTL-6H-4). In addition to affecting RIP and MAT, this QTL also significantly affected TGW and HEI. HvNAM-1 (NO APICAL MERISTEM-1), coding for an NAC transcription factor (Distelfeld ), is a promising candidate for this locus. It plays an important role in inducing senescence (Distelfeld ). Moreover, it has been shown that its wheat homologue NAM-B1 influences the remobilization process of nitrogen, iron and zinc in the developing grain (Uauy ; Waters ). In our study the wild allele accelerated senescence, which shortened RIP by approximately 1 day and consequently led to earlier MAT (–3.4 days). Simultaneously, TGW was reduced by 2.6g and HEI was reduced by 4.5cm. Interestingly, there is a second QTL for RIP (QTL-6H-5) that is closely linked to HvNAM-1. It mainly affects early stages by delaying SHO and HEA and shortening RIP. This behaviour points to the candidate gene HvGR-RBP1 (GLYCINE-RICH RNA BINDING PROTEIN), which is consistent with the findings of other groups (Jukanti ; Streitner ; Mason ). This is the first time in barley genomics that we have been able to clearly distinguish this adjacent locus from HvNAM-1.

Comparison of shared QTLs

In previous sections we mentioned many QTLs that simultaneously affected several developmental traits. Comparing their effects across traits provides insights into the different modes of action of the respective QTLs (Fig. 1 and Table 5).

In general, most QTLs affecting HEA also have a similar effect on SHO. This is also reflected by the high correlation of r=0.92 (Table 3). However, comparing the absolute effects on both stages, we can state that the effect on HEA is generally more pronounced. This is reflected by an increased SEL and indicates that the effect may accumulate over time. However, the denso/sdw1 locus is a striking exception to this rule. Although both time to SHO and HEA are reduced, it simultaneously extends the time needed for SEL. This indicates that the effect on SHO is stronger than on HEA at denso/sdw1. We therefore assume that denso/sdw1 has the biggest impact in the very early developmental stages and that its influence decreases during subsequent developmental phases.

Almost half of all developmental QTLs in our study simultaneously showed significant effects on HEI. However, the direction in which HEI was influenced differed between QTLs, even if they shared the same tendencies with regard to developmental traits. This is interesting when it comes to separating general developmental effects from purely reproduction-promoting effects. In most cases, HEI increased as the duration of SEL increased and vice versa. This is already indicated by the positive correlation of both traits (r=0.45). However, in the case of QTL-4H-1 and Vrn-H2, the prolonged duration of SEL is accompanied by a reduction in HEI. Since these QTLs retard development and simultaneously reduce HEI, we assume that their gene products may generally hamper plant development.

Relevance of developmental subphases with regard to TGW

We observed numerous QTLs that were significantly associated with the developmental subphases SEL and RIP. We observed a general influence on TGW for the duration of grain-filling, i.e. RIP; as RIP increased, TGW increased as well. However, we could not find any direct correlations between SEL and TGW. There were hints that the interplay between RIP and SEL may be of interest. In general, we observed a negative correlation between both traits with r=–0.6. Consequently, a shortened SEL would cause an extended RIP. This would then cause a longer period of grain-filling and thus increase TGW. The time of physiological maturity (MAT) is largely predetermined by environmental factors, as discussed above. Therefore, increasing the duration of RIP leads to earlier flowering, and very likely also to a shorter period of SEL. Thus, if the duration of RIP is extended to increase TGW, this will result in a reduction in SEL. This lowers the amount of biomass that would be able to supply the developing grains with assimilates (Alqudah and Schnurbusch, 2015). Normally this reduction would not have serious consequences since plants produce excess vegetative mass (Egli, 1998). However, if the effect is severe (e.g. in the case of Ppd-H1 in our study), the beneficial effect on TGW may be compensated. We therefore presume that increasing TGW is based on a trade-off between an extended RIP and a shortened SEL. However, one has to keep in mind that yield formation is based on a complex interplay between several yield components. To scrutinize our findings, we therefore suggest studying all yield components simultaneously, i.e. tiller number, number of grains per spike and grain weight, as well as total grain yield, in future HEB-25 field experiments. In addition, spike photosynthesis may also require analysis, since its role in grain-filling has recently been highlighted (Sanchez-Bragado ).

Beneficial wild germplasm present in HEB-25

One major advantage of using HEB-25 is the ability to directly evaluate the value of wild barley alleles in a cultivated background. Thus, we were able to identify 10 QTLs where the wild allele increased TGW by up to 6.6g compared with the elite Barke allele. These wild QTL alleles are promising candidates for improving future barley grain weight once they are introduced into barley breeding programmes. However, the optimum strategy of adjusting plant developmental subphases to increase grain weight may vary depending on the ecogeographic region in which barley is grown. As latitudes increase, alleles conferring a decelerated plant development would be favourable for making use of the longer and cooler season for grain-filling. For barley cultivated in lower latitudes, early flowering would be favourable, for instance, to escape early season terminal drought (Comadran ). Therefore, it is currently not possible to make a general statement about the usefulness of specific plant development QTL effects. However, as we see that there is tremendous effect variation for every trait, we conclude that HEB-25 harbours a multitude of beneficial plant development QTL alleles for different ecogeographic regions.

In this study, the most prominent example for putative advantageous effects of wild alleles is the denso/sdw1 QTL on chromosome 3H. Barke carries the dwarfing denso allele, which was introduced in barley breeding programmes during the ‘Green Revolution’. It reduces plant height and improves resistance to lodging, whereby also yield is increased (Jia ). The wild allele extends SEL and RIP, and increases TGW, compared with the cultivated Barke allele. Unfortunately, these beneficial effects are accompanied by the unfavourable effect of increasing HEI. However, we observed lines (e.g. HEB_25_050) in HEB-25 that carry the wild allele at the denso/sdw1 locus and nevertheless demonstrate high agronomic performance. This is most likely due to introgression of additional wild alleles, which appear to compensate for the effect of increasing HEI (Supplementary Table S6). We thus encourage breeders to integrate wild germplasm into their breeding programmes. This way, the elite barley gene pool can be replenished with new favourable alleles to overcome future agricultural challenges.

Supplementary data

Supplementary data are available at JXB online.

Figure S1. Frequency distributions of BLUEs for all investigated traits, plotted as histograms.

Figure S2. Venn diagrams indicating the number of shared QTLs across traits.

Table S1. Results of BARLEYMAP alignments to detect QTL candidate genes.

Table S2. Descriptive statistics of all investigated traits, grouped by years.

Table S3. Heritabilities (h 2) of all investigated traits, including the comparison of days and GDD for developmental traits.

Table S4. Tabular overview of all results gathered from GWAS, along with family-specific effect estimates for each SNP.

Table S5. Overview of all QTLs and their effects, with position, estimated effect on the different traits, and plausible candidate genes.

Table S6. Raw phenotype data and BLUEs of the investigated traits across years and blocks, along with information about the allelic state at the denso/sdw1 locus.