Genetic control of arbuscular mycorrhizal colonization by Rhizophagus intraradices in Helianthus annuus (L.)

Plant symbiosis with arbuscular mycorrhizal (AM) fungi provides many benefits, including increased nutrient uptake, drought tolerance, and belowground pathogen resistance. To develop a better understanding of the genetic architecture of mycorrhizal symbiosis, we conducted a genome-wide association study (GWAS) of this plant-fungal interaction in cultivated sunflower. A diversity panel of cultivated sunflower (Helianthus annuus L.) was phenotyped for root colonization under inoculation with the AM fungus Rhizophagus intraradices. Using a mixed linear model approach with a high-density genetic map, we identified genomic regions that are likely associated with R. intraradices colonization in sunflower. Additionally, we used a set of twelve diverse lines to assess the effect that inoculation with R. intraradices has on dried shoot biomass and macronutrient uptake. Colonization among lines in the mapping panel ranged from 0–70% and was not correlated with mycorrhizal growth response, shoot phosphorus response, or shoot potassium response among the Core 12 lines. Association mapping yielded three single-nucleotide polymorphisms (SNPs) that were significantly associated with R. intraradices colonization. This is the first study to use GWAS to identify genomic regions associated with AM colonization in an Asterid eudicot species. Three genes of interest identified from the regions containing these SNPs are likely related to plant defense.


