Original Research ARTICLE
Association Mapping Reveals Genetic Loci Associated with Important Agronomic Traits in Lentinula edodes, Shiitake Mushroom
- 1Institute of Applied Mycology, Huazhong Agricultural University, Hubei, China
- 2Institute of Bast Fiber Crops, Chinese Academy of Agricultural Sciences, Changsha, China
- 3College of Informatics, Huazhong Agricultural University, Hubei, China
- 4School of Life Sciences, The Chinese University of Hong Kong, Hong Kong, Hong Kong
Association mapping is a robust approach for the detection of quantitative trait loci (QTLs). Here, by genotyping 297 genome-wide molecular markers of 89 Lentinula edodes cultivars in China, the genetic diversity, population structure and genetic loci associated with 11 agronomic traits were examined. A total of 873 alleles were detected in the tested strains with a mean of 2.939 alleles per locus, and the Shannon's information index was 0.734. Population structure analysis revealed two robustly differentiated groups among the Chinese L. edodes cultivars (FST = 0.247). Using the mixed linear model, a total of 43 markers were detected to be significantly associated with four traits. The number of markers associated with traits ranged from 9 to 26, and the phenotypic variations explained by each marker varied from 12.07% to 31.32%. Apart from five previously reported markers, the remaining 38 markers were newly reported here. Twenty-one markers were identified as simultaneously linked to two to four traits, and five markers were associated with the same traits in cultivation tests performed in two consecutive years. The 43 traits-associated markers were related to 97 genes, and 24 of them were related to 10 traits-associated markers detected in both years or identified previously, 13 of which had a >2-fold expression change between the mycelium and primordium stages. Our study has provided candidate markers for marker-assisted selection (MAS) and useful clues for understanding the genetic architecture of agronomic traits in the shiitake mushroom.
Mushrooms have been consumed for their high nutritional and medicinal values for millennia. The production of mushrooms is environmental friendly, due to its capability of converting lignocellulosic waste materials into food, feed and fertilizers (Fan et al., 2006). Due to their importance in industry, agriculture, medicine, and health, the demand for mushrooms is huge and increasing rapidly worldwide (FAOSTAT, http://faostat.fao.org/).
Lentinula edodes (Berk.) Pegler, also called Xianggu or shiitake, is one of the most popular edible mushrooms worldwide. L. edodes was firstly cultivated at least 900–1,000 years ago in China, and is now the second most widely produced mushroom in the world after Agaricus bisporus (Miles and Chang, 2004). L. edodes is rich in minerals, vitamins, essential amino acids, and lentinan (Chang, 1980). It also has immunomodulatory (Xu et al., 2015), anticancer (Nagashima et al., 2013), and antiviral functions (Di Piero et al., 2010). Moreover, L. edodes could prevent environmental impacts caused by the accumulation of forest and agricultural wastes since it secretes hydrolytic and oxidative enzymes that are responsible for the degradation of organic substrates (Silva et al., 2005).
Breeding elite cultivars is important for sustainable development of the modern mushroom industry. The fruiting body of L. edodes is the main target of breeding schemes. Strains with a fast mycelium growth rate (MGR), high precocity, fine morphological characteristics of fruiting body, and high yield are selected as cultivars. The majority of important agronomic traits of L. edodes are quantitative traits controlled by multiple genes or quantitative trait loci (QTLs), which are highly influenced by the environment and show a continuous variation (Santoyo et al., 2008). Dissecting the genetic basis of important agronomic traits is essential for marker-assisted selection (MAS), which could be integrated into conventional breeding schemes in L. edodes and thus enhancing the breeding efficiency.
Linkage mapping was first introduced to understand the genetic basis of quantitative traits since the 1980s (Lander and Botstein, 1989). Since the 2000s, this approach has been successfully used to map QTLs controlling important traits in several mushroom species, such as A. bisporus (Foulongne-Oriol et al., 2012a,b), Pleurotus ostreatus (Larraya et al., 2002, 2003) and P. eryngii (Im et al., 2016). However, linkage mapping requires the construction of segregating populations, which is time-consuming and laborious, thus limiting its application in breeding programs.
As a complementary method of linkage mapping, association mapping based on linkage disequilibrium (LD) is rapidly becoming a powerful strategy for dissecting the genetic architecture of complex traits in plants (Ingvarsson and Street, 2011). Compared to linkage mapping, association mapping has the following advantages: (1) it employs naturally occurring variations in diverse germplasms; (2) it assesses a wide range of alleles rapidly; and (3) it possesses a higher mapping resolution (Yu and Buckler, 2006). Association mapping have been gradually used for genetic dissection of quantitative traits in fungi since the 2010s. However, such studies were mainly limited to the model fungi and fungal pathogens, such as Saccharomyces cerevisiae (Mehmood et al., 2011; Connelly and Akey, 2012), Neurospora crassa (Palma-Guerrero et al., 2013), Fusarium graminearum (Talas et al., 2012), and Heterobasidion annosum (Dalman et al., 2013).
In L. edodes, the genetic dissection of important agronomic traits through linkage mapping and association mapping has recently been performed in our laboratory (Xiang, 2015; Gong et al., 2016). Utilizing linkage mapping, Gong et al. (2016) mapped QTLs underlying seven morphological traits related to characteristics of fruiting body in L. edodes. Six QTL hotspot regions on six linkage groups (LGs) were detected across two segregating testcross populations. A preliminary study of association mapping in a L. edodes population containing 95 wild and one cultivated strains was conducted by using 95 Insertion-Deletion (InDel) and two simple sequence repeat (SSR) markers (Xiang, 2015). Twenty markers were detected to be significantly associated with 10 traits by the general linear model, with phenotypic variations explained by each marker (R2) ranging from 1.29% to 37.60%. However, in these pioneer studies, the phenotypes were only determined in a single environment, and the reliability of the marker-trait association remains to be confirmed.
In this study, 89 L. edodes cultivars collected from major producing areas in China were used to perform association mapping for 11 important agronomic traits by using 297 genome-wide InDel and SSR markers. Cultivation tests were conducted in 2013 and 2014. The objectives of this study were to (1) evaluate the phenotypic and genotypic diversities of L. edodes cultivars in China, (2) investigate QTLs associated with important agronomic traits, and (3) mine candidate genes associated with these traits. The vast number of polymorphic molecular markers identified here could provide a valuable resource for future studies on genetics and breeding in L. edodes. Findings of this study could also benefit MAS and help illustrate the genetic architecture of important agronomic traits in L. edodes.
Materials and Methods
Lentinula edodes Strains
A total of 89 L. edodes cultivars currently cultivated in main producing areas of China were used in this study (Supplementary Table S1). All the tested strains were either provided by research institutes or collected from different mushroom growing farms, and were preserved in the Institute of Applied Mycology, Huazhong Agricultural University (Wuhan, China, 114.35°E, 30.48°N).
Cultivation tests of the 89 L. edodes cultivars were performed in the Huazhong Agricultural University. All strains were allocated in a mushroom house in accordance with the randomized-block design with two blocks, and seven bags of each strain were contained in each replication. Phenotypes were measured from 10 fruiting bodies of each strain. Data were collected from August 2013 to May 2014, and August 2014 to May 2015. Phenotype evaluation of nine fruiting body-related traits was carried out as previously: PD, pileus diameter (mm); PT, pileus thickness (mm); PW, pileus weight (g); SL, stipe length (mm); SD, stipe diameter (mm); SW, stipe weight (g); NF, number of fruiting bodies (per/bag); WF, weight of a single fruiting body (g); Y, total weight of fruiting bodies per bag (g/bag) (Gong et al., 2014b). Two traits about MGRs on both the MYG medium (2% malt extract; 2% glucose; 2% agar; 0.1% peptone; 0.1% yeast extract) and the mixed sawdust (SD) medium (DGR-myg and DGR-sd) were also measured as previously (Gong et al., 2013, 2014a). A brief description of all the 11 traits is summarized in Table 1.
Descriptive statistics, analysis of variance (ANOVA) and correlation analysis were performed using SPSS 17.0 (SPSS Inc., Chicago, IL, USA). For each trait, the phenotypic data were subjected to normality tests (Kolmogorov-Smirnov test). The coefficient of variation (CV) of each trait was calculated by dividing the standard deviation by the mean. ANOVA was first performed with data for each trait in each year independently for the genotypic variation according to the model: Y = μ + G + ε, where μ is the mean value, G is the genotypic effect and ε represents the residual effect. Data from the two independent year experiments were then combined for analyses according to the model: Y = μ + G + E + G × E + ε, where μ is the mean value, G represents the genotypic effect, E is the year effect, G × E is the genotype and year interaction, and ε is the residual effect. Broad-sense heritability (H2) was assessed with the formula H2 = /[ + (/n)] for the two traits relevant to MGR, and H2 = /[ + (/nr) + (/n)] for the nine fruiting body-related traits, where is the genotypic variance, is the error variance, is the variance of genotype and year interaction, n is the number of replicates in the experiment, and r is the number of years.
Phenotypic data in both years were averaged for each fruiting body-related traits. The means were then combined with the phenotypic data of MGRs to construct an Unweighted Pair Group Method with Arithmetic Averaging (UPGMA) dendrogram using the NTSYS-pc (version 2.10e) program (Rohlf, 2000).
SSR and InDel Genotyping
A total of 389 InDel and 104 SSR markers were developed according to our previous description based on the L54A reference genome of L. edodes (GenBank accession number: LOHM00000000) (Gong et al., 2016). Among these markers, two SSR and 166 InDel markers were employed in our previous studies (Xiang, 2015; Gong et al., 2016). PCR and electrophoresis in the SSR and InDel analyses were conducted as previously, with minor modifications (Xiang et al., 2016). Eight strains (Cr04, L135, S602, 430, 9508, Hunong-3, Qin02, and Shandong-1) were randomly selected to screen for the 389 InDel and 104 SSR markers. A total of 379 polymorphic markers, comprising 328 InDels and 51 SSRs, were then selected. Only 297 markers (249 InDels and 48 SSRs) with a minor allele frequency (MAF) of ≥0.05 and missing data ≤5% in the 89 strains were utilized for further analysis. These 297 markers were distributed in 241 scaffolds of the L54A reference genome (Supplemental Table S2).
Genetic Diversity Analysis
Several genetic parameters were assessed to evaluate the genetic diversity of the tested strains. The observed number of alleles (Na), effective number of alleles (Ne), percentage of polymorphic loci (PPL), observed heterozygosity (Ho), expected heterozygosity (He), and Shannon's information index (I) were calculated by using POPGENE 1.32 (Yeh, 1997). Gene diversity (H) and polymorphism information content (PIC) were calculated via PowerMarker V3.25 (Liu and Muse, 2005).
Analysis of Population Structure and Differentiation
Three methods were used to infer the population structure of the 89 L. edodes cultivars. First, an unrooted neighbor-joining (NJ) tree was constructed using Powermarker V3.25. Second, Principal coordinate analysis (PCA) was conducted with GenAlEx version 6.501 (Peakall and Smouse, 2012). All 297 markers were utilized to compute the NJ tree and PCA. Third, STRUCTURE version 2.3 was employed to perform a Bayesian clustering analysis among our tested strains (Pritchard et al., 2010). To reduce the effect of LD between markers on the STRUCTURE results, only 241 markers that were distributed in different scaffolds were utilized in our analysis. The number of subgroups (K) was set from 1 to 15 to survey the optimal range of K using admixture model assumptions with correlated alleles. For each K, seven independent runs were performed, with a burn-in of 10,000 Markov Chain Monte Carlo (MCMC) iterations and a run length of 100,000. Then the ΔK method was used to identify the optimal value of K (Evanno et al., 2005). Strains assigned to the corresponding genetic groups were graphically visualized by Excel 2013. The genetic differentiation within and among populations and the FST values were measured by analysis of molecular variance (AMOVA) using GenAlex version 6.501 (Peakall and Smouse, 2012).
Marker-Trait Association Analysis
Linkage disequilibrium (LD) between loci and association mapping was analyzed using TASSEL 3.0 (Bradbury et al., 2007). The association between molecular markers and the 11 traits was analyzed with the mixed linear model (MLM). To reduce error from the population structure and kinship relationships, Q-matrix and K-matrix were generated by STRUCTURE version 2.3 and TASSEL 3.0, respectively. To control the false positive rate in association mapping, P-values were corrected with FDR (false discovery rate) using the method described by Benjamini and Hochberg (1995), which is widely used in association mapping (Cook et al., 2012; Slovak et al., 2014; Wen et al., 2014; Celik et al., 2016). A FDR-adjusted P-value cut-off of ≤0.05 was used to detect significant marker-trait associations. Annotated genes of the L54A reference genome within a ± 2 kb scope of the significant trait-associated markers were proposed as putative candidate genes for corresponding traits. Functions of these genes were annotated with the Blast2GO package (Conesa et al., 2005).
Because the 89 shiitake cultivars belong to different temperature types, eight strains were unable to fruit or produce an adequate number of fruiting body (>10) in 2 years' cultivation. Therefore, these eight strains were excluded from further analyses about fruiting body-related traits. Furthermore, 17 strains produced <10 fruiting bodies in 2013 or 2014. As a result, only 64 strains that produced >10 fruiting bodies both in 2013 and 2014 were used for two-way analysis of variance.
Frequency distribution of all surveyed agronomic traits was in the form of continuous variation (Figure 1, Supplementary Figure S1), suggesting that the traits were under quantitative and polygenic control. Most traits, except NF and DGR-myg, were distributed normally. After data transformation, the phenotypic data of NF and DGR-myg were still not fit the normal distribution. Thus, we used the original data of NF and DGR-myg for the statistical analysis.
Figure 1. Histograms showing the frequency distribution of 11 agronomic traits in 2013. The y-axis denotes the number of strains, whereas the x-axis indicates value range of traits.
For each trait, phenotypic variation among strains was confirmed by its range, mean, standard deviation, and CV (Table 1). The CV among strains ranged from 5.08% (DGR-sd) to 88.95% (NF) in 2014. For the nine fruiting body-related traits, the CV among strains ranged from 19.54% (PD) to 78.57% (NF) in 2013, and from 14.96% (PD) to 88.95% (NF) in 2014. In general, the CV of fruiting body-related traits were higher than those in mycelium-related traits. For the nine fruiting body-related traits, significant correlations (P < 0.05) were detected between experiments conducted in the 2 years (Table 1). As revealed by ANOVA analysis, genotype had highly significant effect on DGR-myg and DGR-sd (data not shown). Also, the effects of genotype, year, and interaction of genotype × year on the nine fruiting body-related traits were all significant (P < 0.05) (Supplementary Table S3). Values of the broad-sense heritability (H2) of the 11 agronomic traits were high, and varied from 81.30% (Y) to 98.82% (DGR-myg) (Table 1).
Pearson's correlation was used to uncover the relationship of the nine fruiting body-related traits in both years. There were significant positive correlations (P < 0.01) between six traits related to characteristics of single fruiting body, i.e., PD, PT, PW, SD, SL, and SW. Yield (Y) was significantly positively correlated with the yield-component trait NF but significantly negatively correlated with another yield-component trait WF (P < 0.01); NF was significantly negatively correlated with WF (P < 0.01; Table 2).
Table 2. A matrix of Pearson correlation coefficients (r) for the 11 agronomic traits in Chinese Lentinula edodes cultivars.
Sixty-four strains that we can get at least 10 fruiting bodies for each year were employed in phenotype-based UPGMA clustering. These strains were divided into two groups at the Euclidean distance of 118.19 (Supplementary Figure S2). Strains belonging to different groups in the NJ tree clustered together in the UPGMA dendrogram.
A total of 873 alleles were detected from the 297 loci in all 89 strains, and the number of alleles in each locus varied from two to seven with a mean of 2.939 (Supplementary Table S2). The PIC value varied from 0.022 (S676_SSR2) to 0.731 (S95_ID5) with an average of 0.381. In all 89 strains, the average number of alleles from the 48 SSR markers was 3.229, higher than that from the 249 InDel markers (2.884). However, the PIC value of SSRs (0.379) was comparable to that of InDels (0.381), suggesting a similar genetic variation of both types of markers in the Chinese shiitake cultivars. From all 89 strains, the mean values of Ne, He, I, PPL, H, and PIC were 1.955, 0.456, 0.734, 100%, 0.454, and 0.381, respectively, indicating relatively low genetic variations among Chinese L. edodes cultivars. As for the two groups defined by the NJ tree, all the genetic parameters showed that the genetic diversity in Group A was lower than that in Group B (Table 3). The lower genetic diversity in Group A could be due to the smaller sample size.
Population Structure and Linkage Disequilibrium
According to the genotyping data, the 89 strains represent 89 unique genotypes, and therefore are not clones. In the NJ tree of L. edodes, all strains except Xiangjiu clustered into two distinct groups (Figure 2). Group A consisted of 21 strains and Group B contained 67 strains. PCA also identified two groups congruent with those in the NJ tree (Supplementary Figure S3). The percentages of variation explained by the first 3 axes were 32.7%, 14.1%, and 9.8%.
Figure 2. A neighbor-joining tree of 89 Lentinula edodes cultivars. The ancestors of the strains in the inferred groups are represented by different colors.
Analysis of molecular variance (AMOVA) results suggested that the majority of genetic variation was included within populations (75.33%) (Table 4). The overall FST value across all the strains except Xiangjiu was 0.247, suggesting a great differentiation among L. edodes cultivars in China.
Table 4. Analysis of molecular variance (AMOVA) among and within populations of Lentinula edodes cultivars in China.
Model-based STRUCTURE was also utilized to investigate the population structure of the 89 strains. In the analysis of ΔK, a clear maximum was detected for K = 2 (ΔK = 3015) (Figure 3A). Therefore, two groups were identified in the collection of 89 L. edodes cultivars in China (Figure 3B), which agreed with the results of NJ tree and PCA. None of the 89 strains was assigned exclusively to one group or the other, and all strains shared mixed ancestries from the two groups. These demonstrate that the Chinese shiitake cultivars were genetically closely related.
Figure 3. Results of STRUCTURE analysis. (A) Estimation of the number of populations for K ranging from 2 to 15 by ΔK values; (B) Classification of 89 L. edodes cultivars into two genetic groups. The distribution of the strains assigned to different groups is indicated by the color code (Group A: red, Group B: blue). The y-axis quantifies the cluster membership, and the x-axis lists the different strains. Strains from the different groups defined in the NJ tree are marked in different symbols: ♦, Group A; ■, Group B; ▴, Xiangjiu, excluded from the two groups in the NJ tree.
A total of 19,122 (43.50%) InDel and SSR marker pairs displayed significant LD among all 89 strains (P ≤ 0.001). The r2 values among these marker pairs varied from 0.128 to 1, with an average of 0.316 (Supplementary Table S4). At the highly significant threshold of r2 ≥ 0.2, 30.43% (13,378) of the marker pairs remained in LD. In this study, owing to the fact that only 73 of these 297 makers were used to construct linkage map in our previous study (Gong et al., 2016), the genome-wide LD decay along with the increase of genetic distance were not detected (Supplementary Table S5). The averaged r2 was 0.419 in the region from 0 to 20 cM, and decreased to 0.384 in the region from 20 to 40 cM, then dropped to 0.361 in the region from 40 to 60 cM. However, the averaged r2 increased to 0.381 in the region of >60 cM. Extremely high LD levels were observed across all 89 strains. The average r2 of the unlinked marker pairs was 0.327, which was much lower than that in the linked marker pairs (0.400) (Supplementary Table S5).
Association between Traits and Molecular Markers
The numbers of strains used for association mapping are 75 (in 2013) and 69 (in 2014). Marker-trait associations for the 11 traits were evaluated using the Q + K model in Tassel 3.0 (Table 5). A total of 78 associations covering four traits (PD, PW, SW, and WF) were detected at FDR-adjusted P ≤ 0.05, and the phenotypic variation explained by each marker ranged from 12.07 to 31.32% (Table 5). In 2013 and 2014, the number of significant marker-trait associations were 17 and 61, respectively (Table 5). Eleven markers were associated with two traits (PD and PW) in 2013, and 41 markers were associated with four traits (PD, PW, SW, and WF) in 2014. The number of markers identified to be associated with the four traits varied from 9 (PD) to 26 (PW) with an average of 17.75.
Table 5. Associations between agronomic traits and markers in the Chinese Lentinula edodes cultivars.
The association of a single marker with multiple traits could be the result of pleiotropy. Twenty-one molecular markers were associated with two to four traits (Table 6), five of which (S127_ID1, S328_ID5, S278_ID10, S278_ID41, and S704_inID1) were identified to be associated with the same trait in both years. For instance, S278_ID10 was found to be associated with PD and PW in 2013, and associated with WF, PD, PW, SW in 2014. In particular, associations between S278_ID10 and PD and PW were detected in both years.
Ninety-seven annotated genes were detected within a ± 2 kb scope of the significant trait-associated markers in the L. edodes reference genome, 31 of which contained trait-associated markers. Results from Blast2GO indicated that the 97 genes were involved in a wide range of molecular functions (level 3), including organic cyclic compound binding, small molecule binding, ion binding, transferase activity, carbohydrate derivative binding, hydrolase activity, heterocyclic compound binding, oxidoreductase activity, and protein binding (Supplementary Figure S4).
Genetic Variation of the 11 Agronomic Traits
In this study, ANOVA revealed the significant influences of the genotypes for all the surveyed traits. Heritability was also high for all the 11 traits. These observations indicated that the majority of the investigated traits were highly inheritable, and genotypic effect was the main factor to generate phenotypic variations.
The CV of the two mycelium growth-related traits were smaller than those of the nine fruiting body-related traits. This is because the two traits were determined in incubators and therefore not affected by the changing environmental conditions under which the nine traits were measured. Indeed, ANOVA revealed significant differences of the fruiting body-related traits between the 2 years, suggesting a strong effect of environmental factors (i.e., year) on these traits.
Here, we observed extensive significant correlations between the nine fruiting body-related traits of L. edodes cultivars, in agreement with previous results detected in two segregating populations and one natural population (Gong et al., 2014b). Yield and yield-component traits of L. edodes were found to exhibit the triangular relationship as displayed in our recent report (Gong et al., 2014b) and in A. bisporus (Foulongne-Oriol et al., 2012a), i.e., yield was positively correlated with NF but negatively correlated with WF; and NF was negatively correlated with WF. For nine fruiting body-related traits and two mycelium growth-related traits, no obvious correlation was observed in 2013, which indicated that the growth of mycelium and development of fruiting bodies may be independently controlled.
A total of 21 molecular markers were identified to be associated with two to four traits by association mapping. For instance, S106_inID1 was associated with four traits (WF, PD, PW, and SW), and S560_ID1 associated with three traits (WF, PW, and PD), partly illustrating the genetic basis of phenotypic correlation between these traits. The two major reasons for trait correlations are pleiotropy and close linkage between QTLs controlling different traits (Mackay et al., 2009; Chen and Lübberstedt, 2010). Our recent work on QTL mapping in two segregating populations also suggested that the co-localization of QTLs underlining different traits may be the genetic basis for phenotypic correlation of fruiting body-related traits in L. edodes (Gong et al., 2016). Combining evidences from both association mapping and our recent results from QTL mapping, we postulate that the genetic basis of phenotypic correlation in L. edodes is the tight linkage of QTLs and pleiotropy.
Multigenic effects were also observed in this study (Table 5). For instance, nine markers, including S278_ID10, S278_ID36, S278_ID41, and S328_ID5, were associated with PD. Multigenic effects suggested that these traits in L. edodes were complex quantitative traits that were affected by polygenes.
Understanding the genetic diversity and genetic basis underlying important agronomic traits could improve breeding schemes of L. edodes. In general, the processes of domestication and breeding have a strong impact on the genetic diversity of cultivated species (Font i Forcada et al., 2015). Here, the values of Shannon's information index (I) and polymorphism information content (PIC) revealed by InDel and SSR markers were 0.734 and 0.381, respectively. In a wild population containing 88 Chinese L. edodes strains, the I and PIC values were 0.836 and 0.395, respectively (Xiang et al., 2016). In another study, the PIC value was 0.53 in 89 L. edodes strains from East Asia (Kim et al., 2009). Genetic variation of Chinese L. edodes cultivars is low and was postulated to be derived from a limited number of elite strains (Chiu et al., 1996). Therefore, the wild strains should be introduced into the breeding schemes to diversify the genetic basis of shiitake cultivars in China.
Population Structure and Linkage Disequilibrium
Detailed knowledge on population structure is important to control spurious associations between phenotypes and genotypes in association mapping (Pritchard et al., 2000). Model-based analysis of population structure could provide necessary information in association mapping. Here, population structure analyses based on three methods demonstrated that the Chinese L. edodes cultivars could be divided into two unique groups, with Xiangjiu being the sole exception. This strain was proven to be distinct from other L. edodes cultivars in previous clustering analyses (Fu et al., 2010; Liu et al., 2015). Using strains different from the current study, the L. edodes cultivars in China was also separated into two main groups (Zhang et al., 2007; Fu et al., 2010; Liu et al., 2012). Therefore, it is reasonable to speculate that L. edodes cultivars in China contained two different gene pools, which possibly resulted from domestication and breeding. Great genetic differentiation existed between the two groups as indicated by a FST value of 0.247, which is comparable to that in the Chinese wild L. edodes population (FST = 0.252) (Xiang et al., 2016).
A high level of LD among marker pairs was observed in this study. As mentioned before, the narrow genetic base of Chinese shiitake cultivars revealed here and in previous studies (Chiu et al., 1996; Fu et al., 2010) might be one of the factors that could explain the high level of pairwise LD. Also, a small number of tested strains and molecular markers may cause bias of LD estimates. Furthermore, the population structure contributes to increasing LD level. Population structure could create unexpected LD between unlinked loci across the genome (Yan et al., 2011). The mixing of individuals belonging to different subpopulations with different allele frequencies creates LD, when these subpopulations are admixed to construct a panel of lines for association mapping. Significant LD between unlinked loci results in false-positive associations between a marker and a trait (Soto-Cerda and Cloutier, 2012). In this study, small structured population with narrow genetic base may be the major factors that causes the high level LD and then led to spurious associations between marker alleles and the phenotypes. Therefore, the further association studies require the careful choice of germplasm, as well as a larger number of markers and strains.
Association mapping has been successfully utilized in crop species, such as maize, cotton, wheat, and rice (Abdurakhmonov and Abdukarimov, 2008). Due to their edible and medicinal values, mushrooms have been consumed by humans for a long time. However, linkage and association mapping in mushroom species are still in their infancy, and information is limited to identify QTLs controlling agronomic traits in mushroom species. Here, we detected 78 marker-trait associations covering 43 molecular markers and four traits. Marker-trait associations detected by this method could provide valuable information for MAS in breeding schemes of L. edodes.
Among the 297 markers used here, 73 markers were the same as those used by Gong et al. (2016), while 47 markers were employed by Xiang (2015). Five markers were found to be associated with the same traits as in previous reports (Xiang, 2015; Gong et al., 2016) (Table 5). Gong et al. (2016) revealed that marker S48_ID1 was located in a QTLs-hotspot region in LG2 that was related to PD, PT, PW, SL, SD, SW, and WF. In this study, S48_ID1 was also identified to be associated with WF. S560_ID1 was significantly associated with WF, PW and PD, consistent with previous findings (Gong et al., 2016). S346_ID1 was found to be associated with SW and lie in a QTL hotspot region in LG4 related to SW, SL, and SD (Gong et al., 2016). Three of the six hotspot regions previously identified by Gong et al. (2016) were confirmed by this association mapping, thus suggesting that there are reliable regions harboring QTLs related to fruiting body in L. edodes. Moreover, S255_ID1 was found to be associated with PW and WF, and S163_E1 associated with PW, in agreement with a previous report (Xiang, 2015).
Five markers were identified to be associated with the same trait in both years. They are S127_ID1 and S328_ID5 for PW, S278_ID10 and S704_inID1 for PW and PD, and S278_ID41 for PD. It is worth mentioning that S278_ID41 is located in a <1 kb position to S278-R/F that lay in a QTL hotspot in MG4 special for SL, PD, and PW (Gong et al., 2016). Hence, S278_ID41 is also located in this QTL hotspot in MG4.
The foregoing 10 markers were verified by previous studies or detected in both years, suggesting that association mapping in shiitake cultivars can effectively identify major QTLs underlying agronomic traits. Apart from the five markers that were previously detected, 38 markers related to four traits in L. edodes were newly reported here, demonstrating that association mapping could detect more trait-associated markers than conventional linkage mapping. Moreover, unlike linkage mapping, which is based on a bi-parental population, association mapping investigates genetic variations in a natural population, and thus can evaluate many alleles simultaneously (Abdurakhmonov and Abdukarimov, 2008). Therefore, this study confirmed the feasibility and reliability of association mapping in L. edodes.
The 10 aforementioned markers resided on or were close to 24 annotated genes of the L. edodes reference genome (Supplementary Table S6). These genes could be potential candidate genes related to agronomic traits. By checking the RNA-seq expression levels of the 24 candidate genes in shiitake strain L54 (unpublished data), 13 of them (54.17%) were found to be significantly up-regulated or down-regulated during the transition from mycelium to primordium, with fold changes >2 (Supplementary Figure S5). Hence, it is reasonable to speculate that most of these candidate genes are involved in the development of fruiting body.
Despite the importance of edible mushrooms, research on their breeding and production are still very limited as compared to other crops, which may be partly due to the lack of knowledge of their genetics and breeding system (Chakravarty, 2011). Linkage mapping is usually conducted in purpose-created segregating populations, such as progeny of selected parents. However, the resolution of linkage mapping is hampered by the limited number of recombination in the segregating population. Moreover, there are some particularities in linkage mapping in edible mushrooms. Construction of genetic map is performed using a haploid progeny, whereas the phenotypic evaluation of some traits is only possible at the dikaryotic stage, after crossing haploid progenies with the compatible monokaryotic tester (Foulongne-Oriol, 2012). For each locus, these mushrooms have the same allele from the tester line in one nucleus and contain the allele from one of the parents in the other nucleus (Gao et al., 2015). The identified QTLs reflect the allelic substitution effect of the segregating allele with their interactions with the constant allele from the tester nucleus, and such constraints may lead to inconsistencies in QTL detection (Foulongne-Oriol, 2012). Alternatively, association mapping compares favorably to linkage mapping for the dissection of natural variation by using diverse germplasms, such as those derived from the wild populations, germplasm collections or subsets of breeding germplasm (Zhu et al., 2008). Because many generations have passed and more recombinations occurred, the resolution of association mapping is considerably higher than that in simple bi-parental populations (Rafalski, 2010). However, the complex population structure is one of the sources of false positives in association mapping. For instance, the complex and heterogeneous population structure of a S. cerevisiae population was reported to lead to a high type I error rate in association mapping (Connelly and Akey, 2012). A discernable population structure was also reported in shiitake cultivars here. Moreover, the allele frequency distribution has considerable impact on the detection power of association mapping. Here, only the markers with a (MAF) ≥0.05 were utilized for association analysis.
Overall, the combination of linkage mapping and association mapping could be a more powerful strategy for dissecting the genetic architectures of quantitative traits in edible mushrooms. To begin with, the genomic regions underlying quantitative traits of interest could be defined by linkage mapping. Then, based on the results of linkage mapping, fine mapping could be performed by candidate gene-based association analysis to identify the loci or genes for the traits of interest. In addition, the utilization of nested association mapping populations (Yu et al., 2008) could be a promising approach for genetic dissection of quantitative traits in edible mushrooms.
In summary, we reported here the genetic diversity, population structure and association mapping of agronomic traits in a Chinese L. edodes population containing 89 cultivars by using 297 genome-wide markers. A narrow genetic base with a discernable population structure was observed in the Chinese shiitake cultivars. In association mapping, a total of 43 markers were detected to be significantly associated with four traits. Five of these marker-trait associations were verified by previous studies and another five of them were significantly detected in cultivation tests performed in two consecutive years. Our results have highlighted the significant potential of LD-based association mapping of complex agronomic traits in shiitake with consideration of the population structure. Associations identified here could provide insights into the genetic architecture of important agronomic traits, thus paving a way toward implementation of MAS in L. edodes.
YX and YB conceived and designed the experiments; CL and LZ performed the experiments; WN and HK developed the molecular markers; CL, WG, and ZY analyzed the phenotypic and genotypic data; YX, WG, and MC wrote the paper.
Conflict of Interest Statement
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
This work was financially supported by the National Natural Science Foundation of China (Grant No. 31372117), the National Key Technology Support Program in the 12th 5-Year Plan of China (Grant No. 2013BAD16B02), and Research Grants Council of Hong Kong (Grant No. CUHK466312).
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/article/10.3389/fmicb.2017.00237/full#supplementary-material
Abdurakhmonov, I. Y., and Abdukarimov, A. (2008). Application of association mapping to understanding the genetic diversity of plant germplasm resources. Int. J. Plant Genomics. 2008:574927. doi: 10.1155/2008/574927
Bradbury, P. J., Zhang, Z., Kroon, D. E., Casstevens, T. M., Ramdoss, Y., and Buckler, E. S. (2007). TASSEL: software for association mapping of complex traits in diverse samples. Bioinformatics 23, 2633–2635. doi: 10.1093/bioinformatics/btm308
Celik, I., Camci, H., Kose, A., Kosar, F. C., Doganlar, S., and Frary, A. (2016). Molecular genetic diversity and association mapping of morphine content and agronomic traits in Turkish opium poppy (Papaver somniferum) germplasm. Mol. Breed. 36, 1–13. doi: 10.1007/s11032-016-0469-8
Chakravarty, B. (2011). Trends in mushroom cultivation and breeding. Aust. J. Agric. Eng. 2, 102–109. Available online at: https://search.informit.com.au/documentSummary;dn=683694909801019;res=IELENG
Chiu, S. W., Ma, A. M., Lin, F. C., and Moore, D. (1996). Genetic homogeneity of cultivated strains of shiitake (Lentinula edodes) used in China as revealed by the polymerase chain reaction. Mycol. Res. 100, 1393–1399. doi: 10.1016/S0953-7562(96)80069-4
Conesa, A., Götz, S., García-Gómez, J. M., Terol, J., Talón, M., and Robles, M. (2005). Blast2GO: a universal tool for annotation, visualization and analysis in functional genomics research. Bioinformatics 21, 3674–3676. doi: 10.1093/bioinformatics/bti610
Cook, J. P., McMullen, M. D., Holland, J. B., Tian, F., Bradbury, P., Ross-Ibarra, J., et al. (2012). Genetic architecture of maize kernel composition in the nested association mapping and inbred association panels. Plant Physiol. 158, 824–834. doi: 10.1104/pp.111.185033
Dalman, K., Himmelstrand, K., Olson, Å., Lind, M., Brandström-Durling, M., and Stenlid, J. (2013). A genome-wide association study identifies genomic regions for virulence in the non-model organism Heterobasidion annosum s.s. PLoS ONE 8:e53525. doi: 10.1371/journal.pone.0053525
Di Piero, R. M., de Novaes, Q. S., and Pascholati, S. F. (2010). Effect of Agaricus brasiliensis and Lentinula edodes mushrooms on the infection of passionflower with Cowpea aphid-borne mosaic virus. Braz. Arch. Biol. Technol. 53, 269–278. doi: 10.1590/S1516-89132010000200004
Evanno, G., Regnaut, S., and Goudet, J. (2005). Detecting the number of clusters of individuals using the software STRUCTURE: a simulation study. Mol. Ecol. 14, 2611–2620. doi: 10.1111/j.1365-294X.2005.02553.x
Fan, L., Pan, H., Soccol, A. T., Pandey, A., and Soccol, C. R. (2006). Advances in mushroom research in the last decade. Food Sci. Biotechnol. 44, 303–311. Available online at: http://mychagaworld.com/pdf/Advances%20in%20Mushroom%20Research%20in%20the%20Last%20Decade.pdf
Font i Forcada, C., Oraguzie, N., Reyes-Chin-Wo, S., Espiau, M. T., and Fernández i Martí, A. (2015). Identification of genetic loci associated with quality traits in almond via association mapping. PLoS ONE 10:e0127656. doi: 10.1371/journal.pone.0127656
Foulongne-Oriol, M., Rodier, A., Rousseau, T., and Savoie, J. M. (2012a). QTL mapping of yield-related components and oligogenic control of the cap colour in the button mushroom Agaricus bisporus. Appl. Environ. Microbiol. 78, 2422–2434. doi: 10.1128/AEM.07516-11
Foulongne-Oriol, M., Rodier, A., and Savoie, J. M. (2012b). Relationship between yield components and partial resistance to Lecanicillium fungicola in the button mushroom, Agaricus bisporus, assessed by quantitative trait locus mapping. Appl. Environ. Microbiol. 78, 2435–2442. doi: 10.1128/AEM.07554-11
Fu, L. Z., Zhang, H. Y., Wu, X. Q., Li, H. B., Wei, H. L., Wu, Q. Q., et al. (2010). Evaluation of genetic diversity in Lentinula edodes strains using RAPD, ISSR and SRAP markers. World J. Microb. Biotechnol. 26, 709–716. doi: 10.1007/s11274-009-0227-8
Gao, W., Weijn, A., Baars, J. J., Mes, J. J., Visser, R. G., and Sonnenberg, A. S. (2015). Quantitative trait locus mapping for bruising sensitivity and cap color of Agaricus bisporus (button mushrooms). Fungal Genet. Biol. 77, 69–81. doi: 10.1016/j.fgb.2015.04.003
Gong, W. B., Bian, Y. B., Xu, R., and Xiao, Y. (2013). Genetic relationship between growth rate, biomass, germination period and mating type of monokaryons in Lentinula edodes. J. Pure Appl. Microbiol. 7, 863–869.
Gong, W. B., Li, L., Zhou, Y., Bian, Y. B., Kwan, H. S., Cheung, M. K., et al. (2016). Genetic dissection of fruiting body-related traits using quantitative trait loci mapping in Lentinula edodes. Appl. Microbiol. Biotechnol. 100, 5437–5452. doi: 10.1007/s00253-016-7347-5
Gong, W. B., Liu, W., Lu, Y. Y., Bian, Y. B., Zhou, Y., Kwan, H. S., et al. (2014a). Constructing a new integrated genetic linkage map and mapping quantitative trait loci for vegetative mycelium growth rate in Lentinula edodes. Fungal Biol. 118, 295–308. doi: 10.1016/j.funbio.2014.01.001
Gong, W. B., Xu, R., Xiao, Y., Zhou, Y., and Bian, Y. B. (2014b). Phenotypic evaluation and analysis of important agronomic traits in the hybrid and natural populations of Lentinula edodes. Sci. Hortic. 179, 271–276. doi: 10.1016/j.scienta.2014.09.044
Im, C. H., Park, Y. H., Hammel, K. E., Park, B., Kwon, S. W., Ryu, H., et al. (2016). Construction of a genetic linkage map and analysis of quantitative trait loci associated with the agronomically important traits of Pleurotus eryngii. Fungal Genet. Biol. 92, 50–64. doi: 10.1016/j.fgb.2016.05.002
Kim, K. H., Kim, Y. Y., Ka, K. H., Lee, H. S., Bak, W. C., Jeong, S. J., et al. (2009). Microsatellite markers for population-genetic studies of shiitake (Lentinula edodes) strains. Genes Genomics 31, 403–411. doi: 10.1007/BF03191853
Larraya, L. M., Idareta, E., Arana, D., Ritter, E., Pisabarro, A. G., and Ramírez, L. (2002). Quantitative trait loci controlling vegetative growth rate in the edible basidiomycete Pleurotus ostreatus. Appl. Environ. Microbiol. 68, 1109–1114. doi: 10.1128/AEM.68.3.1109-1114.2002
Larraya, L. M., Mikel, A., Pisabarro, A. G., and Ramírez, L. (2003). Mapping of genomic regions (quantitative trait loci) controlling production and quality in industrial cultures of the edible basidiomycete Pleurotus ostreatus. Appl. Environ. Microbiol. 69, 3617–3625. doi: 10.1128/AEM.69.6.3617-3625.2003
Liu, J., Wang, Z. R., Li, C., Bian, Y. B., and Xiao, Y. (2015). Evaluating genetic diversity and constructing core collections of Chinese Lentinula edodes cultivars using ISSR and SRAP markers. J. Basic. Microbiol. 55, 749–760. doi: 10.1002/jobm.201400774
Liu, J. Y., Ying, Z. H., Liu, F., Liu, X. R., and Xie, B. G. (2012). Evaluation of the use of SCAR markers for screening genetic diversity of Lentinula edodes strains. Curr. Microbiol. 64, 317–325. doi: 10.1007/s00284-011-0069-0
Mehmood, T., Martens, H., Saebø, S., Warringer, J., and Snipen, L. (2011). Mining for genotype-phenotype relations in Saccharomyces using partial least squares. BMC Bioinformatics 12:318. doi: 10.1186/1471-2105-12-318
Nagashima, Y., Maeda, N., Yamamoto, S., Yoshino, S., and Oka, M. (2013). Evaluation of host quality of life and immune function in breast cancer patients treated with combination of adjuvant chemotherapy and oral administration of Lentinula edodes mycelia extract. Onco Targets Ther. 6, 853. doi: 10.2147/OTT.S44169
Palma-Guerrero, J., Hall, C. R., Kowbel, D., Welch, J., Taylor, J. W., Brem, R. B., et al. (2013). Genome wide association identifies novel loci involved in fungal communication. PLoS Genet. 9:e1003669. doi: 10.1371/journal.pgen.1003669
Peakall, R., and Smouse, P. E. (2012). GenAlEx 6.5: genetic analysis in Excel. Population genetic software for teaching and research-an update. Bioinformatics 28, 2537–2539. doi: 10.1093/bioinformatics/bts460
Pritchard, J. K., Stephens, M., and Donnelly, P. (2000). Inference of population structure using multilocus genotype data. Genetics 155, 945–959. Available online at: http://www.genetics.org/content/155/2/945.short
Santoyo, F., González, A. E., Terrón, M. C., Ramírez, L., and Pisabarro, A. G. (2008). Quantitative linkage mapping of lignin-degrading enzymatic activities in Pleurotus ostreatus. Enzyme Microb. Technol. 43, 137–143. doi: 10.1016/j.enzmictec.2007.11.007
Silva, E. M., Machuca, A., and Milagres, A. M. (2005). Effect of cereal brans on Lentinula edodes growth and enzyme activities during cultivation on forestry waste. Lett. Appl. Microbiol. 40, 283–288. doi: 10.1111/j.1472-765X.2005.01669.x
Slovak, R., Göschl, C., Su, X., Shimotani, K., Shiina, T., and Busch, W. (2014). A scalable open-source pipeline for large-scale root phenotyping of arabidopsis. Plant Cell 26, 2390–2403. doi: 10.1105/tpc.114.124032
Talas, F., Wüerschum, T., Reif, J. C., Parzies, H. K., and Miedaner, T. (2012). Association of single nucleotide polymorphic sites in candidate genes with aggressiveness and deoxynivalenol production in Fusarium graminearum causing wheat head blight. BMC Genet. 13:14. doi: 10.1186/1471-2156-13-14
Wen, Z., Tan, R., Yuan, J., Bales, C., Du, W., Zhang, S., et al. (2014). Genome-wide association mapping of quantitative resistance to sudden death syndrome in soybean. BMC Genomics 15:809. doi: 10.1186/1471-2164-15-809
Xiang, X., Li, C., Li, L., Bian, Y., Kwan, H. S., Nong, W., et al. (2016). Genetic diversity and population structure of Chinese Lentinula edodes revealed by InDel and SSR markers. Mycol. Prog. 15, 1–13. doi: 10.1007/s11557-016-1183-y
Xu, X., Yang, J., Ning, Z., and Zhang, X. (2015). Lentinula edodes-derived polysaccharide rejuvenates mice in terms of immune responses and gut microbiota. Food Funct. 6, 2653–2663. doi: 10.1039/C5FO00689A
Yu, J., Holland, J. B., McMullen, M. D., and Buckler, E. S. (2008). Genetic design and statistical power of nested association mapping in maize. Genetics 178, 539–551. doi: 10.1534/genetics.107.074245
Zhang, R., Huang, C., Zheng, S., Zhang, J., Ng, T. B., Jiang, R., et al. (2007). Strain-typing of Lentinula edodes in China with inter simple sequence repeat markers. Appl. Microbiol. Biotechnol. 74, 140–145. doi: 10.1007/s00253-006-0628-7
Keywords: shiitake mushroom, microsatellite, InDel, quantitative trait, association analysis
Citation: Li C, Gong W, Zhang L, Yang Z, Nong W, Bian Y, Kwan H-S, Cheung M-K and Xiao Y (2017) Association Mapping Reveals Genetic Loci Associated with Important Agronomic Traits in Lentinula edodes, Shiitake Mushroom. Front. Microbiol. 8:237. doi: 10.3389/fmicb.2017.00237
Received: 26 August 2016; Accepted: 03 February 2017;
Published: 17 February 2017.
Edited by:Martin G. Klotz, Queens College (CUNY), USA
Reviewed by:Chandra Nayak, University of Mysore, India
Marie Foulongne-Oriol, Institut National de la Recherche Agronomique, France
Copyright © 2017 Li, Gong, Zhang, Yang, Nong, Bian, Kwan, Cheung and Xiao. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) or licensor are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
†These authors have contributed equally to this work.