Multi-Locus Genome-Wide Association Studies for 14 Main Agronomic Traits in Barley

The agronomic traits, including morphological and yield component traits, are important in barley breeding programs. In order to reveal the genetic foundation of agronomic traits of interest, in this study 122 doubled haploid lines from a cross between cultivars “Huaai 11” (six-rowed and dwarf) and “Huadamai 6” (two-rowed) were genotyped by 9680 SNPs and phenotyped 14 agronomic traits in 3 years, and the two datasets were used to conduct multi-locus genome-wide association studies. As a result, 913 quantitative trait nucleotides (QTNs) were identified by five multi-locus GWAS methods to be associated with the above 14 traits and their best linear unbiased predictions. Among these QTNs and their adjacent genes, 39 QTNs (or QTN clusters) were repeatedly detected in various environments and methods, and 10 candidate genes were identified from gene annotation. Nineteen QTNs and two genes (sdw1/denso and Vrs1) were previously reported, and eight candidate genes need to be further validated. The Vrs1 gene, controlling the number of rows in the spike, was found to be associated with spikelet number of main spike, spikelet number per plant, grain number per plant, grain number per spike, and 1,000 grain weight in multiple environments and by multi-locus GWAS methods. Therefore, the above results evidenced the feasibility and reliability of genome-wide association studies in doubled haploid population, and the QTNs and their candidate genes detected in this study are useful for marker-assisted selection breeding, gene cloning, and functional identification in barley.