Introduction
Arbuscular mycorrhizal (AM) symbiosis is formed between plants and members of the Glomeromycota, a monophyletic group of obligate symbiotic fungi (Schüßler et al. 2001;Delaux 2017). Phylogenomic analyses indicate that the Glomeromycota have formed symbioses with plants since the early diversification of embryophytes and likely played an important role in the evolution of land plants (Remy et al. 1994;Delaux et al. 2015;Berbee et al. 2017). Presently, an estimated 65-72% of flowering plant species form symbiotic relationships with AM fungi, as numerous plant lineages have independently lost genes necessary for AM symbiosis (Wang and Qiu 2006;Brundrett and Tedersoo 2018). AM fungi colonize the roots of host plants, forming structures such as intraradical hyphae, vesicles, and arbuscules (Brundrett 2004;Smith and Read 2008). This symbiosis is characterized by the exchange of photosynthetically fixed carbon products for nutrients at arbuscules, the sites of nutrient exchange (Bonfante and Genre 2010). AM fungi are influential in modulating ecosystem structure and productivity and contributing to the coexistence and maintenance of species through plant-soil feedbacks (Hartnett and Wilson 2002;Bever 2003;Rillig 2004).
In recent years, forward and reverse genetics have identified genes that are necessary for the recognition of symbiotic partners, transmission of signal, and transcription of genes necessary for symbiotic interactions with AM fungi (MacLean et al. 2017). Many of the genes belonging to this symbiotic pathway were also recruited for the recognition and initiation of nitrogen-fixing symbioses with rhizobia around 60 Mya, producing "the common symbiosis signaling pathway" (Oldroyd 2013). Given that AM symbiosis is believed to have only evolved once during the evolution of early land plants (Berbee et al. 2017), many components of this pathway are ubiquitous throughout the plant kingdom. Despite the early origin of symbiosis, some plants, including the model plant Arabidopsis thaliana and other members of the Brassicaceae, have lost many genes required for AM symbiosis and are therefore unable to form the symbiosis with AM fungi, resulting in a subset of genes necessary for AM symbiosis that are uniformly found in plants that host AM fungi but absent in all non-hosts (Harrison et al. 2002;Zhang et al. 2010;Gobbato et al. 2012;Wang et al. 2012). Using this common pattern, phylogenomic approaches have been used to identify additional gene families that are conserved in AM symbiosis (Favre et al. 2014;Bravo et al. 2016). While phylogenomic approaches have identified additional genes conserved for AM symbiosis across the plant kingdom, they do not generate an exhaustive list of genes associated with AM symbiosis, as some genes, such as nucleoporins, strigolactone-related genes, and genes involved in symbiosis signaling, have secondary functions not related to symbiosis (Delaux et al. 2013). These genes do not follow the pattern of strict conservation and cannot be identified by such analyses (Bravo et al. 2016). While AM fungi are typically considered symbionts, the effects of AM fungi on plant growth and fitness are varied. Meta-analyses on the effects of AM fungi inoculation in agricultural systems show that there is generally a net benefit for plant responses such as yield, biomass, and nutrient acquisition when AM fungal inoculants are added to the soil (Pellegrino et al. 2015;Schütz et al. 2018;Zhang et al. 2019), although these results are often inconsistent due to the context-dependent nature of plant response to mycorrhizal inoculants (Hoeksema et al. 2010). In addition to benefiting plants through increased nutrient uptake, AM symbiosis has a number of other positive effects on plants, such as increasing drought tolerance (Auge 2001), protecting against belowground plant pathogens (Sikes et al. 2009), ameliorating the harmful effects of heavy metal and salt stress (Audet et al. 2007;Chandrasekaran et al. 2014), and aggregating soil particles around roots (Rillig and Mummey 2006;Wilson et al. 2009).
Phenotypic diversity in AM colonization has been documented within many crop species, including sunflower (Turrini et al. 2016), wheat (Singh et al. 2012;Lehnert et al. 2017), soybean (Pawlowski et al. 2020), and even floriculture crops such as marigold (Linderman and Davis 2004). While there is evidence of variation in AM colonization in crop species, several authors have suggested that positive growth responses to AM colonization in some crop plants may have been reduced in the last century of plant breeding due to inadvertent selection, particularly in high phosphorus conditions (Lehmann et al. 2012;Martin-Robles et al. 2018). There is indeed evidence of reduced variation in mycorrhizal colonization and response in several cultivated species, though not all (Hetrick et al. 1992;Lehmann et al. 2012;Kokkoris et al. 2019). The trend of reduced root colonization with domestication has been documented in wild and cultivated sunflower, both as part of this work (Methods S1, Fig. S1) and in at least one previous study (Turrini et al. 2016).
This study examines the phenotypic and genetic diversity in AM colonization in cultivated sunflower (Helianthus annuus L.) by leveraging a genome-wide association framework to identify the genetic architecture contributing to differences in AM colonization within this valuable crop species. Sunflower is a globally important oilseed crop that was domesticated approximately 4000 years ago from the conspecific wild sunflower by Native Americans and, like other members of the Asteraceae family, generally forms Arum-type AM morphologies (Blackman 2011;Turrini et al. 2016). In sunflower, GWAS has been used to elucidate the genetic architecture of numerous traits such as Sclerotinia head rot resistance (Fusari et al. 2012), flowering time and branching (Mandel et al. 2013;Nambeesan et al. 2015), root and seedling morphology (Masalia et al. 2018), salt tolerance and ionomics (Temme et al. 2020), and floral architecture (Dowell et al. 2019). Although identification of genomic regions underlying AM colonization in sunflower has not yet been attempted, association mapping of this trait has been attempted in several other species, including wheat (Lehnert et al. 2017), sorghum (Leiser et al. 2016), and soybean (Pawlowski et al. 2020). The large continuous variation in AM colonization in diversity panels of crop species suggests that AM colonization as a trait is highly polygenic, with many genes of small effect Leiser et al. 2016;Lehnert et al. 2017;Davidson et al. 2019). The objectives of this study are to quantify the variation in AM colonization in cultivated sunflower, identify genomic regions associated with increased AM colonization, and identify candidate genes involved with AM colonization.

