Contributed by Joseph R. Ecker, June 5, 2012 (sent for review April 30, 2012)
Regulation of gene expression by DNA methylation is crucial for defining cellular identities and coordinating organism-wide developmental programs in many organisms. In plants, modulation of DNA methylation in response to environmental conditions represents a potentially robust mechanism to regulate gene expression networks; however, examples of dynamic DNA methylation are largely limited to gene imprinting. Here we report an unexpected role for DNA methylation in regulation of the Arabidopsis thaliana immune system. Profiling the DNA methylomes of plants exposed to bacterial pathogen, avirulent bacteria, or salicylic acid (SA) hormone revealed numerous stress-induced differentially methylated regions, many of which were intimately associated with differentially expressed genes. In response to SA, transposon-associated differentially methylated regions, which were accompanied by up-regulation of 21-nt siRNAs, were often coupled to transcriptional changes of the transposon and/or the proximal gene. Thus, dynamic DNA methylation changes within repetitive sequences or transposons can regulate neighboring genes in response to SA stress.
Cytosine methylation is a stable and heritable modification of the DNA that imparts epigenetic control throughout the genome, including regulation of coding and noncoding elements. In many eukaryotes, regulation of gene expression by DNA methylation is crucial for defining cellular identities and coordinating organism-wide developmental programs. At the molecular level, DNA methylation is coupled to nucleosome positioning (1), specific histone modifications (2), and transcriptional activity (3, 4). For example, in Arabidopsis thaliana, simultaneous disruption of the chromatin remodeling enzymes KYP, SUVH5, and SUVH6 results in a concomitant decrease in cytosine methylation and H3K9me2 levels, and consequently transcriptional reactivation of heterochromatic transposons (5).
In Arabidopsis, DNA methylation is deposited at CG, CHG, and CHH sequences (where H is A, C, or T) through three genetically separable pathways to regulate transposon silencing (6), genomic imprinting (7⇓–9), and stable gene silencing (3, 4, 10, 11). Cytosine methylation is established in all sequence contexts by de novo methyltransferases (DRM1/2) through a small RNA-directed DNA methylation (RdDM) pathway (12). Here DICER-dependent 21- to 24-nt siRNAs guide Argonaute proteins (AGO4/AGO6) to complementary sequences within the genome, likely through a siRNA:nascent RNA base pairing mechanism, to direct cytosine methylation (13, 14). Methylation of CGs and CHGs are maintained through DNA replication by MET1, a homologue of the mammalian DNA methyltransferase DNMT1, and the plant-specific CMT3 methyltransferase, respectively (15⇓–17). Conversely, active demethylation of methylcytosines is catalyzed by the DEMETER (DME) family of DNA glycosylases (18⇓–20). How these proteins and others collaborate to shape the epigenetic landscape has been intensely studied in recent years; however, the potential for these pathways to be dynamically regulated during nondevelopmental processes has yet to be thoroughly investigated.
In the seed, hypermethylation of transposable elements (TEs) in the embryo genome, along with pervasive demethylation of the embryo-nourishing endosperm genome by DME, function to silence TEs and coordinate gene imprinting, respectively, and are hallmarks of early development (7, 8). For example, proper gametophyte development requires temporal expression of the maternally imprinted FWA gene in the endosperm, a locus that is dynamically regulated by active demethylation of methylcytosines within upstream repetitive sequences (9). Although it remains unclear if widespread alterations in DNA methylation analogous to gene imprinting can be elicited by stress, plants subjected to heat stress display transient changes in nucleosome density, as well as transcriptional derepression, at some repetitive elements (21, 22), indicating that the epigenetic landscape can be dynamically modified.
Here, we demonstrate that DNA methylation imparts control over the Arabidopsis defense response against the biotrophic pathogen Pseudomonas syringae pv. tomato DC3000 (Pst). To expand on these observations, we generated genome-wide, single-base resolution DNA methylomes of plants exposed to different biotic stresses. Early studies examining DNA methylation alterations in response to stress have been limited by low-resolution, nonquantitative, or noncomprehensive approaches. More recently, immunoprecipitation of methylated DNA followed by a partial-genome analysis with tiling arrays has provided increased resolution of DNA methylation changes in the progeny of abiotically stressed plants (23, 24). Although these initial studies have supported the potential for localized methylation changes in response to some abiotic stresses, a high-resolution quantitative analysis of DNA methylation alterations in response to different biotic stresses would provide unprecedented insight into the dynamics of cytosine methylation during infection. Our unbiased, genome-wide approach uncovered unique aspects of stress-induced dynamic DNA methylation changes, including a striking relationship among hypomethylation, biogenesis of specific siRNAs, and transcriptional derepression at some transposons.
Mutants defective in CG or non-CG methylation exhibit numerous developmental pleiotropies arising from organism-wide gene derepression (15, 25) and transposon mobilization (6) (stable or heritable methylation), as well as endosperm-specific defects in gene imprinting (dynamic methylation). To investigate if DNA methylation regulates nondevelopmental processes, we exposed mutant plants globally defective in maintenance of CG methylation (met1-3) or non-CG methylation (ddc, drm1-2 drm2-2 cmt3-11) to the bacterial pathogen Pst. Surprisingly, both mutants were markedly resistant to bacterial colonization and failed to develop the characteristic chlorotic phenotype associated with pathogen infection (Fig. 1 A and B). Similarly, mutants partially impaired in CG (ddm1) or non-CG (rdr1, rdr2, rdr6, drd1, nrpd1a, and dcl2/3/4) methylation displayed modest increases in Pst resistance (SI Appendix, Fig. S1). Bacterial growth of avirulent or nonpathogenic P. syringae strains was also restricted in the met1-3 and ddc mutants (SI Appendix, Fig. S2), indicating that loss of DNA methylation nonspecifically enhances resistance to bacteria.
Mutants defective in DNA methylation are resistant to pathogenic P. syringae. Adult plants were infected with Pst (1 × 105 cfu·mL−1). At the indicated time points, leaf tissue was harvested and (A) representative leaves were photographed or (B) the bacteria were extracted, the colony forming units were quantified, and the data were plotted as the mean ± SEM. (C) A heat map representation of a one-dimensional hierarchical clustering of genome-wide differential expression levels as determined by mRNA-seq for the mutants relative to their WT Col-0 controls during untreated and Pst-infected (5 dpi) conditions.
To investigate the transcriptional consequence of widespread hypomethylation on pathogen-responsive genes, we performed mRNA-seq on untreated or Pst-infected met1-3 and ddc plants. Consistent with previous studies (3), steady-state transcript levels for both mutants were generally similar to WT levels, although we did identify a number of genes that were constitutively misexpressed (>10-fold change, met1-3, 2,034 genes or 8.5%; ddc, 1,286 genes or 5.4%; Fig. 1C and SI Appendix, Fig. S3). More strikingly, exposure of both mutants to Pst resulted in up- and down-regulation of numerous pathogen-responsive genes relative to the infected WT control, indicating that Pst-induced transcriptional networks interface with DNA methylation to regulate gene expression (Fig. 1C). Constitutive and inducible misregulation of gene expression is exemplified at the pathogen-responsive gene PR1 (SI Appendix, Fig. S4), which is not directly targeted by DNA methylation, but rather is misexpressed as a consequence of methylation-dependent alterations in transcriptional networks. Intriguingly, we also observed that met1-3 plants display a fully penetrant dwarf phenotype in the F2 generation (Fig. 1A), as well as constitutively overexpress PR1, defects that may be attributable to constitutive overactivation of plant defenses (26⇓–28). Together, these data are consistent with the observation that global demethylation of the rice genome with 5-azadeoxycytidine enhances bacterial resistance to virulent strains of Xanthomonas (29), as well as supports the hypothesis that reprogramming the DNA methylation states at some loci is an important mechanism for plant defense.
To further investigate the possibility that DNA methylation dynamically responds to biotic stress, we performed genome-wide methylation profiling (MethylC-seq) (3), on leaf tissue from populations of Pst-infected plants at 5 days postinfection (dpi). By using a rigorous purification strategy to enrich for Arabidopsis nuclear DNA (>99.8% bacterial-free), we generated two independent biological replicates of high-coverage (95–98% of mappable cytosines; 93–95% of all cytosines), high-depth (8.2–12.1× per cytosine) DNA methylomes for untreated and Pst-treated plants (SI Appendix, Table S1). A pairwise comparison of the methylation levels at 4,046,648 high-confidence methylcytosines (mCs; SI Appendix) revealed that all four individual methylomes were well correlated (Pearson r = 0.83–0.88), as were the replicate means (Pearson r = 0.92; SI Appendix, Table S2). Additionally, the distribution of genome-wide methylation levels at mCGs, mCHGs, mCHHs were similar across samples (SI Appendix, Fig. S5), further indicating that complete, pancellular reprogramming of DNA methylation does not occur after 5 d of Pst exposure.
To examine the possibility that dynamic regulation of DNA methylation occurs locally in response to stress, we applied a linear mixed-effect model to identify differentially methylated regions (DMRs; Fig. 2A and SI Appendix, Figs. S6–S9). By using this approach, we found that differentially methylated cytosines (DmCs) were enriched in gene-rich and depleted in gene-poor regions of the genome (Fig. 2B and SI Appendix, Fig. S7), suggesting a role for Pst-induced methylation changes in transcriptional control. Therefore, we focused our subsequent analysis on the 8,360 DmCs, or 0.39% of the methylated cytosines, in gene-rich regions (1,513 DMRs; SI Appendix, Fig. S7 and Table S3). Moreover, we observed similar levels of differential methylation at these DMRs at an earlier 3 dpi time point (Fig. 2C), suggesting that targeted methylation changes persist between different infection time points (i.e., 3–5 dpi). Generally, Pst-induced DMRs were small in size (>90% were eight DmCs or less; SI Appendix, Fig. S9), and largely composed of cytosines in the CG and CHH sequence context (Fig. 2D), whereas DmCHGs were heavily underrepresented (P = 9.8 × 10−324). This was exemplified at the At3g50480 locus (Fig. 2A), which encodes a homologue of RPW8, a canonical disease-resistance protein that confers resistance to powdery mildew (30). Finally, we observed hypo- and hypermethylation of mCs in every sequence context (Fig. 2E), indicating that active remethylation or demethylation of cytosines can occur dynamically in response to pathogen. Importantly, methylation changes at most DMRs were relatively modest (SI Appendix, Table S3), likely because the dynamic methylation response occurs in a localized or cell type-specific manner, which may in turn collaborate with other cell-autonomous or -nonautonomous transcriptional programs to control bacterial growth.
Characterization of genome-wide DNA methylation changes in response to bacterial pathogen. (A) AnnoJ genome browser depicts a region of demethylation induced by Pst infection. (B) Enrichment of DmCs in gene-rich or gene-poor regions relative to mCs. (C) A dendrogram comparing the untreated, 3 dpi, and 5 dpi samples based on Pearson correlation of the methylation levels (i.e., mC/C) at Pst-induced DmCs. (D) The sequence context breakdown for DmCs in gene-rich regions and (E) the fraction that are comprised of hypomethylated (hypo-DmC) or hypermethylated (hyper-DmC) cytosines. (F) Enrichment of DmCGs and DmCHHs at the indicated genome features or (G) in 300-bp bins upstream (−1.5 kb) and downstream (+1.5 kb), or five equally sized bins within (green) all protein-coding genes. TAIR9 genome annotations were used for all analyses. P values were derived from a hypergeometric test between the number of DmCs and mCs within each genome feature (F) or bin (G).
To determine if certain genomic regions are preferentially targeted for differential methylation during stress, we examined several types of genomic and genic features for DmC enrichment, focusing exclusively on DmCGs and DmCHHs. Consistent with a function in proximal cis-acting gene regulation, we observed a strong enrichment of DmCs in intergenic regions (Fig. 2F), as well as a peak of DmCG enrichment ∼1 kb upstream of the transcriptional start site of protein-coding genes (Fig. 2G). Furthermore, we found a considerable amount of DmCHH enrichment flanking protein-coding genes (Fig. 2G), much of which is a consequence of DmCHH enrichment within TEs (Fig. 2F and SI Appendix, Fig. S10), suggesting that dynamic regulation of methylation within transposons may contribute to short-range gene regulation. Importantly, DmCs were also enriched at the 5′ and 3′ ends of protein-coding genes (Fig. 2 F and G), consistent with the observation that demethylation of imprinted genes in the endosperm preferentially occurs at the 5′ and 3′ ends (7, 8), indicating that pathogen-induced dynamic methylation may use similar mechanisms. Notably, loss of the DME family of demethylases results in preferential hypermethylation of mCs flanking gene bodies, including upstream promoter regions and 3′ UTRs (3), further suggesting that active modulation of methylation at gene boundaries contributes to gene regulation.
In Arabidopsis, Pst infection elicits a basal defense response that, in turn, is suppressed by an arsenal of bacterial virulence factors to enhance colonization; however, the plant immune system can recognize some virulence factors through a sophisticated surveillance system to robustly activate defense and hormonal pathways, including up-regulation of salicylic acid (SA) signaling, to render the bacterium avirulent (31). To examine dynamic regulation of cytosine methylation during different types of stress conditions, we generated two replicates of high-coverage (93–98% of mappable cytosines; 90–94% of all cytosines), high-depth (8.7–16.3× per cytosine) DNA methylomes for plants treated with avirulent bacteria [Pst(avrPphB), 5 dpi] or exogenous SA (SI Appendix, Table S1). Although methylation levels were globally correlated across all the DNA methylomes generated (Pearson r = 0.83–0.88; SI Appendix, Table S2), a hierarchical clustering analysis of the methylation levels at Pst-induced DmC positions revealed that methylation at CGs and CHGs were similarly altered in plants exposed to SA or avirulent or virulent Pst strains (Fig. 3A). In contrast, methylation levels at Pst-responsive DmCHHs were unique to Pst infection, suggesting that dynamic regulation of CHH methylation, but not necessarily CG or CHG methylation, may respond differently to a particular stress; or, alternatively, nearby DmCHHs within the same RdDM-targeted region may impart similar transcriptional control during different stress conditions.
Profiling DNA methylation during different stress conditions reveals unique aspects of dynamic differential methylation. (A) A dendrogram comparing untreated (red) or stress-treated (black) plants based on Pearson correlation of the methylation levels (i.e., mC/C) at the Pst-induced DmC positions. (B) The total number and sequence context breakdown of DmCs in response to each stress. (C) The number, direction, and degree of methylation changes at DmCGs, DmCHGs, and DmCHHs for each treatment (Hypo, hypomethylation, blue; Hyper, hypermethylation, red). (D) The average differential methylation at mCGs, mCHGs, mCHHs within 1-kb sliding windows of chromosome 1. The signal-to-noise (STN) ratio is plotted as μtreated − μuntreated/σtreated + σuntreated, where μ and σ are the mean and SD of the replicate mC/C values within a window, respectively. (E) A comparison of the differential methylation values (treated − untreated) at gene-rich DmCGs or DmCHHs that overlap between the indicated stresses (scatter-plot), and the number of observed (“O”) or expected (“E”) overlapping DmCs for each comparison (bar graph).
To further compare methylation levels between differently stressed plants, we identified DMRs in the methylomes of plants treated with SA and avirulent bacteria (Fig. 3B and SI Appendix, Table S3). Plants exposed to virulent or avirulent bacteria exhibited similar numbers of DmCGs and DmCHGs; however, SA induced substantially more DmCs in both sequence contexts (4.4–7.4 fold; Fig. 3 B and C), most likely a result of the high concentrations of hormone applied to the plants and pancellular induction of stress responses. The number of DmCHHs, in contrast, varied across samples, providing additional support for differential regulation of CHH dynamic methylation pathways in response to distinct stresses. Unexpectedly, we observed a preponderance of hypomethylated (77%) or hypermethylated (89%) DmCs in response to SA or avirulent bacteria, respectively (Fig. 3C). This prompted us to investigate if similar, but more widespread, differential methylation can be seen within these methylomes by using a low-resolution (1-kb bins) approach rather than the single-nucleotide pairwise comparisons. Indeed, we detected megabase-scale, sequence context-specific hypomethylation (with SA) or hypermethylation (with avirulent bacteria) of pericentromeric regions that was consistent with the levels of differential methylation at DmCs in gene-rich regions (Fig. 3D). Demethylation of centromeric repeats, as well as decondensation of heterochromatin, has been previously observed 24 h after Pst infection (32); however, at 5 dpi, we were unable to detect widespread pericentromeric hypomethylation in response to Pst. It is possible that methylation levels in these heterochromatic regions are highly sensitive to stress: an acute or persistent stress may initiate robust demethylation [e.g., 24 h after Pst infection (32) or chronic SA exposure], which undergoes remethylation as the stress condition is removed (e.g., 5 dpi Pst) or transitions into a protective hypermethylated state to guard against future stress (e.g., 5 dpi avirulent Pst), as has been described in the progeny of salt-stressed plants (24).
Despite these differences at heterochromatic regions, we observed numerous DmCGs and DmCHHs that were similarly targeted for differential methylation during multiple stress conditions (Fig. 3E and SI Appendix, Fig. S11). Furthermore, the levels of differential methylation at these DmCs were well correlated (Pearson r = 0.71–0.88), suggesting that a common mechanism may be responsible for modifying these cytosines during a variety of stressful conditions. This phenomenon was particularly evident at At1g13470, an uncharacterized gene that is normally repressed by RdDM within the gene body in a manner analogous to transposon silencing, yet exhibits robust stress-induced demethylation, as well as a concomitant increase in gene expression (SI Appendix, Fig. S12).
To further expand on these observations, we examined the short-range relationship between DMRs and transcript abundance on a genome-wide scale (strand-specific mRNA-seq; SI Appendix, Table S4). First, we assigned every DMR to a proximal protein-coding gene based on the relative positioning of each feature (SI Appendix). A Gene Ontology analysis (33) of DMR-associated genes in response to Pst or SA revealed a strong enrichment of genes that function in plant defense (Fig. 4A). Infection with avirulent bacteria, in contrast, induced DMRs proximal to genes largely involved in transcriptional regulation and, to a lesser extent, the defense response (Fig. 4A and SI Appendix, Table S5), indicating that DMRs at late stages of avirulent Pst infection (5 dpi) likely influence different aspects of the plant defense system. Next, we determined the transcript levels of all annotated protein-encoding genes (34) for each stress treatment (SI Appendix, Fig. S13). Differential expression of DMR-associated genes was dependent on the direction of methylation change, as genes proximal to hypomethylated DMRs exhibited higher levels of differential transcript abundance relative to all genes (Wilcoxon P = 7.4 × 10−4 to 2.4 × 10−7), whereas genes near hypermethylated regions displayed similar, or only modestly higher (avirulent Pst), levels of differential transcript abundance compared with all genes (Wilcoxon P = 0.40–2.0 × 10−3; Fig. 4B). Furthermore, an enrichment analysis of DMR-associated differentially expressed genes (DEGs) (35) revealed a strong correlation between gene body demethylation and increased transcript abundance, as well as a weaker relationship between hypomethylation of nongenic regions and up-regulation of proximal genes (Fig. 4C). Finally, we determined the transcript abundance of all of the stress-induced hypomethylated DMR-associated genes in the met1-3 and ddc mutants. In the methylation mutants, these genes were highly misexpressed relative to all genes (Fig. 4D), and many were misregulated in the absence or presence of pathogen (SI Appendix, Fig. S14). Moreover, a significant number of these genes were similarly misexpressed in both mutants, implicating CG and non-CG methylation in their regulation (SI Appendix, Fig. S14). Together, these data indicate that DNA methylation is at least partially responsible for transcriptional control of these genes. To our knowledge, these data represent the first widespread identification of stress-responsive genes whose transcript abundance is coupled to dynamic changes in DNA methylation (virulent, 148 genes; avirulent, 380; SA, 753; SI Appendix, Table S3).
Induction of differential methylation is coupled to differential gene expression. (A) Gene Ontology enrichment analysis of DMR-associated genes (SI Appendix, Table S5). (B) Differential expression levels of all genes (white) or genes associated with hypomethylated (blue) or hypermethylated (red) DMRs are displayed as box-plots (boxes represent the quartiles, whiskers mark data within 1.5 interquartile ranges of the quartile, and outliers are suppressed; Wilcoxon P values are reported). (C) Enrichment of DEGs proximal to hypomethylated regions (Hypo-DMRs) positioned within protein-coding genes (genic) or all other regions (nongenic). (D) Histograms showing that hypomethylated DMR-associated genes are more misexpressed in the DNA methylation mutants than expected (relative to all genes). The enrichment is calculated independently for each differential gene expression range. (C and D) A hypergeometric test was used to calculate enrichment P values (*P < 0.05, **P < 1 × 10−5, and ***P < 1 × 10−10).
Disruption of DNA methylation pathways results in genome-wide changes in small RNA (smRNA) abundance (3). Thus, we sequenced smRNAs from untreated and SA-challenged plants and found that globally smRNA levels remained constant upon SA treatment; however, we observed a marked increase in 21-nt, but not 24-nt, siRNAs specifically at transposons (Fig. 5A). This phenomenon was particularly evident at TE-associated DMRs (Fig. 5B), suggesting that dynamic methylation may function with 21-nt siRNAs to shape the epigenetic landscape at transposons in response to SA. Curiously, loss of MET1 or down-regulation of the chromatin remodeling ATPase DDM1 results in transcriptional activation of some transposons, as well as an increase in 21-nt siRNAs (3, 36). Notably, MET1 and DDM1 transcripts were both down-regulated in response to biotic or SA stress (SI Appendix, Fig. S15), further implicating these proteins and 21-nt siRNAs in TE regulation.
Transposon-associated DMRs are associated with siRNA and transcriptional changes. (A) SA-induced siRNA levels of the indicated sizes were calculated at the indicated genome features (31,189 TEs) as a ratio relative to untreated levels. Reads corresponding to miRNAs were discarded. (B) Normalized 21- and 24-nt siRNA levels in 100-bp nonoverlapping bins upstream (−1 kb) and downstream (+1 kb) of TE-associated DMRs. (C) Enrichment of DETs at TE-associated hypo- or hypermethylated regions. A hypergeometric test was used to calculate the P value for each type of association. (D) Top: Enrichment of transposons (left y-axis) in 250-bp nonoverlapping bins upstream (−5 kb) and downstream (+5 kb) of the transcriptional start sites (TSS) of nongenic DMR-associated DEGs or non-DEGs calculated relative to all DEGs or all non-DEGs, respectively. The broken gray line indicates the normalized gene density as a function of length (right y-axis). Bottom: Average number of TEs per gene (within 5 kb of the TSS) for the indicated gene classes. (E) Heat map representation of the differential methylation (DM) levels and differential expression (DE) of 21-nt siRNAs across all DETs, as well as the differential expression levels of neighboring DEGs (<5 kb between TSSs). TE-associated DMRs are also reported.
DNA methylation at TEs contributes to transcriptional silencing of transposons and proximal genes within the heterochromatin, as well as dynamic regulation of some imprinted euchromatic genes (37). To examine the possibility that methylation within TEs is targeted during stress to regulate gene expression, we first examined the transcript levels of DMR-associated transposons and found a strong relationship between transcriptional up-regulation of the TEs and demethylation of the underlying DNA (Fig. 5C). Notably, methylation-independent transcriptional activation of the ONSEN transposon is triggered by heat stress and results in coinduction of the proximal protein-coding gene (21). Therefore, we next determined the TE density at DMR-associated DEGs or non-DEGs relative to all DEGs or all non-DEGs, respectively. We identified peaks of TE enrichment ∼2 kb upstream of the transcriptional start site and near the 3′ end of these genes (Fig. 5D). The TE enrichment was substantially more pronounced for DMR-associated DEGs compared with non-DEGs (P = 2.9 × 10−43 at −2 kb), suggesting that transcription of these genes is coupled to the methylation state of the proximal TEs. Importantly, TE enrichment flanking DMR-associated genes was not simply an artifact of SA-induced DEGs or stress-responsive genes generally exhibiting enhanced transposon densities, as TE densities are similar across all gene classes (Fig. 5D, Bottom).
To further investigate whether dynamic alterations in methylation and/or 21-nt siRNAs at TEs drive expression changes of the transposon and the neighboring gene, we profiled SA-induced alterations in DNA methylation and 21-nt siRNAs at all differentially expressed transposons (DETs). We found that transcriptional up-regulation of TEs correlated with biogenesis of 21-nt siRNAs and demethylation of the underlying DNA, whereas transposon down-regulation did not show a similar correlation (Fig. 5E). In many cases, but not all, TE expression changes correlated with similar alterations in the expression of the neighboring protein-coding genes, suggesting that their transcriptional regulation may be coupled. Notably, this correlation was not dependent on 21-nt siRNA biogenesis, as several transcriptionally repressed TEs were associated with down-regulated DEGs. Together, our data suggest that methylation levels within TEs may dynamically control expression of the transposon and, in some cases, the proximal gene in response to stress, as exemplified at the At4g39860 locus (SI Appendix, Fig. S16).
Plants use a sophisticated series of defense mechanisms to restrict the growth of biotrophic bacteria upon infection, characterized by a coordinated induction of various hormonal signals that trigger alterations in gene expression networks, which, together, restrict bacterial growth. Here, we demonstrate that some defense genes are modulated by DNA methylation, as mutants globally lacking CG or non-CG methylation display constitutive and inducible misexpression of pathogen-responsive genes. These transcriptional alterations correlate with a dramatic enhanced Pst resistance phenotype in the met1-3 and ddc mutant plants. Interestingly, ddm1 mutants, which display CG methylation defects at repeats and transposons but not genes (37), only partially phenocopy the met1 mutants, suggesting that CG methylation represses plant defenses by regulating repeats/transposons and protein-coding genes.
Loss of many RdDM pathway components, including simultaneous deletion of DRM1/2 and CMT3, results in measurable differences in Pst growth upon infection. In contrast, global depletion of CHG methylation (cmt3-7) or CHH methylation (drm1-2 drm2-2) alone fails to impart Pst resistance (38); however, this may be a consequence of compensatory methylation by the remaining DNA methyltransferases in these mutant backgrounds (3, 39). Surprisingly, ago4 mutants display increased susceptibility to Pst even though loss of RNA polymerase V, which recruits AGO4 to chromatin and functions upstream of DRM2, results in enhanced resistance (38, 40, 41). These conflicting phenotypes may be an indication that AGO4 has targets outside of the canonical RdDM pathway. Although it remains unclear whether these pathogen phenotypes are a direct or indirect consequence of transcriptional changes at improperly methylated defense gene(s), our data clearly indicate that the defense response against Pst as a whole is negatively regulated by DNA methylation.
The striking phenotypes that we observed in the DNA methylation mutants prompted us to profile the methylome in response to biotic stress. Our genome-wide, high-resolution analysis of stress-induced methylation changes has provided unique insight into the malleability of DNA methylation in response to the environment. Notably, our approach provides only a snapshot of the methylation changes at a specific time point during infection. Thus, a detailed temporal analysis of the methylation dynamics from the onset of infection, for which our study provides a framework, will be necessary to fully understand dynamic methylation in the context of disease progression. Furthermore, it is important to note that our MethylC-seq approach is unable to uncouple cell-autonomous from non–cell-autonomous methylation dynamics or resolve what cells and/or cell types are responding to the stress. The subtle methylation changes that we observe suggest that alterations are likely to occur locally during infection, potentially in cells that are directly contacting bacteria. Moreover, the biotic stresses examined in our study are expected to trigger widespread cell death. It is unlikely that dynamic methylation is restricted to dying cells because it is often associated with active processes (e.g., increased gene expression); however, we cannot eliminate the possibility that cell death pathways directly or indirectly contribute to methylation alterations. Thus, the scope and specificity of stress-induced methylation changes is unlikely to be fully appreciated until the development of single-cell DNA methylation profiling.
Our data support a model whereby DNA methylation imparts persistent control over some defense genes during nonstressful conditions, but, in response to environmental stimuli, can change dynamically to alter gene expression. Importantly, our genome-wide observations are consistent with previous reports linking stress-induced methylation changes to transcriptional control of specific loci in other plant systems (29, 42⇓⇓⇓–46). It remains unclear, however, why some defense genes are held under the control of DNA methylation. Notably, hyperactivation of plant defenses through the SA signaling pathway results not only in pathogen resistance but also dwarfism (26⇓–28). Thus, DNA methylation may repress some SA pathway genes to avoid these developmental defects until their expression is required for defense. This is consistent with our observation that genes targeted by stress-induced dynamic demethylation tend to be constitutively misexpressed in the met1-3 and ddc mutants.
Although many of the DMRs we identified were associated with protein coding genes, there were also numerous transposons that were targeted by dynamic methylation in response to SA. Reprogramming methylation at TEs correlated with changes in their expression and/or the expression of proximal genes, as well as the biogenesis of TE-associated 21-nt siRNAs. Curiously, these smRNAs are indistinguishable from the mobile epigenetically activated siRNAs in pollen (36). In this case, TEs are activated in the pollen vegetative nucleus thereby triggering the biogenesis of epigenetically activated siRNAs, which in turn accumulate in the gametes, where they enforce transposon silencing. It is intriguing to speculate that SA-induced hypomethylation at TEs triggers 21-nt siRNA biogenesis that may function in a similar non–cell-autonomous manner to shape the epigenetic landscape of surrounding vegetative cells, or potentially, reproductive tissues. This could serve as a conduit for transgenerational memory, whereby enhanced stress resistance encoded in epigenetic information is conferred upon the progeny. Indeed, some initial studies have supported the potential for transgenerational memory of stress in several systems, including plants, nematodes, and mammals (23, 24, 46⇓⇓⇓–50). This process appears to engage DNA methylation and/or smRNA pathways to induce expression changes at specific loci in the progeny of stressed individuals. A careful experimental design that takes advantage of high-resolution DNA methylation profiling will be essential to determine if biotic stress in vegetative cells of Arabidopsis can trigger alterations in methylation in reproductive tissues, which, in turn, can be propagated into the progeny.
A. thaliana plants, ecotype Col-0, were grown in a Promix-HP:vermiculite (2:1) soil mix at 22 °C. The mutant plant lines met1-3 (17), met1-9 (51), ddm1-8 (52), ddm1-10 (53), drm1-2 drm2-2 cmt3-11 (54), ago4-2 (38), ago6-2 (13), rdr1-1 (55), rdr2-1 (55), rdr6-15 (56), drd1-5 (57), nrpd1a-4 (58), nrpd1b-11 (59), nrpd2a-2 nrpd2b-1 (60), and dcl2-1 dcl3-1 dcl4-2 (61) have been previously described. The pathogen strains used in this study were Pst, Pst(avrPphB) (62), and Pst(hrcC−) (63). All Pseudomonas strains were grown on King's B agar plates at 28 °C and antibiotic selection was carried out by using the following concentrations (in μg·mL−1): kanamycin, 50; rifampicin, 100.
All experiments were performed in adult WT Col-0, met1-3 (F2 generation), and drm1-2 drm2-2 cmt3-11 (ddc, F1 generation) plants (SI Appendix). Bacterial growth assays were performed in plants infected at 1 × 105 cfu·mL−1 (OD600, 0.0002) by vacuum infiltration. Bacterial isolation and quantification of the associated colony forming units was performed as previously described (64), and data were plotted as the mean and SE of the decimal logarithm [log(cfu·cm−2)] of approximately eight replicate plants. For qualitative measurement of disease symptoms, plants were infected with Pst at 1 × 105 cfu·mL−1 by vacuum infiltration, and representative photographs were taken of uninfected leaves or leaves at 1, 3, or 5 dpi. All infection experiments were performed at least twice.
A single Col-0, met1-3, or ddc seed stock was used to grow plants for all of the genomic studies to minimize the possibility of variability in DNA methylation levels that may have been induced in previous generations. Short-day–grown adult plants (∼5 wk old) were infected with virulent Pst or avirulent Pst(avrPphB) bacteria at 1 × 105 cfu·mL−1 by vacuum infiltration and maintained at high humidity, and leaf tissue was harvested at 5 dpi. A 3 dpi Pst sample was also generated (one biological replicate). For SA treatments, adult plants (∼5 wk old, short-day–grown) were sprayed with 1 mM SA (Fisher Scientific) containing 0.01% Silwet l-77 (Lehle Seeds) every day for five consecutive days, and tissue was collected on day 6 (i.e., 5 d of exposure). The untreated control plants were grown under identical conditions as the stressed plants; however, tissue was collected before any stress treatment. Two individual populations of plants, representing two biological replicates, were grown and treated independently for each condition. Preparation of nuclei from leaf tissue was performed as previously described (65, 66), with only minor modifications, as detailed in SI Appendix. Genomic DNA was extracted from purified nuclei using the Plant DNeasy kit (Qiagen) according to the manufacturer’s instructions, and MethylC-seq libraries were prepared as previously described (67).
Total RNA was isolated from ∼250 mg of frozen leaf powder using the mirVana miRNA Isolation kit (Ambion) or the Plant RNeasy kit (Qiagen) according to the manufacturer’s instructions for isolation of total RNA from plant tissue. Non–strand-specific mRNA-seq libraries were generated from 4 μg of total RNA (Plant RNeasy kit) and prepared by using the TruSeq RNA Sample Prep kit (Illumina) according to the manufacturer’s instructions. For the strand-specific mRNA-seq libraries, polyadenylated mRNA was purified from total RNA (40–80 mg; mirVana miRNA Isolation kit; Ambion) with two sequential poly(A) selections using the Oligotex mRNA Mini kit (Qiagen) per the manufacturer’s instructions. The RNA was ethanol-precipitated following each poly(A) selection. Sequencing libraries were prepared as described in the Directional mRNA-seq Sample Preparation kit protocol (Illumina), with modifications that are detailed in SI Appendix.
RNAs enriched for smRNAs were isolated from frozen leaf powder (∼200 mg) using the mirVana miRNA Isolation kit (Ambion) according to the manufacturer’s instructions. Following ethanol precipitation, the 5′ ends of the RNAs were dephosphorylated with RNA 5′ polyphosphatase (Epicentre) per the manufacturer’s instructions, and the smRNAs were resolved by electrophoresis on a 15% (wt/vol) TBE-urea gel. RNAs corresponding to ∼15 to 50 nt in length were excised, eluted from the gel, and ethanol-precipitated. smRNA-seq libraries (two biological replicates per condition) were generated by using the smRNA Sample Prep kit (version 1.5; Illumina) according to the manufacturer’s instructions.
Details of Illumina high-throughput sequencing, data processing, and statistical analyses are described in SI Appendix.
We thank Dr. Mathew Lewsey for critical review of this manuscript. This work was supported by a Catharina Foundation postdoctoral fellowship (to M.P.), National Institutes of Health National Research Service Award Postdoctoral Fellowship F32-HG004830 (to R.J.S.), a California Institute for Regenerative Medicine Training Grant (to R.L.), National Institutes of Health Grant AI060662 (to J.E.D.), the Mary K. Chapman Foundation (J.R.E.), Howard Hughes Medical Institute (J.R.E.), and Gordon and Betty Moore Foundation (J.R.E.). J.R.E. is a Howard Hughes Medical Institute/Gordon and Betty Moore Foundation Investigator.
↵1Present address: Department of Molecular Biology, Massachusetts General Hospital, Boston, MA 02114.
↵2Present address: Center for Genomic Science of IIT@SEMM, Istituto Italiano di Tecnologia at the IFOM-IEO Campus, 20139 Milan, Italy.
↵3Present address: Whitehead Institute for Biomedical Research, Cambridge, MA 02142.
Author contributions: R.H.D., J.E.D., and J.R.E. designed research; R.H.D., R.J.S., J.M.D., and J.R.N. performed research; R.H.D., M.P., R.J.S., and R.L. analyzed data; and R.H.D. and J.R.E. wrote the paper.
The authors declare no conflict of interest.
Data deposition: The sequences reported in this paper have been deposited in the National Center for Biotechnology Information Sequence Read Archive (accession no. SRA036590).
See Commentary on page 12842.
See Author Summary on page 12858 (volume 109, number 32).
This article contains supporting information online at www.pnas.org/lookup/suppl/doi:10.1073/pnas.1209329109/-/DCSupplemental.
Freely available online through the PNAS open access option.