INTRODUCTION
Barley (Hordeum vulgare L, 2n = 2x = 14), one of the first domesticated grains in the Fertile Crescent (Zohary et al., 2012), has been used widely as animal feed, human health foods, and a source of beer. Its yield and quality are the most important breeding objectives in crop breeding programs.
Most agronomic traits, such as plant height-and yield-related traits, are controlled by quantitative trait loci (QTLs) in barley, so it is difficult to obtain their genetic foundation and molecular mechanism. Plant height and its component traits serve as major plant morphological traits affecting barley seed yield. An appropriate plant height is a prerequisite for obtaining the desired yield in barley-breeding programs. To date, more than 30 types of dwarfing or semidwarfing genes have been detected, while only a few have been successfully used in barley breeding program, such as uzu and sdw1/denso (Jia et al., 2009;Ren et al., 2010Ren et al., , 2013. Moreover, a large number of QTLs for plant height related traits were reported to be located on all the seven chromosomes (Sameri et al., 2006;Baghizadeh et al., 2007;Wang et al., 2010;Ren et al., 2014). Grain yield is the key trait for the breeder in barley breeding program, therefore, the yield related traits including spike number per plant (SP), grain number per plant (GP), grain weight per plant (GWP), and 1,000 grain weight (TGW), have gained more attentions in the genetic dissection of yield related traits. A huge number of QTLs for yield related traits were detected to be located across all the chromosomes (Li et al., 2005;Sameri et al., 2006;Baghizadeh et al., 2007;Wang et al., 2010Wang et al., , 2016aRen et al., 2013).
Traditionally, QTL mapping has been widely applied in the genetic dissection of quantitative traits in barley (Zhuang et al., 1997;Li et al., 2005Li et al., , 2007Peng et al., 2011). As the development of DNA sequencing technologies, it is relatively easy to obtain high-density SNP genotypes for association mapping (AM) population, which offers a huge convenience for genomic and genetic research in different species. Therefore, genomewide association studies (GWAS) present a powerful tool to reconnect the complex quantitative traits with their genes. Due to the development of cheaper, faster and higher-throughput molecular markers, AM has been widely used for mapping QTLs and genes in many crops, such as maize, soybean, rice, barley and wheat (Huang et al., 2010;Yang et al., 2010;Pasam et al., 2012;Hu et al., 2015). In comparison with traditional QTL mapping, AM has three obvious advantages, including shorter construction time, much higher mapping resolution and a greater number of alleles (Zhang et al., 2005;Yu and Buckler, 2006). In barley, AM has been widely applied for complex traits including disease resistance (Massman et al., 2011), drought tolerance (Varshney et al., 2012;Wójcik-Jagła et al., 2018), salinity tolerance (Fan et al., 2016) and especially agronomic traits (Gawenda et al., 2015;Xu et al., 2018).
Beside for natural population, nowadays GWAS have been widely applied to the genetic analysis for complex traits in familybased populations, such as nested association mapping (NAM) and multi-parent advanced generation intercross (MAGIC) populations, and proved to be powerful tool for uncovering the basis of key agronomic traits in maize and barley (Tian et al., 2011;Cook et al., 2012;Maurer et al., 2015Maurer et al., , 2016. Moreover, there are also successful cases that combine linkage analysis with GWAS in several bi-parental segregation populations, such as recombinant inbred line (RIL) population (Lu et al., 2010;Reif et al., 2010). For single segregating population, successful but fewer cases were performed using GWAS (Gao et al., 2015;Henning et al., 2016;Liu et al., 2018). Henning et al. (2016) conducted GWAS for downy mildew resistance in a segregating population of "Teamaker" × USDA 21422M in hop (Humulus lupulus L.), Gao et al. (2015) determined the location of TTKSK resistance in the 108 doubled haploid (DH) lines using a GWAS method implemented by R package rrBLUP (Endelman, 2011), andLiu et al. (2018) used two strategies (QTL mapping and GWAS) to reveal the genetic bases of fiber quality traits and yield components in 231 RILs. Therefore, it is feasible to use GWAS to dissect the genetic foundations of complex traits in single bi-parental segregating population. A combination of linkage and association methodologies should provide the more accurate and powerful approach for revealing the genetic bases of complex traits (Ott et al., 2011). Association mapping of drought tolerance-related traits was performed in barley to complement a traditional bi-parental QTL mapping study in Wójcik-Jagła et al. (2018).
The objectives of this study were to: (a) use GWAS to further dissect the genetic foundations for main agronomic traits in our previous studies of Ren et al. (2013Ren et al. ( , 2014 and Wang et al. (2016a), and compare the quantitative trait nucleotide (QTN) results with those in previous studies, (b) evaluate if GWAS is feasible and reliable in the genetic dissection of complex traits in DH population, and (c) mine the candidate genes in the regions of the QTNs. The outcome of this study will provide more precise and complete information for further gene cloning, and marker-assisted selection in barley breeding.

Plant Materials and Field Experiments
One hundred and twenty-two DH lines, derived from a cross between barley cultivar "Huaai 11" (six-rowed and dwarfing) and barley cultivar "Huadamai 6" (two-rowed), was used in this study. The details of the materials and field experiments were described in the previous studies of Ren et al. (2010Ren et al. ( , 2013Ren et al. ( , 2014 and Wang et al. (2016a).

Phenotyping Data
Fourteen agronomic traits for the above DH population were measured in 2008-2009, 2009-2010, and 2011-2012, and all the three datasets had been reported by Ren et al. (2013Ren et al. ( , 2014 and Wang et al. (2016a). These traits included plant height (PH), first internode length (IL1), second internode length (IL2), third internode length (IL3), fourth internode length (IL4), main spike length (MSL), spike number per plant (SP), spikelet number of main spike (SMS), spikelet number per plant (SLP), grain number per plant (GP), grain number per spike (GS), grain weight per plant (GWP), grain weight per spike (GWS), and 1,000 grain weight (TGW). All the details have been described in Ren et al. (2013Ren et al. ( , 2014. The best linear unbiased predictions (BLUPs) for each trait of 3 years were calculated using the R package Lme4 (Bates et al., 2014) with the following model: y = lmer (Trait ∼ (1|Genetype) + (1|Year)). The three single-year phenotypic values (Ren et al., 2013(Ren et al., , 2014Wang et al., 2016a) and their BLUP values were used for GWAS. The results of phenotype statistics of BLUP for each trait were summarized in Table S1. used in this study. Base on the recent genome sequence release in barley (Ibsc, 2016;Beier et al., 2017;Mascher et al., 2017), all the SNP markers were aligned to the most reliable genome position (http://webblast.ipk-gatersleben.de/barley_ibsc/). The Vrs1 locus controlling row number of barley was integrated with SNP markers for GWAS. All the above information has been described in Ren et al. (2016).

GWAS
Q matrix was calculated by STRUCTURE software (Falush et al., 2003), and the optimal K was inferred in Figure S1. The kinship (K) matrix between the lines was calculated as previously described in Wang et al. (2016b). All the 9680 SNPs for the above 122 DH lines were used to conduct GWAS for the above 14 traits in 3 years and their BLUP values using five multilocus GWAS methods, including mrMLM (Wang et al., 2016b), FASTmrMLM (Zhang and Tamba, 2018), FASTmrEMMA , pLARmEB  and ISIS EM-BLASSO , which were included in the R package mrMLM v3.1 (https://cran.r-project.org/web/packages/ mrMLM/index.html). All parameters in GWAS were set at default values. The critical thresholds of significant association for the five methods were set as LOD = 3 (or P-value = 2 × 10 −4 ; Wang et al., 2016b).
The significant QTNs, repeatedly detected in at least two environments or methods, were viewed as reliable. The associated regions on chromosomes, repeatedly located on same or similar traits in at least 2 years or methods, were viewed as reliable QTN clusters. These QTNs (or clusters) were named as "qtn (qtnc)" + trait name abbreviation + chromosome + detected QTL order on chromosome.

Phenotypes Difference Corresponding to QTNs
For each QTN, all the DH lines were firstly divided into two groups based on their QTN genotypes, then t-test was used to test the phenotypic difference between the two genotypes.

Identification of Candidate Genes
According to the recent genome sequence release of barley (Ibsc, 2016;Beier et al., 2017;Mascher et al., 2017) and the gene annotation information (http://plants.ensembl.org/Hordeum_ vulgare/Info/Index and https://www.uniprot.org/uniprot), some genes around reliable QTNs (or clusters) were selected for each trait. By combining gene annotation information, protein domain function in database and previous reports, and expressional information (http://barlex.barleysequence.org), then, candidate genes for each trait were mined.

GWAS for 14 Agronomic Traits
Using five multi-locus GWAS methods in the R package mrMLM v3.1, GWAS for 14 agronomic traits were performed. A total of 913 significant QTNs were found to be associated with the 14 agronomic traits in 3 years and their BLUP values ( Table S2). The number of significant QTNs varied across various traits, ranging from 4 for 2012_IL4 to 37 for BLUP_PH (Figure 1; Table S2), the chromosomal distribution of all identified QTNs revealed that 2H had the maximum number of significant QTNs, which weren't evenly distributed on the genome, and five QTN hotspots on chromosomes 2H, 3H, 6H, and 7H were observed (Figures 1, 2).

The QTNs for Spike and Yield Related Traits
Main spike length (MSL): Three reliable QTNs and one QTN cluster were detected for MSL. Among which, qtnMSL-2H-7 FIGURE 1 | Chromosomal distribution of QTNs identified in this study. The x-axis indicates genomic locations by chromosomal order, and the significant QTNs are plotted against genome location. Each row represents one QTN identified by a different method. The red arrows show the QTN hotspots. located on chromosome 2H: 727,985,438 bp, was repeatedly detected not only in three environments and BLUP values but also by multiple methods to be significantly associated with MSL, and explained phenotypic variation about 2.09-18.16% (LOD score: 3.29-22.92; Tables 1, S3; Figure 2).
Spikelet number of main spike (SMS): qtnSMS-2H-9 referred to Vrs1 (the morphological markers for row number of barley), located at 2H (position: 652,030,802 bp), was identified in all the situations and methods to be significantly associated with SMS (Tables 1, S3), accounting for the largest phenotypic variation (65.07-90.41%). As already known, Vrs1, a gene controlling row number of barley, was validated controlling row number in the DH population derived from a cross between the sixrowed barley cultivar "Huaai 11" and the two-rowed barley cultivar "Huadamai 6." Thus, Vrs1 should control the spike related traits, such as SMS and SLP. Moreover, qtnSMS-2H-8 (2H: 535,680,815 bp) with minor effect was associated with SMS ( Table 1).
Spikelet number per plant (SLP): One reliable QTN qtnSLP-2H-10 (2H: 649,558,019 bp) and one reliable QTN clusters qtncSLP-2H-11 (2H: 652,030,802-653,982,961 bp) close to Vrs1 (2H: 652,030,802 bp) were identified in multiple environments and by multiple methods to be significantly associated with SLP, explaining high proportions of total phenotypic variation, 14.87-55.05 and 8.97-80.30%, respectively. The reliable QTN cluster qtncSLP-4H-1 with a minor effect (2.82-4.14%), mapped on the region 15,498,372-16,168,735 bp of chromosome 4H, was detected in 2 years and the BLUP (Tables 1, S3; Figure 2). FIGURE 2 | Chromosomes location of reliable QTLs for 14 agronomic traits in both previous studies (Ren et al., 2013(Ren et al., , 2014Wang et al., 2016a) and the current studies. The peak positions of previous QTLs were used for mapping, Genetic distance scale in physic position (Mb) is placed at left margin. Green is for the QTLs detected in Ren et al. (2013Ren et al. ( , 2014, black is for the QTLs Wang et al. (2016a), red is for the QTNs and QTNs clusters of current study, the region of QTNs clusters was marked with red on the bar, the cyan is for the candidate genes.
Spike number per plant (SP): Two reliable QTN clusters for SLP were detected. qtncSP-2H-12, located on 2H: 648,821,931-652,030,802 bp overlapping with Vrs1, was detected in three environments and BLUP value and by multiple methods to be associated with SP, accounting for 7.61-52.32% of the phenotypic variation. qtncSP-2H-13, located on chromosome 2H (Position: 662,335,248-663,628,734 bp), was detected in two environments and by multiple methods to be associated with SP, explaining 6.49-11.95% of the phenotypic variation (Tables 1, S3; Figure 2).
Grain number per plant (GP) and grain number per spike (GS): The QTN cluster qtncGP-2H-14 with a high proportion of total phenotypic variation (19.56-37.44%), located at 2H: 649,558,019-650,438,830 bp, was detected in 2010 and BLUP value to be associated with GP. The qtncGP-2H-15 with a high explanation (3.41-34.59%), close to Vrs1 (2H: 652,030,802-652,604,015 bp), was detected in three environments and by multiple methods to be associated with GP (Table 1). One small-effect QTN on 4H (Position: 596,447,744 bp), was detected in 2009 to be significantly associated with GP. Meanwhile, this QTN was also detected for GWP and GWS ( Table 1). One QTN cluster and two QTNs were detected for GS in at least two environments and by multiple methods. The qtncGS-2H-16 (2H: 649,558,019-651,399,477 bp), close to Vrs1, showed a high proportion of phenotypic variation (9.77-71.5%) for GS. The qtnGS-2H-17 referred to Vrs1 (2H: 652,030,802 bp) was significantly associated with GS, explaining a high percentage (49.99-72.25%) of the phenotypic variation.
Grain weight per plant (GWP) and grain weight per spike (GWS): One QTN on 4H (Position: 596,447,744 bp) was found to be associated with GP, GWP and GWS with 3.97-9.17% proportions for the phenotypic variation. The QTN qtnGWS-2H-19, derived from the same associated SNP 2_625783669 (2H: 764,361,924 bp) as GS, was detected for GWS in 2009 and 2010 and by multiple GWAS methods with 0.74-6.25% proportion of phenotypic variation (Tables 1, S3; Figure 2).

Phenotypic Difference Corresponding to QTNs
According to the QTN genotypes, all the DH lines were divided into two different groups to test whether the significant difference  of corresponding phenotypes of the QTN genotypes exist using t-test. Here, six reliable QTNs were used to underlying the phenotypes difference as an example. The details were showed in Figure 3. Among six QTNs, three yield related QTNs (qtnSMS-2H-9 for SMS, qtnGS-2H-17 for GS, and qtnTGW-2H-20 for TGW) had the significant differences of phenotypic averages between their two genotypes in all four environments (Figures 3A-C), and three PH related QTNs (qtnPH-7H-2, qtnIL2-7H-6, and qtnMSL-2H-7) had the significant differences in all four environments (Figures 3D-F), indicating their reliability. It was worth noting that the qtnMSL-2H-7 and qtnSMS-2H-9 was detected in all the four environments, the qtnPH-7H-2, qtnTGW-2H-20, and qtnGS-2H-17 was only detected in three environments, and the qtnIL2-7H-6 was detected only in 2012 and the BLUP ( Table 1).

Identification of Candidate Genes Around Reliable QTNs (or Clusters)
According to the recently released genome sequence of barley (Ibsc, 2016;Beier et al., 2017;Mascher et al., 2017) and the gene annotation information, ten candidate genes for the traits of interest were detected around the reliable QTNs and QTN clusters ( Table 2). The chromosomal distribution of candidate genes was showed in Figure 2. Among which, two genes correspond to the previously reported genes, such as sdw1/denso and Vrs1, while eight candidate genes were new and their functions were derived from the annotated information, which need to be further validated (

DISCUSSION Previously Reported and Novel QTNs Detected With Multi-Locus GWAS Analysis
The comparison between the reliable QTNs (or clusters) for the main agronomic traits and the reliable QTLs in previous studies (Ren et al., 2010(Ren et al., , 2013(Ren et al., , 2014Wang et al., 2016a) were conducted (Table 1; Figure 2). According to the physic positions of associated markers, the reliable QTNs (or clusters) in this study were integrated to the physic map with the reliable QTLs using MapChart 2.32 (Voorrips, 2002 ; Figure 2). Among 39 reliable QTNs (QTN clusters) detected by GWAS, 19 were located on the same regions of QTLs in previous studies (Sameri et al., 2009;Ren et al., 2014;Wang et al., 2016a), while 20 reliable QTNs (QTN clusters) including some minor effect QTNs were novel (Table 1). Totally, 8 of 18 QTNs (QTN clusters) associated with plant height related traits, were same as those in previous studies (Sameri et al., 2009;Ren et al., 2014), and 10 were new. For spike and yield related traits, 11 of 21 QTNs (QTN clusters) were the same as the QTLs in Wang et al. (2016a), the others were new in the current study (Table 1).
The comparison between QTLs in previous studies (Sameri et al., 2009;Ren et al., 2013Ren et al., , 2014Wang et al., 2016a) and in current study indicated that the consistent results should be much more reliable, which is valuable for further gene cloning and molecular marker assistant selection for breeding. Meanwhile, it illustrated that GWAS is feasible and reliable to detect significant associations for complex quantitative traits in DH population. Moreover, some new QTNs with minor effects were detected, suggesting that GWAS should be a good complementary to traditional QTL mapping. The combination of linkage and association analysis should provide the more accurate and powerful approach for reveling the genetic base of complex quantitative traits (Ott et al., 2011).
HORVU7Hr1G008720, was identified at 7H:108,834,015-108,839,990 bp, which was 0.15Mb from the QTN qtnPH-7H-2 for PH, encoding Alpha-mannosidase like AMS1p (Table 2; Figure 2). Alpha-mannosidase is the component of cell wall, involving in cell wall biosynthesis or modification, which participated in the cell growth of internodes with pectinesterase and alpha-xylosidase in plant (Wu and Cao, 2008). Moreover, Alpha-mannosidase is the member of cytoplasm-to-vacuole targeting (Cvt) pathway with AuTophaGy8 (ATG8) gene, and soybean transgenic lines over-expressed GmATG8c showed higher plant height than the wild type (Xia et al., 2012). Therefore, HORVU7Hr1G008720 should be a reliable candidate gene, which affected PH by regulating the cell growth as AMS1p does.
HORVU7Hr1G058360 at 7H: 258,860,422-258,866,854 bp was close to qtnIL2-7H-6 (7H: 258,071,311 bp) for IL2, encoding a S-acyltransferase with DHHC-cysteine-rich domain (Table 2; Figure 2). DHHC-cysteine-rich domain S-acyltransferase proteins are involved in plant development and stress responses in Arabidopsis . AtPAT10 is an S-acyl transferase, which affects the vascular development through controlling the cell division and expansion in Arabidopsis. AtPAT10 mutants are semi-dwarfed, and the reduction of plant height is due to the reduced length of the internodes, which appears to be the result of reduction in both cell number and cell size in these tissues (Qi et al., 2013). Therefore, HORVU7Hr1G058360 is a reliable candidate gene regulating the IL2 as AtPAT10 dose.
HORVU3Hr1G096010 ( Figure 2). Homeobox gene was reported to be involved in the regulation of morphological development in plants, homeobox gene OSH15 affects the architecture of internodes resulting in d6 dwarf plants (Sato et al., 1999). Moreover, the RAD gene in Arabidopsis, encoding small plant-specific single SANT/MYB domain protein, affects the growth and development of Arabidopsis. Overexpression of the RAD gene can repress Arabidopsis growth, resulting in dwarfing and delaying flowering (Baxter et al., 2007;Zhang et al., 2011). Thus, HORVU3Hr1G096010 is a reliable candidate gene regulating plant height as the function of homeobox gene or RAD gene.
For  Figure 2). HvCul4 gene encodes a BLADE-ON-PETIOLE-like (BOP-Like) protein containing BTB/POZ domain, which shares high similarity with Arabidopsis BOP1 and BOP2 (Tavakol et al., 2015;Jost et al., 2016). It was reported that HvCul4 controlled the tiller and leaf pattern in barley (Tavakol et al., 2015;Jost et al., 2016), and Arabidopsis BOP1 and BOP2 acted at boundary regions to regulate axillary development and leaf morphogenesis (Ha et al., 2004). In addition, according to the barley expression database from Barlex (http://barlex.barleysequence.org), HORVU2Hr1G094080.1 showed highest level of expression in developing inflorescences. Therefore, candidate gene HORVU2Hr1G094080 performed the similar function as HvCul4, BOP1, and BOP2 to control the spike number per plant (SP).
For grain number per spike (GS) and grain weight per spike (GWS), marker 2_625783669 (2H: 764,361,924 bp) was detected significantly associated with these 2 traits. HORVU2Hr1G126690 (2H: 764,279,329-764,290,102 bp), encoding a protein with N-acetyltransferase domain, was detected about 83 Kb from the marker 2_625783669 (Table 2; Figure 2). OsSNAT1 encodes N-acetyltransferase1, it was reported that overexpression of T2 homozygous OsSNAT1 in rice increased panicle number and seed weight per plant, while decreased spikelet numbers per panicle under paddy field conditions (Lee and Back, 2017). Moreover, the expression of HORVU2Hr1G126690.4 is much higher in developing inflorescences than in other tissues according to the barley expression database from Barlex (http://barlex.barleysequence. org). Thus, the candidate gene HORVU2Hr1G126690 may affect the GS and GWS through the similar function of OsSNAT1.
Marker 4_497278091 (4H: 596,447,744 bp) was significant associated with GP, GWP, and GWS. HORVU4Hr1G075070 (4H: 596,446,043-596,448,382 bp), encoding Patatin, was detected at 4_497278091 ( Table 2). Overexpression of a patatin-like protein in Camelina sativa (Li et al., 2015) or in Arabidopsis  reduced growth and overall seed production, but increased seed oil content. Therefore, HORVU4Hr1G075070 is a reliable candidate gene which might affect GP, GWP and GWS as the function of patatin-like gene.
Among the above ten candidate genes, two were previously reported, such as sdw1/denso and Vrs1, eight were new, which were derived from the annotated information. Based on the annotations of these candidate genes, homologous genes or proteins with same function or function domain were reported to be regulated the corresponding traits in barley, Arabidopsis and rice. The reliable QTNs and QTN clusters for these traits may be the effect of the candidate genes with similar function as the homologous genes or proteins does. The functions of eight reliable candidate genes need to be further validated. In summary, it is feasible and reliable to use multi-locus GWAS in bi-parental segregation populations.

The New Multi-Locus GWAS for Bi-Parental Segregation Population
Traditionally, segregation populations were used for QTL analysis, and GWAS are commonly used in natural populations. Nowadays, as the development of high-throughput SNP markers and high-throughput phenotypes, GWAS have been widely applied to the genetic analysis for complex traits in familybased populations (such as NAM and MAGIC populations) and proved to be powerful tool for uncovering the basis of key agronomic traits in maize and barley (Tian et al., 2011;Cook et al., 2012;Maurer et al., 2015Maurer et al., , 2016. However, for single segregating population, successful but fewer cases were performed using GWAS (Gao et al., 2015;Henning et al., 2016;Liu et al., 2018). It indicated that GWAS for segregating population are feasible. However, high false positive rate is an obvious problem in the traditional single-locus GWAS using general linear models (GLMs) and mixed linear models (MLMs) Pace et al., 2015). And the P threshold (P = 0.05/n, n is the number of SNPs) leads to missing many significant QTNs, particularly small-effect QTNs (Wang et al., 2016b). Some multi-locus GWAS methodologies, such as mrMLM (Wang et al., 2016b), FASTmrMLM (Zhang and Tamba, 2018), FASTmrEMMA , ISIS EM-BLASSO , pLARmEB , and pKWmEB (Ren et al., 2018) have been developed to remedy the shortcomings mentioned above. These multi-locus GWAS methods have been used to analyze the published data, indicated that these methods constituted effective approaches with high detection power and less stringent criteria (Wang et al., 2016b;Tamba et al., 2017;Zhang et al., 2017;Wen et al., 2018). Totally, five multi-locus GWAS methods were used in our study, which improved the detection power and accuracy of QTNs for interesting traits. Moreover, the QTNs (QTN clusters), which were repeatedly detected in multiple environments and GWAS methods, were selected as reliable QTNs (QTN clusters). This greatly improved the accuracy of the association results and reduced its false positive, and more small-effect QTNs were detected within a certain rate of false positive. In addition, the t-test results of the phenotypic difference corresponding to QTNs, demonstrated that GWAS have the more stringent threshold of significance than t-test, and the advantages of accuracy and false positive controlling. In our study, 39 reliable QTNs and QTN clusters were detected, among which 10 reliable candidate gene were identified. Meanwhile, several new reliable QTNs with small-effect were also detected, which were different from the previous reports (Table 1). There was a limitation to identify candidate genes for all the reliable QTNs and QTN clusters, especially the small-effect ones, based on the imperfect annotation database of barley. However, these results indicated that multi-locus GWAS methods are feasible and reliable for DH population, and good complementary to traditional QTL mapping for the detection of new reliable QTNs even with smalleffect. which will provide more useful information for future works.

CONCLUSIONS
Available online at: In this study, five multi-locus GWAS methods were performed for 14 main agronomic traits in 122 doubled haploid (DH) lines. Thirty-nine reliable QTNs and/or QTN clusters were repeatedly detected in multiple environments and methods 10 candidate genes for the interest traits were detected, 19 QTNs and two genes (sdw1/denso and Vrs1) were previously reported, and eight candidate genes need to be further validated.
The results validated the feasibility and reliability of GWAS in DH population and the good complementary to traditional QTL analysis. All the results will facilitate elucidating genetic basis of agronomic traits and improving marker-assisted selection breeding in barley.

AUTHOR CONTRIBUTIONS
DS and XR conceived and designed the experiments. XH, XR, JW, and LL conducted the experiments and phenotyping measurements. XH and JZ performed the analysis. XH and JZ wrote the paper. GS, XR, and DS modified the manuscript. CL produced the Huaai 11 and Huadamai 6 DH population.
All the authors read and approved the final version of this manuscript.

FUNDING
This project was supported by the earmarked fund for China Agriculture Research System (CARS-5).