Genome‐wide association study reveals novel players in defense hormone crosstalk in Arabidopsis

Abstract Jasmonic acid (JA) regulates plant defenses against necrotrophic pathogens and insect herbivores. Salicylic acid (SA) and abscisic acid (ABA) can antagonize JA‐regulated defenses, thereby modulating pathogen or insect resistance. We performed a genome‐wide association (GWA) study on natural genetic variation in Arabidopsis thaliana for the effect of SA and ABA on the JA pathway. We treated 349 Arabidopsis accessions with methyl JA (MeJA), or a combination of MeJA and either SA or ABA, after which expression of the JA‐responsive marker gene PLANT DEFENSIN1.2 (PDF1.2) was quantified as a readout for GWA analysis. Both hormones antagonized MeJA‐induced PDF1.2 in the majority of the accessions but with a large variation in magnitude. GWA mapping of the SA‐ and ABA‐affected PDF1.2 expression data revealed loci associated with crosstalk. GLYI4 (encoding a glyoxalase) and ARR11 (encoding an Arabidopsis response regulator involved in cytokinin signalling) were confirmed by T‐DNA insertion mutant analysis to affect SA–JA crosstalk and resistance against the necrotroph Botrytis cinerea. In addition, At1g16310 (encoding a cation efflux family protein) was confirmed to affect ABA–JA crosstalk and susceptibility to Mamestra brassicae herbivory. Collectively, this GWA study identified novel players in JA hormone crosstalk with potential roles in the regulation of pathogen or insect resistance.

respond effectively to each attacker or to multiple attackers at the same time, hormonal signalling pathways cross communicate in antagonistic or synergistic manners. In particular, SA and ABA have been shown to interact with the JA pathway, thereby strongly modulating the JA-induced defense output (Pieterse et al., 2012).
In Arabidopsis thaliana (hereafter Arabidopsis), two distinct branches of the JA pathway have been shown to antagonize each other: the ethylene response factor (ERF) branch and the MYC branch, which are coregulated by ET and ABA, respectively Pieterse et al., 2012;Robert-Seilaniantz et al., 2011).
SA has been reported to have a major impact on JA-induced defenses in both the ERF and the MYC branch of the JA pathway (Bostock, 2005;Pieterse et al., 2012;Stout, Thaler, & Thomma, 2006). Although the effect of SA on the JA pathway can be antagonistic, synergistic, or neutral, in Arabidopsis, antagonistic interactions seem to prevail (Pieterse et al., 2012;Tsuda, Sato, Stoddard, Glazebrook, & Katagiri, 2009)