Germplasm resources
We used the sunflower association mapping (SAM) panel developed and described by Mandel et al. (2011Mandel et al. ( , 2013. This mapping panel consists of 288 lines that are maintained by the USDA North Central Regional Plant Introduction Station (NCRPIS) and the French National Institute for Agricultural Research (INRA). The panel is estimated to capture ca. 90% of the allelic diversity in the cultivated sunflower gene pool, which makes it ideal for the robust mapping of phenotypic traits in crop sunflower (Mandel et al. 2011(Mandel et al. , 2013. The mapping panel includes lines from oilseed and confectionary market types, as well as both major heterotic groups, relevant to the use of cytoplasmic male sterility in hybrid seed production. A subset of twelve lines, referred to as the Core 12 (Mandel et al. 2011), capture ca. 50% of the allelic diversity within the SAM panel and thus provide a broad cross section of genetic diversity in cultivated sunflower.

Experimental design
In the fall of 2018, we performed phenotyping of the degree of colonization exhibited across the SAM population by Rhizophagus intraradices (N.C. Schenck and G.S. Smith) C. Walker and A. Schüßler. We selected R. intraradices to inoculate plants because it is a generalist species of AM fungus commonly used in research (Stockinger et al. 2009). This species was first described by Schenk and Smith (1982) after being recovered from crops around the state of Florida. The fungal inoculum R. intraradices UT126 was originally cultured on peanuts and was donated to the International Culture Collection of (Vesicular) Arbuscular Mycorrhizal Fungi (INVAM) at the University of West Virginia (https:// invam. wvu. edu/) who provided the culture used for this experiment. The inoculum provided included roots, spores, hyphae, and growth medium, and was proliferated on sudangrass (Sorghum sudanense) at INVAM. We grew four replicate plants of each of the 288 lines in 5-L tree pots (Model CP512, Stuewe and Sons, Tangent, Oregon) in an evaporatively cooled greenhouse with ambient lighting on the University of Central Florida campus (28.591297, − 81.192183). Four greenhouse benches were used, constituting spatial blocks, with one replicate plant per line randomized spatially within each block. Pots were filled with uniformly screened coarse silica sand (low baseline fertility) which was sterilized using an electric soil sterilizer (Model SS-15, Pro-Grow, Phoenix, Arizona). A mesh screen was placed at the bottom of each pot to stabilize and retain the sand. For efficiency and space, we decided against growing control plants for each line because non-inoculated controls should experience no root colonization and would require the growth of more than a thousand additional plants to obtain no additional phenotypic data. Instead, to ensure our soil sterilization treatment was successful, we grew an additional four replicate plants of each Core 12 line to serve as non-inoculated controls and distributed these pots randomly throughout the four blocks. The non-inoculated controls also allowed for the examination of the relative effects of AM colonization on plant growth from a representative subsection of the SAM panel. Controls were flagged with red tape to avoid any accidental contamination during the experiment but were otherwise treated identically to all other plants apart from fungal inoculation. With four replicates of all 288 lines, plus four additional replicates per line for the Core 12 serving as controls, a total of 1200 pots were prepared.
Before the seeds were planted, each non-control pot was inoculated with 10andnbsp;mL of a homogenized mixture containing 90% moist sand and 10% fungal inoculum by volume, added directly under the seed in each planting hole. After inoculation, two seeds were planted in each pot and were thinned to one plant per pot 10andnbsp;days after planting. Pots were watered every other day until most plants had established one true leaf pair and then watered to field capacity biweekly. Plants were fertilized once per week with 200andnbsp;mL of modified half-strength Hoagland solution that was determined in an earlier trial to be optimal for AM colonization (Methods S1). After 6andnbsp;weeks of growth, the root system of each plant was gently removed from the sand, rinsed, and transported to the lab for processing and phenotyping.

Analysis of aboveground biomass
For the Core 12 lines, the shoots of both the inoculated replicates and non-inoculated controls were harvested and dried at 60 °C in a forced-air drying oven until they reached a constant mass. They were then weighed to obtain aboveground biomass. The dried plant tissue samples were sent to the Louisiana State University Extension Soil Testing and Plant Analysis Laboratory for nutrient analysis of the tissue using inductively coupled plasma mass spectrometry (ICP-MS) (Benton 1991). Mycorrhizal growth response for each line was calculated as follows (Hetrick et al. 1992): Nutrient content was used to measure crop agronomic performance (Pellegrino et al. 2015). Total grams of phosphorus and potassium in the dried shoots of each plant were calculated by multiplying the total biomass by the concentration of the phosphorus or potassium per gram of plant tissue, giving the estimated total of phosphorus or potassium in the shoot tissue. Shoot P and K responses for each line were then calculated as follows:

Determination of AM colonization
Roots were stained following the traditional protocol from Phillips and Hayman (1970) with a few modifications. Roots from the harvested plants were cut into 1-2-cm segments and 0.10-0.15 g of root sample was placed in a histology cassette. The root samples were then moved to a boiling 2% KOH (w/v) solution. After 15 min, the cassettes were rinsed three times and placed in 2% HCl (v/v) solution at room temperature. After another 15 min, the cassettes were moved to a 1:1:1 solution of 85% (v/v) lactic acid, glycerol, and distilled water with 0.05% (v/v) trypan blue stain. The cassettes were left to stain overnight, then rinsed, and frozen at -18 °C for later microscopic examination. The most common way to score AM colonization is with the gridline intersect method (Giovanetti and Mosse 1980), but this method is prone to errors due to misinterpretation by observers, particularly those inexperienced in scoring roots (Dodd and Jeffries 1986;McGonigle 1990). Because the scale of phenotyping required several individuals of different experience levels to assist with scoring roots, a new method was devised to control for observer bias. Ten randomly selected root pieces from each replicate plant, representing about 2-50% of each plant's total root system, were examined under 100 × magnification. A single random image was taken of each root piece on the focal plane that showed the entirety of the colonized parts of the root (Fig. S2). Ten evenly spaced, vertical lines were superimposed on each image and a single experienced observer counted the number of lines that were intersected by an AM structure. This method made it possible for several people to prepare samples and create images of the roots, with no error introduced due to variation among observers. The single observer validated this method against the traditional gridline intersect method using the Core 12 samples (n = 47), with an R 2 value of 0.839 (Fig. S3).
Colonization is highly variable, even in a controlled greenhouse setting (Davidson et al. 2019). In any greenhouse experiment, environmental gradients can influence plant phenotypes. To account for microclimate variation within the greenhouse, the spatial layout of pots was converted into a grid system, where x and y positions for each pot were used to fit a standard least-squares linear model, generating least-squares means to account for possible variation in airflow, temperature, or light levels. These least-squares means reflect the normalized average for colonization of each line and this is the primary response variable of interest, referred to hereafter as mean colonization. Two other complementary metrics were calculated-the maximum colonization among the four replicates (representing the highest potential colonization observed for each line) and the range of colonization within each line (as a metric of variability in colonization under otherwise relatively uniform environmental conditions).

Genome-wide association studies
GWAS was performed using the methods of Temme et al. (2020). This pipeline uses the Genome-wide Efficient Mixed Model Association (GEMMA) algorithm to account for population structure as identified by principal component analysis and familial relatedness (Zhou and Stephens 2012).
Twenty-seven lines exhibited higher than expected rates of residual heterozygosity and are not included in these analyses, leaving a total of 261 lines that were used for GWAS (Masalia et al. 2018;Temme et al. 2020). Lines that were not retained for mapping analysis were disproportionately landraces and open-pollinated varieties, so these were included in phenotyping to provide colonization data on these important sunflower germplasm pools.
A traditional Bonferroni correction is generally used to prevent false positives, but due to the non-independent nature of the markers in this dataset, this correction would result in a significance threshold that is far too conservative, resulting in very low statistical power, so a Gao correction which adjusts the significance threshold based on the effective number of independent tests was used instead (Gao et al. 2008;Gao et al. 2010). Temme et al. (2020) used the PLINK v1.9 indep-pairwise linkage disequilibrium (LD) variant pruning function to estimate the effective number of independent tests (Chang et al. 2015). This analysis used a sliding window of 100 kbp, a SNP step size of 10, and an R 2 threshold of 0.8. This approach identified LD blocks across the sunflower genome, and the total number of blocks was used to calculate the Gao significance threshold (α = 0.05/estimated number of independent tests). Similarly, instead of classifying all SNPs as independent associations, SNPs were binned into significantly associated regions using observed patterns of linkage disequilibrium with LDSelect v1.0 (Carlson et al. 2004). SNPs were collapsed into regions based on an R 2 threshold of 0.8, a threshold observed to be sufficient for defining independent blocks in association studies (Masalia et al. 2018). Heritability estimates were obtained using the R-package heritability (Kruijer et al. 2015).
A list of genes within significant regions was compiled using the HA412-HO genome assembly v.2 (Temme et al. 2020;Todesco et al. 2020). A list of genes within "suggestive" Fig. 1 Frequency distribution of mean root colonization by Rhizophagus intraradices for 288 lines from the mapping panel, 6 weeks after planting. Least-squares means are reported to account for environmental gradients in the greenhouse. Mean colonization for lines ranged from 0.0 to 70.1% regions that contained the top 0.01% of SNPs across the genome was also compiled. Both sets of genes were compared with a list of genes that are differentially expressed in sunflowers when colonized by AM fungi (Vangelisti et al. 2018) as well as with a list of conserved genes for AM symbiosis across land plants (Bravo et al. 2016).

Phenotypic variation of AM symbiosis
Light microscopy revealed the presence of mycorrhizal structures in most, but not all, inoculated plants (Fig. S2). Mycorrhizal structures were not present in the non-inoculated control plants, indicating that soil sterilization was successful. The SAM panel had considerable phenotypic variation in mean colonization, with line means ranging from 0.0-70.1% (Fig. 1, Data S1). The mean colonization among lines in the panel was 33.2%. The line means of AM colonization across the panel formed a normal distribution (Shapiro-Wilk, P = 0.1315). Within the spatial least-squares linear model, the spatial position of the pots was statistically significant but explained only a very small portion of the variation. Genotype arithmetic and least-squares means differed on average by only 0.04% (± 0.72% S.D.) and the largest observed difference was only 3%, demonstrating that this spatial correction of phenotypic data would not likely alter the results of subsequent analyses. Within the Core 12 lines, the variation in mean colonization ranged 11.7-43.3% and averaged 28.8%. Heritabilities for the mean colonization, maximum colonization, and range of colonization within the SAM panel were calculated as 0.227, 0.152, and 0.088 respectively.
Lines in the Core 12 (Data S2) exhibited high variability in mycorrhizal growth response (MGR). Core 12 lines had anywhere from a 29% decrease in biomass to an 88% increase in biomass when inoculated with AM fungi or MGR values from − 0.29 to 0.88 (Figs. S4 and S5). MGR

Fig. 2 Manhattan plots for
A mean root colonization, B maximum root colonization, and C range of root colonization. The genomic position of each SNP is represented on the x-axis, colored by chromosome, and the negative logarithm of the association p-value for each SNP is presented on each y-axis. Red lines represent the significance threshold (p = 2.43 × 10 −6 ) using a Gao correction with the effective number of tests (n = 20,562). Blue lines represent the threshold for suggestive SNPs, which were the top 0.01% of all SNPs was significantly correlated with elevated shoot P content (P = 0.0012) and shoot K content (P = 0.0002) in colonized plants (Fig. S4). In contrast, mean colonization was not correlated with MGR, shoot P, or K response variables, and lines with highly colonized roots among the Core 12 had both high and low MGR and shoot P and K responses.

GWAS results
As previously described elsewhere (Temme et al. 2020;Todesco et al. 2020), a total of 1.81 million SNPs were identified from the 261 lines included in this work. These SNPs were distributed across all 17 chromosomes of the 3.6 Gbp genome of cultivated sunflower. Filtering to remove markers with > 30% missing data or heterozygous loci and markers with minor allele frequency (MAF) < 5% resulted in a set of 1.47 million SNPs that were used in this analysis. The total number of LD blocks, and therefore the effective number of tests, was calculated to be 20,562, resulting in an adjusted significance threshold of 2.43 × 10 −6 .
For the three focal traits used in the mapping (mean root colonization, maximum root colonization, and range of root colonization), one significant SNP was identified for each trait (Fig. 2). Mean root colonization mapped to a single SNP (HA412HOChr07:148,640,639) on chromosome 7. This single-SNP block explained 19.7% of the variance in mean root colonization. Maximum colonization mapped to a single SNP (HA412HOChr02:183,246,676) on chromosome 2 which explained 14.4% of the variance in the trait. The range of colonization mapped to a single SNP (HA412HOChr05:179,283,961) on chromosome 5 which explained 27.9% of the variance in the trait. Sixteen genes were located in the three associated regions for the traits tested (Table 1). One of the two genes located on the significant region on linkage group 2 was Ha412HOChr02g0090921, a member of a large clade of ankyrins that have a PGG domain. The gene annotation indicates that it is an integral part of the cellular membrane (GO:0,016,021) and may be important for membrane transport. Transcripts of this gene were also found to be differentially regulated in the transcriptome of a sunflower when colonized with AM fungi (Vangelisti et al. 2018). Two genes of interest in significant regions on linkage groups 5 and 7, Ha412HOChr05g0240851 and Ha412HOChr07g0318991, are involved with lipid metabolic processes (GO:0,016,298, GO:0,006,629, GO:0,004,806) and their predicted products were annotated as putative triacylglycerol lipases (Table 1). A clade of this subfamily containing one of these genes, Ha412HOChr05g0240851, is also conserved for AM symbiosis (Bravo et al. 2016). The genes within this clade were identified as GDSL lipases-esterases.
An additional eight genes were identified by crossreferencing the top 0.01% of SNPs with compiled lists of differentially expressed and conserved genes for AM symbiosis in sunflower ( Table 2). Three of these genes, HA412HOChr08g0327781, HA412HOChr08g0327791, and HA412HOChr08g0327801, were identified as a single gene in the sunflower genome annotation used for the classification of differentially expressed genes as well as in the conserved gene analysis. This gene is in the conserved  (Vangelisti et al. 2018). All significant and suggestive genes are reported in Data S3.

Discussion
One objective of this study was to understand the phenotypic variation in AM symbiosis present within the cultivated sunflower germplasm using a large, representative population. Previous studies have shown that there is a normal distribution of variation in AM colonization within sunflower, wheat, and other species (Linderman et al. 2004;Singh et al. 2012;Lehnert et al. 2017;Salloum et al. 2016;Turrini et al. 2016;Davidson et al. 2019). The distribution in mean colonization observed in the SAM mapping panel approximately matched the distribution of colonization in a small set of 11 cultivars chosen by Turrini et al. (2016) which exhibited colonization of 8.6-78.7%, very similar to the 0.0-70.1% observed in this experiment (Fig. 1). The distribution of phenotypic variance in root colonization seen among genotypes in this mapping panel is unimodal and symmetric. This suggests that colonization is highly polygenic, with likely additive effects among many genes of varying effect size, whereas a bimodal or strongly skewed genotypic distribution would indicate alternative genetic architectures, for example, oligogenic architectures or strongly epistatic interactions among a few major genes (Mackay 2001;Barghi et al. 2020). Our finding of few loci of relatively large effect size is consistent with reported trait distributions, as many of the loci of smaller effect cannot be detected by GWAS alone, with polygenic architectures expected to follow an exponential distribution of allelic effects (Mackay 2001). A polygenic trait with largely additive variance is also likely to respond in a predictable manner to directional natural selection or artificial selection, which is a relevant finding for understanding how mycorrhizal colonization has evolved under domestication and crop improvement, and how it may respond to targeted breeding efforts. The results presented by Turrini et al. (2016) suggest that the susceptibility to AM colonization has on average decreased due to genetic bottlenecks that occurred during the domestication and improvement of cultivated sunflower (Liu and Burke 2006;Mandel et al. 2011). The results of this study, however, indicate that although domestication may have decreased susceptibility to AM colonization, there is still a large amount of phenotypic variation in AM colonization within the cultivated sunflower germplasm, although the heritability for root colonization is on the lower end of estimates found in studies of other crops (Leiser et al. 2016;Lehnert et al. 2017;Davidson et al. 2019). Additionally, the substantial amount of phenotypic variability that still exists for AM symbiosis, despite reduced genetic diversity in cultivated lines, is consistent with the variation that remains in other floral, stem, leaf, and root traits in this same mapping panel (e.g., Masalia et al. 2018;Dowell et al. 2019;Temme et al. 2020). Although colonization was highly variable, it was not correlated with other metrics used to measure differences in plant growth and nutrient response to AM fungi. When compared against other root traits measured in the same mapping panel, the only statistically significant correlations were between specific root length and maximum colonization and range of colonization (R 2 < 0.15; Fig. S6; Masalia et al. 2018). This suggests that colonization may be somewhat affected by root tissue density, but the weak correlation between these traits suggests that this was not the primary driver of colonization traits. In some previous studies, high AM colonization has been associated with a large mycorrhizal growth response (Lehmann et al. 2012;Treseder 2013), but others have not found this correlation  (Vangelisti et al. 2018) or identified as a conserved gene (CON) for AM symbiosis (Bravo et al. 2016). Genes within significant regions (Table 1) are not included in  this table   Gene Predicted product Support  (Kaeppler 2000;Hohmann and Messer 2017). We found that there was a high variation in mycorrhizal growth response in the Core 12 lines and that this variation could not be explained by variation in colonization. These results also indicate that depending on the plant genotype, the AM fungus tested (R. intraradices) can cause reduced growth in some plant genotypes when colonized, whereas other genotypes exhibit large increases in biomass after inoculation (Fig. S4). A different study that tested the effects of inoculation of 18 diverse sorghum genotypes with four different fungal species found a similar trend of differential responses in plant genotypes for all four species of AM fungi tested (Watts-Williams et al. 2019). Those authors also found strong correlations between shoot P response and MGR in diverse genotypes, similar to the relationship between MGR and shoot nutrient responses in the present study. These findings have interesting implications because while the variability in MGR and shoot P and K response has been studied in the context of interspecific variation (van der Heidjen et al. 1998;Bever 2003), intraspecific variation might also be a critical component to understanding the evolution of AM symbiosis among plant lineages and for understanding the structure and productivity of plant communities. Large variation in MGR and shoot nutrient responses occur in other species when genotypes are grown under the same environmental conditions, including sorghum and wheat (Ellouze et al. 2016;Watts-Williams et al. 2019). Some authors have suggested that mapping MGR may provide more potential SNPs that would be useful for breeding AM traits in plants (Kaeppler et al. 2000;Ellouze et al. 2016;Watts-Williams et al. 2019). The results from this analysis further indicate that MGR and shoot nutrient response are likely more useful for estimating the positive benefits of AM symbiosis than root colonization alone and that MGR may be an important variable for further investigation as it could reveal a role for different genomic regions than those detected when mapping AM colonization alone. Despite its lack of application to growth benefits received from AM symbiosis, uncovering the genomic regions associated with root colonization may still be an important step for understanding control of colonization in plant roots and should be considered in future studies assessing MGR. Plant genotypes with genomic regions associated with neutral responses and low colonization under high nutrient conditions may be good at modulating nutrient trade when conditions are not favorable. The degree of AM colonization may be associated with other plant benefits that are not related to nutrient uptakes such as increased performance in water limiting conditions and pathogen resistance (Auge 2001;Sikes et al. 2009). While this experiment only included nutrient stress, it is possible that under other environmental conditions there might be a stronger correlation between colonization and measures of plant benefit.
Previous attempts to use GWAS to determine the genetic basis of AM colonization have indicated that the trait is likely controlled by many genes of small effect size. In sorghum, no SNPs above the Bonferroni significance threshold were detected, and the three SNPs with the lowest p-values could each only explain ca. 8% of the observed variation (Leiser et al. 2016). The phenotypic variation explained by a given marker in winter wheat was even lower, with no marker explaining more than 1.14% (Lehnert et al. 2017). In a study on rice cultivars, no SNPs were found above the most stringent significance threshold, but 23 putative SNPs of relatively small effect were revealed when suggestive regions were compared with transcriptomic data (Davidson et al. 2019). In soybean, six significant loci explained 24% of the variation in root colonization, but none of the loci explained more than 7.1% of this variation (Pawlowski et al. 2020). High amounts of variation in colonization as well as low amounts of phenotypic variation explained by each genomic region in these studies suggests that there are likely many genes of small effect controlling colonization that have gone undetected due to conservative significance thresholds and low effect sizes. Surprisingly, the variation explained by the three significant loci in sunflower was considerably higher than expected, with each locus explaining 14.4-27.9% of the variation in their respective traits. To the best of our knowledge, this is the first attempt to use GWAS to identify genomic regions associated with AM colonization in an Asterid eudicot species, so it is possible that there are genes of larger effect outside the Poaceae crops (i.e., rice, wheat, and sorghum) that have primarily been studied to date. This would be consistent with recent findings in soybean, to date the only other non-Poaceae crop examined by GWAS (Pawlowski et al. 2020).
An important caveat is that these results may be specific to the particular fungal culture that we tested (i.e., R. intraradices UT126). Some plant genotypes may have higher or lower colonization when colonized by other AM fungus cultures. It is difficult to know without testing additional fungus species (and genotypes within species) whether the significant regions found in this study are associated with differences in the ability of a plant genotype to form AM symbiosis generally or if they are indicative of a specific response to this strain of R. intraradices. When inoculated with many AM fungus species, plant host genotype has been found to affect AM fungal community composition, while not affecting the overall rate of colonization, indicating that AM colonization may be constant in diverse AM fungal communities and that differences between plant cultivars are due to their differences in specificity to AM fungi (Mao et al. 2014). It is thus possible that observed differences in colonization in this study are not indicative of each line's ability to form AM symbiosis but rather its specific response to R. intraradices. If true, these results nevertheless still provide valuable information about genes that may be responsible for partner specificity in AM symbiosis.
The three significant SNPs identified for root colonization traits are in genomic regions with genes that are known to be important for plant immunity and defense. The significant region on chromosome 2 contains a gene that is a member of a large clade of ankyrins that have a PGG domain (Ha412HOChr02g0090921). While the function of these genes is largely unknown, this gene family has been classified as a mycorrhizal-expanded gene family (Kramer et al. 2019). Other members of this gene family are ACCELER-ATED CELL DEATH 6 (ACD6), which plays a role in systemic acquired resistance and immunity, and INEFFECTIVE GREENISH NODULES (IGN1), which prevents premature senescence of cells infected with nitrogen-fixing bacteria at the onset of infection (Rate et al. 1999;Kumagai et al. 2007). Transcripts of this sunflower ankyrin gene were also found to be differentially regulated in the transcriptome of a sunflower when colonized with AM fungi (Vangelisti et al. 2018), which further supports the validity of this association. Within the other two significant regions are two genes encoding putative triacylglycerol lipases (Ha412HOChr05g0240851, Ha412HOChr07g0318991). When subjected to a BLAST search against the Medicago truncatula genome, these two genes were identified as GDSL lipases, a conserved gene family in plants that form AM symbiosis (Bravo et al. 2016). GDSL lipases are a large subfamily of lipases that shows very broad substrate specificity (Akoh et al. 2004). These enzymes are involved in lipid metabolism and are believed to play roles in modulating defense responses, notably against pathogenic fungi (Lee et al. 2009;Gao et al. 2017;Lai et al. 2017). The importance of lipid metabolism in mycorrhizal symbiosis has been increasingly recognized in the last several years, as lipids produced in the plant are believed to be sent to the fungus using an array of uniquely conserved genes, such as RAM1, RAM2, and STR (Bravo et al. 2016, Keymer et al. 2017Luginbuehl et al. 2017).
Analysis of suggestive genes using published transcriptomics data also provided another subset of genes that may be important for AM colonization, even though they did not surpass the stringent significance threshold. The suggestive genes located on chromosome 8 (HA412HOChr08g0327781, HA412HOChr08g0327791, and HA412HOChr08g0327801) are annotated as a single gene that is in the conserved HYP5c gene family identified by Bravo et al. (2016). While the function of this gene has not been characterized, this gene may be involved in the glycosylation of cell membrane components (Konečný et al. 2019). HA412HOChr05g0204541 and HA412HOChr05g02044891 are both transporters that are differentially transcribed under AM symbiosis and are also found to be differentially transcribed in Medicago truncatula during rhizobial symbiosis (Boscari et al. 2013), indicating that these two genes belong to the common symbiont signaling pathway. HA412HOChr07g0289501 and HAHO412Chr07g0289521 are both differentially transcribed in AM symbiosis in sunflower and appear to be gene copies with predicted products of salutaridine reductase, which is potentially involved in pathogen defense (Ziegler et al. 2009). Finally, HA412HOChr13g0600911 is a cyclic nucleotide-gated channel that activates the Ca 2+ oscillations that are necessary for the establishment of both rhizobial and AM symbiosis (Charpentier et al. 2016). Given the conservative significance threshold of the Gao correction, it is possible that there are more genes of smaller effect that we lacked the power to detect.
In summary, we characterized variation in AM fungal colonization in a cultivated sunflower diversity panel and GWAS to study its genetic basis. We also used a subset of lines to determine the relationship between AM colonization, mycorrhizal growth response, and shoot nutrient response. There was a large amount of variation in mean colonization across lines, but this trait did not correlate with mycorrhizal growth response or shoot nutrient response. GWAS yielded three significant associations on linkage groups 2, 5, and 7. Two of these significant regions contained genes belonging to a subfamily known to be partially conserved for AM symbiosis, and the other significant region included a gene that belongs to a gene family that has expanded in the lineages of AM hosts. Using the top 0.01% of SNPs, we identified eight additional colonization-associated candidate genes that were also found to be differentially regulated during AM symbiosis or conserved for AM symbiosis across land plants. Some of the candidate genes identified have been studied in the context of AM symbiosis and are involved with defense, signaling, and transport, but more research needs to be done on other candidates to understand their potential function during AM symbiosis.

Competing interests
The authors declare no competing interests.