. Experiments performed with
Arabidopsis revealed that the JA-responsive genes PDF1.2 and VSP2 are highly sensitive to suppression by SA. In many cases, this antagonism between the SA and JA pathways affects plant resistance against necrotrophs or insect herbivores (Caarls, Pieterse, & Wees, 2015). Suppression of the JA pathway by SA is predominantly regulated at the level of gene transcription (Caarls et al., 2015;Van der Does et al., 2013). Important regulators of the interaction between the SA and JA pathways have been identified, such as the redox sensitive transcriptional coregulator NONEXPRESSOR OF PATHOGENESIS-RELATED PROTEINS1 (NPR1; Spoel et al., 2003) and several WRKY and TGA transcription factors (Caarls et al., 2015).
Moreover, SA was shown to promote degradation of the transcription factor ORA59 (Van der Does et al., 2013) and to inhibit ORA59 gene expression (Zander, Thurow, & Gatz, 2014), providing a mechanistic explanation of how SA suppresses the ERF branch of the JA pathway.
Like SA, ABA is also a strong modulator of JA-induced defenses.
Using a small set of 18 Arabidopsis accessions, we previously demonstrated that all tested accessions were sensitive to SAmediated suppression of the JA-responsive marker gene PDF1.2, albeit to different extents (Koornneef et al., 2008), highlighting the potential significance of crosstalk between the SA and JA pathway in induced plant defenses in nature (Thaler, Humphrey, & Whiteman, 2012). We reasoned that the observed natural genetic variation in the level of hormone crosstalk would provide a so far unexplored resource from which genes encoding novel players in hormone crosstalk could be identified. Recent advances in genotyping and sequencing technology have made genome-wide association (GWA) mapping a good approach to mine natural genetic variation and detect molecular markers linked to, for example, stress resistance traits (Assmann, 2013;Atwell et al., 2010;Meijon, Satbhai, Tsuchimatsu, & Busch, 2014). GWA mapping is a method initially utilized in human population studies to identify the genetic basis of complex traits (Hirschhorn & Daly, 2005). GWA mapping has also been successfully utilized in plant studies (Aranzana et al., 2005;Atwell et al., 2010;Bac-Molenaar et al., 2015;Broekgaarden et al., 2015;Chan, Rowe, & Kliebenstein, 2010;Kloth, Thoen, Bouwmeester, Jongsma, & Dicke, 2012 of GWA mapping is that natural variation in the phenotype of a given trait is caused by genetic differences in the population under study. Via single nucleotide polymorphisms (SNPs) in the genomes of the phenotyped genotypes, these genetic differences can be linked to genetic loci and ideally to the responsible genes. In this study, we used the GWA mapping resource of the Arabidopsis Haplotype Map (HapMap) collection that is comprised of over 360 natural accessions collected globally and genotyped for~250k SNPs relative to the Col-0 accession Kim et al., 2007;Nordborg et al., 2005;Weigel & Mott, 2009). We observed a large genetic variation in the magnitude by which SA and ABA affect the expression of the JA-responsive marker gene PDF1.2 in the different Arabidopsis accessions. GWA mapping revealed multiple loci associated with either SA-JA or ABA-JA crosstalk. Several candidate genes were validated using T-DNA insertion mutant analysis, which yielded a number of novel players in SA-JA and ABA-JA crosstalk with effects on the level of resistance against the necrotrophic pathogen B. cinerea or the herbivorous insect Mamestra brassicae.  Table S1). Figure 1a shows that the basal level of 2.2 | GWA mapping reveals loci associated with the effect of SA and ABA on JA-responsive gene expression To identify novel players involved in the effect of SA and ABA on JAresponsive gene expression, we performed GWA mapping using thẽ 214k SNP set that is commonly used for GWA studies in Arabidopsis (Atwell et al., 2010;Bac-Molenaar et al., 2015;Horton et al., 2012;Kim et al., 2007;Li, et al., 2010). To this end, we used normally distrib-   Arabidopsis, the linkage disequilibrium of the population of Arabidopsis accessions used for GWA mapping is estimated 10-50 kb (Kim et al., 2007;Nordborg et al., 2005). Therefore, we considered all genes within 15 kb upstream and downstream of each significant SNP to be candidates for the observed SNP-trait associations. This yielded a list of 60 candidate genes in nine loci for the MeJA + SA/MeJA dataset and 98 candidate genes in eight loci for the MeJA + ABA/ MeJA dataset (Table S2). Their putative functions based on the Gene Ontology tool of The Arabidopsis Information Resource (Lamesch et al., 2012) are given in Table S2.

| Fine mapping of GWA-identified loci
To validate the identified SNP-trait associations in the GWA mapping analyses, we performed locus-specific mapping ( Table S1). In the 15 kb upstream and downstream of each significant SNP identified by both AMM and FAST-LMM (-log 10 (p) > 4; MAF > 5%), SNP-trait associations were identified by LSM. The significance of associations between traits and SNP markers was evaluated using a nonparametric Kruskal-Wallis test (Filiault & Maloof, 2012   Note. Shown are Arabidopsis gene identifier (AGI) numbers of candidate genes that are located in the closest proximity of the identified highly associated  homozygous T-DNA insertion lines of the respective candidate genes (Table S3)  Shown are Arabidopsis AGI numbers of candidate genes that are located in the closest proximity of the identified highly-associated SNPs [log 10 ( In Arabidopsis, B. cinerea induces JA-dependent defenses in the plant that effectively suppress disease , Windram et al., 2012. In turn, B. cinerea has been shown hijack the SA pathway to suppress effective JA-dependent defenses via SA-JA crosstalk (El Oirdi et al., 2011). Therefore, we next tested glyI4 and arr11 for their level of resistance against B. cinerea. To this end, 5-week-old plants were inoculated with B. cinerea spores, and disease symptoms were scored 3 days later. Figure 4c and 4d show that glyI4 developed significantly less-severe disease symptoms than Col-8. Conversely, arr11 developed significantly more-severe disease symptoms than Col-8.
Together, these results indicate that reduced SA-JA crosstalk in glyI4 and enhanced SA-JA crosstalk in arr11, contrastingly affect the level of resistance against B. cinerea. Loss of SA-JA crosstalk correlates with enhanced resistance, while increased SA-JA crosstalk is associated with enhanced susceptibility to this necrotrophic pathogen.

| T-DNA insertion lines analysis of candidate genes associated with ABA-JA crosstalk
To investigate whether the selected candidate genes (Table 2)  In Arabidopsis, the synergistic interaction of ABA on the MYC branch of the JA pathway is associated with increased resistance to herbivory (Vos et al., 2015). Feeding by the leaf-chewing insect M.
brassicae induces the MYC-branch and enhances the expression of the ABA-JA responsive gene VSP2 (Pangesti et al., 2016). To test whether the impaired ABA-JA crosstalk phenotype of mutant At1g16310 is associated with changes in the level of resistance against M. brassicae feeding, we performed an insect-resistance bioassay with this herbivore. One first-instar M. brassicae caterpillar was placed on each plant and allowed to feed for 14 days, after which the caterpillar was weighed. Figure 5d shows that the caterpillars were significantly heavier when they fed from At1g16310 plants than when they fed from Col-8 plants. These results indicating that reduced ABA-JA crosstalk in mutant At1g16310 is associated with enhanced susceptibility to M. brassicae feeding.

| DISCUSSION
Plant hormones have pivotal roles in the regulation of plant defense responses. Their signalling pathways cross communicate, which provides plants with an enormous regulatory potential to rapidly adapt to their biotic and abiotic environment (Reymond & Farmer, 1998). Hormonal crosstalk is thought to be a cost-saving strategy and may have evolved as a means of the plant to reduce allocation costs by repression of unnecessary defenses that are ineffective against the attacker that is encountered (Thaler et al., 2012;Vos et al., 2015). In Arabidopsis, the JA response pathway is particularly sensitive to antagonism by the plant hormones SA and ABA, which potentially impacts the level of JA-dependent resistance against necrotrophic pathogen and herbivorous insects (Pieterse, et al., 2012. In this study,  (Figure 1d and 1e).
Of the 349 accessions tested, 266 accessions displayed both SAand ABA-mediated antagonistic effects on PDF1.2 expression (Figure 1f), which suggests that the antagonistic effects of both SA and ABA on the JA pathway must have important functions in nature. We observed only a weak correlation between the levels of SA-mediated and ABA-mediated antagonism on the JA pathway.
Moreover, the genomic regions that our GWA mapping found to be associated with SA-JA and ABA-JA crosstalk did not overlap.
Together, these findings suggest that the negative effects of SA and ABA on JA-responsive gene expression are based on distinct molecular mechanisms and can function independently with magnitudes that depend on the genetic background.
GWA mapping is a useful tool to study the genetic architecture of traits but has not often been used for the detection of genetic variants associated with plant defense (Bartoli & Roux, 2017). Using a selection pipeline of GWA and LSM, we identified six genomic regions with sig-  (Caarls et al., 2015;Koornneef et al., 2008). It is, therefore, tempting to speculate that GLYI4 modulates SA-JA crosstalk by interfering with this process.
T-DNA insertion mutant analysis pointed to the type-B ARR (Arabidopsis Response Regulators) ARR11 as a second novel player in SA-JA crosstalk. Mutant arr11 displayed hypersensitivity to SAmediated suppression of PDF1.2 and enhanced susceptibility to B. cinerea infection. The type-B family of ARRs play a pivotal role in the early transcriptional response of plants to cytokinin (Argyros et al., 2008). Cytokinin is involved in many plant developmental processes, but it also functions as part of the hormonal network that regulates the balance between plant growth and adaptation to stress (Giron, Frago, Glevarec, Pieterse, & Dicke, 2013;O'Brien & Benkova, 2013).
For instance, another member of the ARR family, the cytokinin-activated transcription factor ARR2, has been shown to bind to the SA response factor TGA3, therewith enhancing SA/NPR1-mediated defense gene expression and plant immunity to the biotrophic pathogen Pseudomonas syringae (Choi et al., 2010). Our data support a role for ARR11 as negative regulator of SA-JA crosstalk, thus positively affecting resistance against B. cinerea. However, the mode of action of ARR11 in this process remains elusive.
The T-DNA insertion mutant for At1g16310 was the only mutant with an altered phenotype related to ABA-JA crosstalk. It showed a reduced negative effect of ABA on MeJA-induced PDF1.2 and, reciprocally, a reduced positive effect of ABA on MeJA-induced VSP2. In line with this, the level of resistance against the insect herbivore M.
brassicae was also reduced in this mutant. Mutant At1g16310 plants also appeared to be less sensitive to ABA, as exemplified by a reduced expression of the ABA-responsive gene RAB18. In Arabidopsis, the synergistic interaction of ABA on the MYC branch of the JA pathway is associated with enhanced expression of VSP2 and increased resistance to herbivory, while it antagonizes the ERF branch of the JA pathway, resulting in suppression of JA-induced PDF1.2 (Bodenhausen & Reymond, 2007;Pangesti et al., 2016;Vos, et al., 2015). Our data support a role for At1g16310 in enhancing the sensitivity of the plant to ABA, therewith stimulating the level of ABA-JA crosstalk and the level of resistance against insect herbivory. The protein encoded by At1g16310 is member of the cation diffusion facilitator family of proteins, which are important for the maintenance of cation homeostasis in bacteria, yeast, plants, and mammals. The role of cation diffusion facilitators in modulating cellular cation concentrations can impact diverse processes, including cation tolerance, oxidative stress resistance, and protein turnover (Delhaize et al., 2007). However, the molecular mechanism by which At1g16310 influences ABA sensitivity and as such ABA-JA crosstalk and herbivore resistance, is yet unknown.
In this GWA study, we pinpointed several loci in the Arabidopsis genome that are associated with variation amongst Arabidopsis accessions in the antagonistic effect of SA and ABA on the JA pathway. It is tempting to speculate that the underlying genes have been subject to evolutionary selection to shape the output of the JA signalling pathway and maximize survival under the prevailing environmental conditions. Future research will be focused on unraveling the mode of action of the identified genes and their corresponding proteins in hormonal interplay, which will increase our understanding of the ingenious ways by which plants adapt to their often hostile environment.

| Plant material and growth conditions
In this study, a total of 349 natural Arabidopsis thaliana accessions of the HapMap collection (Table S1) were used, which are genotyped for~250k bi-allelic SNPs (Baxter et al., 2010;Chao et al., 2012;Platt et al., 2010). After quality control and imputation, this SNP-set was reduced to a set of 214.051 SNPs (Thoen et al., 2017). Seeds of the Arabidopsis accessions were sown in cultivation containers filled with autoclaved river sand. Sand was supplied with half-strength Hoagland solution containing 10 μM Sequestreen (CIBA-Geigy, Basel, Switzerland) as described (Van Wees, Van Pelt, Bakker, & Pieterse, 2013). To attain a high relative humidity (RH) for germination, cultivation containers were enclosed in a tray with water and covered with a transparent lid. Seeds were stratified for 2 days at 4°C in the dark to ensure a homogeneous germination after which the tray was moved to a growth chamber with an 8-hr day/16-hr night rhythm, a temperature of 21°C, and a light intensity of 100 μmol m -2 sec -1 . After 8 days, the lids of the trays were slightly opened and gradually removed over a 2day period to adjust to the 70% RH present in the growth chamber.   (Livak & Schmittgen, 2001). Fold change in gene expression was calculated relative to the mock treatment in wild-type plants.
AMM first performs a genome-wide scan using the approximate inference proposed by Zhang et al. (2010) and Kang et al. (2010) and then updates the smallest 100 p values using an exact mixed model inference (Kang et al., 2008). This algorithm closely resembles the commonly used Efficient Mixed-Model Association eXpedited (Kang et al., 2010). Additionally, GWA analysis was performed using Fast-LMM, as described by Cao et al. (2011). FAST-LMM captures all cofounders in the population structure simultaneously as LMM with the advantage to process larger dataset, making the analysis faster (Lippert et al., 2011). SNPs with a MAF < 5% were not considered in both models because of possibly elevated false-discovery rates (Atwell et al., 2010). The GWAPP Geneviewer was used to zoom in on trait-associated SNPs and reveal their position in the genome to pinpoint candidate genes within 15 kb upstream and downstream of the identified SNP. For each of the candidate genes, the annotations were retrieved from The Arabidopsis Information Resource10 (arabidopsis.org).

| DNA isolation and genotyping
T-DNA insertion mutant lines were routinely genotyped using genomic DNA isolated with the sucrose method (Berendzen et al., 2005).
Primers for T-DNA insertion mutant genotyping were designed using the SIGnAL T-DNA verification primer design tool (http://signal.salk. edu/cgi-bin/tdnaexpress) and are listed in Table S4. DNA amplification was performed over 34 cycles in a Biorad Thermal cycler using Phire Hot Start II DNA Polymerase (Life Technologies, Bleiswijk, the Netherlands) and the following PCR conditions: denaturation at 98°C for 5 s, annealing at 60°C for 10 s, and elongation at 72°C for 20 s. PCR products were separated by agarose gel electrophoresis.

| Pathogen and insect bioassays
Botrytis cinerea strain B05.10 (Van Kan, Van 't Klooster, Wagemakers, Dees, & Van der Vlugt-Bergmans, 1997) was grown for 2 weeks on half-strength potato dextrose agar (PDA; Difco Laboratories, Leeuwarden, the Netherlands) plates containing penicillin (100 μg ml -1 ) and streptomycin (200 μg ml -1 ) at room temperature as described previously . B. cinerea spores were subsequently collected, filtered through glass wool, and resuspended in half-strength potato dextrose broth (PDB; Difco Laboratories, Leeuwarden, the Netherlands) to a final density of 1 × 10 5 spores ml -1 . After a 3-hr incubation period, 5-week-old plants were inoculated by applying 5-μl droplets of the spore suspension to six leaves of each plant. Plants were placed under a lid to increase RH to 100% to stimulate the infection.
Three days after B. cinerea inoculation, lids were removed and the symptoms were scored in four disease severity classes ranging from no symptoms (Class I), nonspreading lesion (Class II), spreading lesion (Class III), up to severe spreading lesions with tissue maceration (Class IV) .
Growth of M. brassicae larvae was assessed over a period of 14 days.
To this end, a single freshly-hatched first-instar (L1) larvae was placed on each plant. After 14 days of growth, larval weight was measured as described .