Gene Body Methylation in Plants: Mechanisms, Functions, and Important Implications for Understanding Evolutionary Processes

1Gene body methylation (hereafter gbM) is a common phenomenon in plants and can affect up to 60% of genes in some species. It has been controversial whether gbM has any function in plants, but recent findings suggest it is under selection and correlated with fitness. Here, we review the scientific literature, include novel analyses, and discuss the potential role of gbM in rapid evolutionary change.

2Epigenetics is the study of changes in gene expression that can be inherited through cell divisions (either mitotic or meiotic) that are not due to modifications in the DNA sequence. A longstanding question is whether epigenetics can play a role in adaptation. stated that “a direct demonstration that other molecules, in addition to DNA [sequence], carry substantial heritable information would represent an important conceptual change in evolutionary biology.” Theoretically, epigenetic marks may be the basis for this conceptual change. If epigenetic marks affect fitness, if they are inherited through generations, and if they epimutate over time, then they can be the target of selection and facilitate adaptation.

3Cytosine methylation is one common epigenetic mark that is generally found in eukaryotes, including vertebrates, insects, fungi, and plants. In some of these groups, cytosines are methylated only in a single context when they are part of a CG dinucleotide. In plants, however, cytosine methylation occurs in three sequence contexts—CG, CHG, and CHH (where H stands for A, T, or C). Methylation marks in these three contexts are produced by different biochemical pathways and have different patterns of inheritance. For example, epimutation accumulation lines in Arabidopsis thaliana have demonstrated that genome-wide methylation divergence at CG dinucleotides increases throughout >30 generations, illustrating that plant CG DNA methylation is transmitted from generation to generation and epimutates over time. In contrast, CHH methylation is mostly erased by demethylation in the A. thaliana male germline and later reset during embryonic development. Therefore, CHH methylation is only transmitted partially over, at most, one or a few generations (with the interesting exception of some asexual plants without meiosis; ). The transgenerational inheritance of the third context—CHG methylation—remains unclear. Although CHG methylation is retained during gametogenesis, epimutation accumulation lines in A. thaliana do not diverge for CHG methylation over generations, suggesting that CHG methylation is not inherited at a genome-wide scale. It is possible, however, that some genomic sites inherit CHG methylation over a few generations, especially in some asexual species. To summarize, of the three methylation contexts in plants, methylation in CG dinucleotides is most prone to transgenerational inheritance and is therefore the best candidate for epigenetic adaptation.

4To consider the possibility of epigenetic adaptation, it is also important to know where these marks reside in the genome. In flowering plants, patterns of DNA methylation vary among genomic regions. Methylation in all three contexts silences transposable elements (TEs) and prevents activity at regulatory elements. Both CHG and CHH genic methylation types are associated with reduced expression levels, as is CG methylation in promoter regions. In contrast, the exons of some genes (∼20% of A. thaliana genes; ) are methylated only in the CG context, a phenomenon called gene body methylation (hereafter gbM). gbM is mostly found in moderately and constitutively expressed housekeeping genes. However, since its initial discovery, the topic of gbM function has been controversial. If it has no function, it is obviously unlikely to contribute to adaptive processes.

5Here, we review our current knowledge of gbM in plants, with the ultimate goal to critically evaluate whether it has a function and may be a target for natural selection. We start by presenting the mechanisms of plant gbM establishment and maintenance because these mechanisms are crucial for understanding how selection could act on this epigenetic state. We then consider the evolution of gbM, specifically whether gbM is a neutral manifestation of epigenomic dynamics or whether there is evidence that it can be advantageous. Adaptive arguments presume that gbM has a phenotype on which natural selection can act. Some, but certainly not all, recent work has established a connection between gbM and gene expression, but questions about generality and mechanisms remain. To address these questions, we review functional analyses of mutants and also comparative epigenomic approaches that have studied hypothetical functions of gbM, namely, its potential role in regulating and stabilizing expression, preventing aberrant transcription, and improving the fidelity of intron splicing. Finally, we present a model synthesizing the prevalence, distribution, and effect of gbM with its potential evolutionary significance.

6CG methylation is maintained during plant cell division by Methyltransferase 1 (MET1), which adds a methyl group on the symmetrical CG dinucleotide of a complementary DNA strand. Epimutation accumulation lines in A. thaliana have shown that the maintenance of CG methylation by MET1 is an inherently error-prone process, with the epimutation rate estimated to be ∼10−3 per generation per haploid epigenome for the loss of CG methylation in genes. Without methylation maintenance mechanisms, CG methylation is quickly diluted and lost over cell divisions, as demonstrated by the absence of CG methylation in A. thaliana met1 mutants.

7Studies on Eutrema salsugineum, a close relative of A. thaliana, have recently clarified the mechanisms responsible for the establishment of gbM in plants. Eutrema salsugineum lacks both gbM and the Chromomethylase 3 (CMT3) gene. The link between CMT3 loss and the absence of gbM was at first enigmatic because, until recently, CMT3 was not known to methylate CG sites. It is now known that the CMT3 protein is involved in a self-reinforcing feedback loop: CMT3 recognizes the histone mark H3K9me2 (histone H3 lysine 9 dimethylation) and then de novo methylates nearby cytosines predominantly in the CHG context but also occasionally in the CG context. CHG DNA methylation in turn leads to H3K9 methylation by SU(VAR) Homologue 4 (SUVH4), leading to a positive feedback loop between CHG methylation and H3K9 methylation. CHG methylation typically suppresses transcription; however, in A. thaliana, CHG methylation is removed in transcribed genes due to active demethylation of H3K9 by Increased in Bonsai Methylation 1 (IBM1).

8The joint loss of CMT3 and gbM evolved independently in two Brassicaceae species, corroborating their association. However, cmt3 mutants in A. thaliana have shown that CMT3 does not affect the maintenance of gbM once it is established, suggesting the action of CMT3 is limited to gbM establishment. Interestingly, transgenic reinsertion of CMT3 into E. salsugineum re-established genic methylation in all three contexts in a subset of genes. This subset of genes has been called “CHG-gain” genes, and these genes tend to be orthologous to gbM genes in A. thaliana. After the CMT3 transgene was lost, CHG-gain genes only maintained methylation in the CG context, presumably due to maintenance by MET1. On average, CHG-gain genes are longer, contain more exons, and exhibit a moderate—but on average higher—level of expression than non-CHG-gain genes. CHG-gain genes are also enriched for CWG trinucleotides (CAG and CTG) as opposed to CCG trinucleotides, consistent with the preferred substrate of CMT3. Finally, CHG-gain genes have a higher frequency of CHG cytosines compared to non-CHG-gain genes.

9The role of CMT3 in genic de novo methylation was recently confirmed in A. thaliana mutants that hyperexpress CMT3 during late embryonic development. CMT3 hyperexpression induces embryonic hypermethylation predominantly in the CWG context, but hypermethylation is also found in other contexts, including CG dinucleotides. These findings confirm that CMT3 is sloppy and can methylate contexts other than CHG. Methylation changes caused by embryonic CMT3 hyperexpression were maintained over cell divisions and still observed in 3-week-old plants, consistent with the model that CMT3-induced epimutations give rise to gbM that can be maintained by MET1 across cell divisions and generations. The same gene patterns were repeatedly observed in independent transgenic lines, confirming that CMT3 hypermethylation is not stochastic and tends to target a specific gene set. CMT3-induced hypermethylation was also enriched in genes characterized by inaccessible chromatin marks and heterochromatin histone variants. Altogether, these observations lead to a model in which gbM establishment is caused by the recruitment of CMT3, the formation of a feedback loop that ultimately produces CHG, CHH, and CG methylation, the eventual removal of CHG and CHH methylation, and the maintenance of the remaining CG methylation by MET1 (fig. 1).

10gbM genes are typically defined statistically as being significantly more methylated than the genic average in the CG context and significantly less methylated than the genic average in the CHG and CHH contexts. Once defined, the proportion of gbM genes varies greatly across species, with as many as ∼60% of genes in Mimulus guttatus but 0% in Marchantia polymorpha, Physcomitrella patens, and E. salsugineum. The lack of gbM in a few species has been used to argue that gbM is dispensable and thus has no function. However, the loss of gbM in a few species does not imply that it is nonfunctional in all plants.

11A remarkable feature of gbM is that it is enriched over a conserved set of orthologs among species as distantly related as ferns and angiosperms. Two alternative hypotheses can explain the remarkable conservation of gbM. The first is the biased establishment of gbM in a subset of specific genes with inaccessible chromatin marks and heterochromatic (H3K9me2) histone variants. If these biases are conserved across species, they could explain the distribution of gbM across both genes and species. This first hypothesis is neutral with respect to selection because it does not assume that gbM has any effect on fitness. Instead, in this scenario, gbM is a consequence of CMT3 activity that is retained and transmitted over generations by MET1.

12An observation in favor of the neutral hypothesis is that gbM genes share many characteristics of the CHG-gain genes described previously. That is, the genes targeted by CMT3 are like gbM genes, in that they are generally characterized as being constitutively expressed at moderate levels and tend to be longer than unmethylated genes, with more exons and a higher frequency of CAG and CTG (as opposed to CCG) trinucleotides. Moreover, CHG-gain genes in E. salsugineum tend to be orthologous to gbM genes in A. thaliana. This observation suggests gbM establishment is biased toward specific genes, potentially explaining the conservation of gbM between orthologs. The overlap between CHG-gain and gbM is not complete because ∼40% of CHG-gain genes in E. salsugineum are orthologous to a gbM gene in A. thaliana. One explanation for the imperfect overlap between CHG-gain genes and gbM genes is that they were defined in different species – that is, CHG-gain genes were defined in E. salsugineum and gbM genes in A. thaliana. Moreover, these two species diverged 47 Ma, which may be ample time for the targets of CMT3 to diverge. Finally, another plausible explanation is temporal. The CHG-gain genes in E. salsugineum were established experimentally over only a few generations; the continuation of this experiment over a much longer timeframe could lead to the establishment of methylation within more genes, potentially increasing the 40% overlap.

13An alternative hypothesis is that the conservation of gbM across genes and species is shaped in part by the action of natural selection. Under this scenario, specific subsets of genes are targeted for de novo establishment of gbM, but selection on or against gbM removes or maintains CG methylation in different gene sets. At least three observations support the hypothesis that some gbM is under selection. First, DNA methylation is mutagenic and elevates C to T substitutions. Therefore, the conservation of gbM in a specific set of orthologous genes is surprising, especially because gbM genes are generally enriched for housekeeping and other important functions and evolve more slowly than unmethylated genes. This suggests the possibility that the mutagenic nature of methylation is compensated by an advantageous effect that maintains gbM in specific genes. The second observation in favor of the selective hypothesis comes from the comparison of the gbM status in orthologous genes of eight grass species, where shifts in the gbM status of genes are almost exclusive to the tips of the phylogeny (i.e., in a single species). This pattern suggests that shifts in gbM are deleterious and generally not favored over evolutionary time (table 1); however, it is also possible that the pattern is driven by epimutational biases.

14The third observation is based on population genetic analyses because selection acting on gbM can be explicitly measured using DNA methylation variation among natural populations of a species. Indeed, if gbM is advantageous within a given gene, an unmethylated allele will be disadvantageous and removed by selection. In such a gene, only a small proportion of individuals should be observed with an unmethylated allele. To infer the intensity of selection, constructed a population model that relies on the site frequency spectrum (SFS) of epigenetic states. In an inspired analysis, applied this model to the SFS of CG sites within all genes of a sample of 92 A. thaliana individuals. They did not detect a deviation from neutrality (table 1), but this result came with two important caveats. The first is that the test is unlikely to be powerful with a small sample, particularly if gbM has a small impact on fitness. used the sample of 92 individuals that was available at the time, but larger samples now exist. The second is that CG methylation within genes is not limited to gbM genes but can also be found in genes that are methylated in all three contexts (i.e., TE-like methylation; ). Methylation in all three contexts within a gene can be caused by a nearby TE insertion, is known to suppress expression, and may be an indication of pseudogenization. Therefore, most genes with TE-like methylation are likely to be under different evolutionary pressures than gbM genes, such that analyzing both gbM and TE-like methylated genes together, as done by, is likely to confound opposing selection pressures.

15For all these reasons, we recently repeated the analyses of with the important difference that we separated gbM gene sets from any genes with TE-like methylation. We also relied on larger data sets—that is, two distinct subsets of 876 and 120 individuals that originated from different sources—from the 1001 methylomes project in A. thaliana. To assess whether selection acts on the gbM state, we characterized the population frequency of methylation at the gene level to estimate the SFS of gene allelic states. Using the population genetic model of, we inferred that genes with ancestral gbM in Brassicaceae were under significant selection to remain CG methylated in A. thaliana (table 1; ) based on the larger data set. Conversely, ancestrally unmethylated genes in Brassicaceae were under selection to remain unmethylated in A. thaliana. We repeated the analyses on the smaller data set and also on an SFS drawn at the level of individual cytosines. The former had similar trends as the larger data set but without a significant effect of gbM, and the latter corroborated our gene-level analyses. That is, the overall impression is that CG sites within ancestrally gbM genes in Brassicaceae have been under selection to remain methylated in A. thaliana, while CG sites within ancestrally unmethylated genes have been under selection to be unmethylated in A. thaliana. Importantly, the results were also confirmed after splitting the gene sets into CHG-gain and non-CHG-gain genes, as characterized by, showing that biases in epimutation rates between gene sets were not completely responsible for the inferred selection acting on gbM. In other words, this control using CHG-gain genes shows that cis effects (either genetic or epigenetic) that locally influence epimutation rates do not explain the inferred selective pressures.

16Like all evolutionary analyses, there are caveats to this analysis, too. First, it relies on a model that simplifies the evolutionary process and includes assumptions that do not strictly fit the study organism (e.g., the model assumes random mating, but A. thaliana is self-fertilizing). Second, there is always the possibility that results are driven by sampling effects, including demographic history, although using two data sets and separate partitions of those data sets somewhat discounts that notion here. Third, and perhaps most importantly, it is difficult to disentangle genetic from epigenetic effects. Overall, however, this work suggests that gbM has a measurable effect on fitness. The estimated selection coefficients were small (4Nes = 1.4) but nonetheless similar to the magnitude of selection acting on codon usage that has been measured in A. thaliana, A. lyrata, and Capsella rubella.

17One interesting feature of codon bias, a phenomenon widely accepted as a genomic feature that is under weak selection, is that it varies among species, with selection detectable in species with large historical population sizes but not detectable in small Ne species. An overarching feature of gbM is that much of the experimental and comparative work on gbM has focused on A. thaliana. It is worth noting that this species may be atypical in at least three respects. First, two independent studies relying on different data sets and approaches have inferred that A. thaliana has lost gbM three times faster than gaining it relative to closely related outcrossing species. Second, the recent shift of A. thaliana to an inbreeding mating system reduced its effective population size, which is likely to have weakened the efficacy of selection on gbM in that species. Finally, methylation mutants usually have little phenotypic effect in A. thaliana, whereas they are often lethal in taxa with higher TE load, such as maize. Together these observations suggest that A. thaliana may not be the best model for measuring the evolutionary effects of methylation, and yet there is still some evidence that selection acts on gbM in that species, raising the possibility that the effects of gbM may be more pronounced in other species. For this reason, we advocate that similar analyses are extended to other taxa when large methylation data sets become available.

18Given that there are some indications that gbM may be under weak selection, one naturally wonders what its function might be. One consistent hypothesis has been that gbM affects gene expression levels. This hypothesis first came from the observation that genic methylation levels across genes within A. thaliana are associated with expression levels: methylated genes tend to be intermediately to highly expressed, with lower expression variance among tissues. These patterns have been interpreted in two ways: either gbM might affect expression patterns (fig. 2a) or, conversely, active transcription might drive gbM. Many highly expressed genes do not have gbM in A. thaliana, an observation that discounts the second hypothesis or at least suggests that the relationship is not completely straightforward. Moreover, it is now known that CMT3 does not depend on gene expression to methylate genes but instead on inaccessible chromatin marks and heterochromatin histone variants, although it remains possible that the initial recruitment of CMT3 requires or depends on gene expression.

19One difficulty in assessing the effect of gbM on gene expression comes from the possible confusion between genetic and epigenetic effects. Indeed, variation in gene expression can be caused by numerous factors—for example, by nearby single nucleotide polymorphisms (SNPs) in regulatory sequences, by a nearby TE insertion (genetic cis effects), by a change in a transcription factor (trans effects), or by a change in the gene DNA methylation level (epigenetic cis effects). Genome-wide association studies (GWASs) and epigenome-wide association studies have shown that DNA methylation variants associated with expression variation are often in linkage disequilibrium with nearby SNPs, making it difficult to disentangle the respective contribution of SNPs and methylation variation on gene expression. However, found a significant association between cis-methylation and gene expression in hundreds of genes across 135 A. thaliana accessions. Interestingly, gbM was positively correlated with gene expression, and the effect remained significant after controlling for SNPs. Overall, the number and magnitude of the affected loci by DNA methylation were smaller than the effect of SNPs, and hence, the authors concluded that DNA methylation has limited effects on expression variation (table 1; ).

20The association between gbM and expression was further tested experimentally in epigenetic recombinant inbred lines (epiRILs) obtained through the cross of a met1 mutant and wild-type (WT) A. thaliana, followed by eight generations of inbreeding. The resulting epiRILs have a mosaic methylome, with some regions derived from the met1 mutant that originally lacked gbM and other regions containing CG methylation derived from the WT parent. inferred differentially expressed genes between the met1-derived regions of epiRILs and their WT homologs. They found only 6 out of 3,471 genes that were gbM in WT plants and became differentially expressed when located in met1-derived regions in epiRILs. On the other hand, they found significantly more genes (46 out of 3,124, P-value = 2.55 × 10−9) that were unmethylated in WT plants and became differentially expressed when located in met1-derived regions in epiRILs. Taken together, these results suggest that gbM loss has little, if any, effect on gene expression (table 1). However, if gbM has a small effect on expression, it is likely to have been missed by differential expression analyses that typically detect individual genes with twofold or more expression differences. More subtle effects may be statistically detectable only by approaches that summarize trends across multiple genes. More importantly, the met1 mutant might be a poor system to study the association between gene methylation and expression level because both methylated and unmethylated genes were upregulated in met1 mutants using microarray data.

21Another approach to test for associations between gbM and expression has been to use comparative and evolutionary rather than experimental approaches. For example, several studies have compared expression between gbM-deprived E. salsugineum with its close relative A. thaliana, but the results have been controversial. In the first study, estimated the expression of unmethylated E. salsugineum genes that are orthologous to gbM genes in A. thaliana and found no difference in expression between species (table 1). reanalyzed the data from and used genes that were unmethylated in both A. thaliana and E. salsugineum as a negative control to measure the average difference in expression between the two species. When taking into account this species effect in a linear model, gbM loss in E. salsugineum was associated with a small but significant decrease in expression (table 1). In the third study using the same data, disagreed on the use of unmethylated genes as a negative control because they have been shown to have more variable expression levels over evolutionary time and again found no effect of gbM loss on gene expression.

22Another effort compared gbM and unmethylated genes between A. thaliana and A. lyrata. Methylated genes were expressed at significantly higher levels, on average, and with less variation between species than non-CG-methylated genes. The authors identified genes that changed methylation status between A. thaliana and Arabidopsis lyrata to examine whether the shift in methylation correlated with gene expression. They found that genes that had gained gbM in one of the two species also tended to shift toward higher expression levels, but these results were not statistically significant (table 1). However, genes that differed in gbM status between A. thaliana and A. lyrata exhibited significantly higher variance in expression between species than genes that were gbM in both species, consistent with previous studies suggesting gbM modulates expression variability. Another comparative study compared the methylomes of eight grass species and found that genes that were gbM in all eight species tended to have higher and less variable expression compared to genes that varied in their methylation state across species. Although the effect was very small, the results suggested a positive effect of gbM on expression level and expression stabilization (fig. 2c). It is worth emphasizing, however, that this approach, like most comparative approaches, cannot determine causality. More recently, we used the 876 A. thaliana methylomes to study the association between gbM and gene expression within a species by comparing the methylation state of alleles both to their expression level and to the variability in expression across the larger data subset of 876 A. thaliana methylomes. Across genes with polymorphic methylation states, the expression of gbM alleles was consistently and significantly higher than unmethylated alleles (table 1). Taken across the entire genome, gbM alleles also had a significantly less variable expression level compared to unmethylated alleles of the same gene. Although consistent across the thousands of genes in the data set, the effect was quite small: on average, a methylated allele had ∼1 more RNA-seq read than an unmethylated allele. A weakness of this work is that it did not disentangle potential genetic effects from epigenetic effects; however, the gbM effect did remain consistent when models included a proxy for genetic variation, by including the number of CG dinucleotides in statistical analyses. Consistent with our A. thaliana results, work on the outcrossing crucifer Capsella grandiflora has revealed that the presence of gbM is a major predictor of cis-regulatory constraint. GbM lowers the probability of allele-specific expression via cis-regulation, again suggesting a stabilizing effect of gbM on expression level.

23As we have just reviewed, several studies have established an association between gbM and stable expression level (fig. 2b and c), which complement the proposal that gbM has a homeostatic effect on expression. This phenomenon has been further investigated by, who studied gene expression levels in A. thaliana roots via single-cell RNA-seq. They found no significant correlation between gbM and gene expression noise (as measured by the variation in the expression level among single cells). However, gbM was significantly positively correlated with gene expression consistency, which they measured as the number of single-cell RNA-seq replicates in which the gene was expressed (fig. 2b). This effect remained after correcting other genomic features such as gene expression, gene length, gene conservation, and gene duplication status. Therefore, found that gbM genes are more consistently expressed than unmethylated genes across cells of a tissue, which can be interpreted as implying that gbM is involved in the maintenance of a consistent gene expression (fig. 2b). If this is true, the mechanism by which this happens remains unknown. One hypothesis comes from the anticorrelation observed between genome-wide distributions of the histone variant H2A.Z and DNA methylation in A. thaliana. H2A.Z is typically associated with transiently expressed response genes, such as immune-responsive or environmental stimulus-responsive genes. In met1 mutants, the loss of DNA methylation was accompanied by a gain in H2A.Z deposition. However, in an h2a.z mutant, DNA methylation patterns were only minimally affected, suggesting that DNA methylation prevents H2A.Z incorporation but not the converse. Based on these observations, it has been proposed that gbM serves to stabilize transcription by preventing deposition of the histone variant H2A.Z.

24Finally, have used quantitative trait loci (QTL) mapping to identify ∼1,000 genes for which the proportion of methylated CG sites significantly correlates with expression level across a sample of over 900 natural A. thaliana accessions. The variance in expression explained by CG methylation is modest for most genes, but for some genes, it reaches levels comparable to the effect of SNPs on expression. gbM is mostly positively correlated with expression; in contrast, TE-like methylation (i.e., in all three contexts) is, as expected, negatively correlated with expression. In a clever extension to control the effect of linked SNPs, identified SNPs with significant effects on expression using GWAS. They then repeated the analysis linking CG methylation and expression within nested sets of accessions that carry the same GWAS allele. For the vast majority of genes, this approach confirmed the significant positive correlation between gbM and expression level, either because there was no GWAS SNP or because at least one nested sample had a significant correlation. A second control analyzing haplogroups, which corrects for cis genetic variation, led separately to the same conclusion. They then studied gene expression in the met1 A. thaliana mutant without gbM, and they found as expected a reduced expression level in genes with these three characteristics: (1) genes with a significant positive correlation between CG methylation and expression level, (2) genes that were methylated in the WT Col-0 accession (which is the accession used for the met1 mutant), and (3) genes for which the effect of gbM was not confounded by linked genetic variants. The authors concluded that gbM is positively correlated with gene expression in hundreds of genes, independently of local genetic variants.

25Although did attempt to account for trans effects as well as cis effects, a shortcoming of most of the comparative studies referred to in this section is that they do not account for possible trans genetic effects on gene expression, which may result in overestimated cis epigenetic effects. Altogether, however, evolutionary and comparative studies tend to find small but detectable relationships between gbM and either gene expression levels or gene expression variation. These results contrast with many direct experimental measurements of gbM based on A. thaliana mutants. If the comparative conclusions are correct, they are important because they suggest that gbM has a phenotype that may be the target of natural selection.

26To date, studies have been inconsistent as to whether gbM associates with gene expression (table 1). When it does associate with expression, it is also difficult to disentangle cause from effect. If, however, we assume there is a real relationship between gbM and gene expression, there remains an open question: what is the mechanism(s) by which gbM affects expression? One hypothesis is that gbM improves transcription by regulating alternative promoters within gene bodies, thereby potentially preventing aberrant internal and/or antisense transcription (fig. 2d;; ). This hypothesis stems from the observation that CG methylation is typically depleted within active promoter regions and also that genes with CG-methylated promoters are silenced. Aberrant transcription, whether in the sense or antisense orientation, is expected to be deleterious because it is energetically costly and leads to the accumulation of both unnecessary transcripts and truncated proteins that can be toxic for the cell. Aberrant antisense transcription is expected to disturb gene expression because RNA polymerases coming from both directions may collide. Moreover, the RNA-directed DNA methylation pathway can be activated by the pairing of sense and antisense transcripts into double-stranded RNA, which may further prevent gene expression. Hence, if gbM prevents aberrant reverse transcription, it could explain the aforementioned association between gbM and gene expression.

27However, the results of tests for this effect have been inconsistent (table 1). Some of these tests have taken place in mammalian systems because they too exhibit CG methylation within genes, even though they mostly do not methylate in the CHG and CHH contexts. For example, studied mouse embryonic stem cells in DNA methyltransferase 3b (Dnmt3b) mutants that lack gbM and compared them to WT. To quantify internal transcription, the authors used a ratio of the number of RNA-seq reads that map to the third exon divided by the number of reads mapping to the first exon (hereafter exon3/exon1). This ratio is expected to increase when there is cryptic intragenic initiation of transcription. found that the loss of gbM was accompanied by higher exon3/exon1 ratios than WT mice, suggesting an increase of spurious internal transcription between exons 1 and 3 when gbM is lost (table 1). The results were confirmed by a sequencing technique that allows characterization of the exact position of the 5′ end of mRNAs by targeting the mRNA cap (DECAP-seq). However, performed similar analyses in Dnmt triple mutants on highly expressed genes, and they were unable to corroborate the findings of (table 1). Results from suggest that the role of gbM in suppressing spurious transcription initiation may be specific to the lack of DNMT3B, but only while other DNMTs are still present.

28Similar work has sought evidence for an effect of gbM on aberrant transcription in plants. For example, compared met1-derived regions of A. thaliana epiRILs with orthologous WT regions. They quantified antisense transcription and found that gbM loss did not lead to an increase in differentially expressed antisense transcripts (table 1). However, detected that the expression of antisense transcripts was activated in 938 genes in h1, met1 double mutants compared to WT. The number of upregulated antisense transcripts was comparatively low in single mutants when compared with WT (145 and 34 for met1 and h1, respectively), suggesting redundancy in H1 and MET1 repression of antisense transcription. This finding demonstrates that, at least for some genes, gbM may repress antisense transcription in A. thaliana jointly with histone H1 (table 1). This study also again exemplifies redundancy among DNA methylation, histone variants, and histone marks. These different epigenetic marks are interdependent and play overlapping roles in the cell, complicating the characterization and inference of potential gbM effects.

29More recently, used RNA long reads sequenced by Oxford Nanopore Technology Direct Sequencing (ONT DRS) to characterize transcription start sites (TSSs) in A. thaliana. They found that the met1–3 mutant, which lacks CG methylation, has significantly more unique TSSs compared to WT, and these unique TSSs occurred in regions where mutant methylation was lower than WT. These results suggest that gbM can prevent the initiation of aberrant transcription. The transcription termination site (TTS) was also affected by DNA methylation. Indeed, the met1–3 mutant had a higher number of unique TTS than WT, indicating that CG methylation also inhibits aberrant transcription termination. Altogether, this work suggests that gbM could ensure proper transcription of genes from start to end (fig. 2d).

30Here, we revisited this issue by analyzing Isoseq (PacBio RNA long read) data in maize and A. thaliana (Supplementary Materials and Methods, Supplementary Material online). We included maize because this is (to our knowledge) the first such attempt to examine this question in a plant other than A. thaliana. We focus on Isoseq data because it can represent full-length mRNA, thanks to the selection of mRNAs that contain a 3′ poly-A tail and, in some cases, a 5′ cap (when sequencing is done with a cap-trap step). The A. thaliana Isoseq data set we analyzed has a 5′ cap so that most Isoseq reads likely represent full-length mRNAs (Supplementary Materials and Methods). Some of the A. thaliana data set was publicly available, and we also generated new Isoseq data for this study; in both cases, the data were generated from Col-0 inflorescences. In contrast, the maize Isoseq data, which was generated on pooled RNA extracted from six tissues at different developmental stages of the B73 inbred line, were not generated with a cap-trap step. With both the maize and A. thaliana data sets, we considered aberrant internal transcription to be reflected in Isoseq reads that begin after the start of exon 1 (fig. 2d). For each gene, the proportion of full-length Isoseq reads with a “conventional” TSS (i.e., that begin prior to the start of exon 1) was computed and compared between gbM and unmethylated (UM) genes. gbM and UM genes were categorized from publicly available methylation data (see Supplemental Materials and Methods, Supplementary Material online).

31In maize, we found that gbM genes had a significantly higher proportion of conventional TSSs (average 0.81) compared to UM genes (average 0.78, Wilcoxon test P-value = 5.63 × 10−11fig. 3a); superficially, this observation complies with the prediction that gbM genes have less aberrant transcription. However, gbM is known to be associated with more highly expressed and longer genes, and these covariates must be taken into account. Genes with higher expression had a significantly higher proportion of conventional TSSs (generalized linear model contrast estimate = 0.042, Z-ratio = 69.15, P-value <2 × 10−16, Supplementary table S1, Supplementary Material online), perhaps reflecting higher selective pressures to remove aberrant transcription for highly expressed genes. Longer genes had significantly fewer conventional TSSs (generalized linear model contrast estimate = −0.178, Z-ratio = −130.7, P-value <2 × 10−16Supplementary table S1, Supplementary Material online), which could be attributable to the higher probability of a long gene harboring an aberrant internal promoter or an experimental artifact (i.e., longer genes may have a higher chance of having their mRNA not fully reversed transcribed during sequencing, leading to 5′ truncated transcripts and wrongly inferred aberrant TSS). Notably, however, only 321 of 2059 detected nonconventional TSSs occurred in introns. After taking gene length and gene expression into account in a generalized linear model (Supplementary Materials and Methods, eq. 6, Supplementary Material online), UM genes had a significantly higher proportion of conventional TSSs compared to gbM genes (generalized linear model contrast estimate = 0.214, Z-ratio = 33.4, P-value <2 × 10−16, Supplementary table S1, Supplementary Material online). This result is not consistent with the expectation that gbM prevents aberrant TSS (table 1).

32The A. thaliana results complement the maize results in some ways but not others. Since the data were generated with a 5′ cap, the overall proportion of conventional TSSs in the A. thaliana Isoseq data set was higher than in maize (fig. 3a). gbM genes had a lower proportion of conventional TSSs (mean 0.90) compared to UM genes (0.96). However, after taking gene length and gene expression into account in a generalized linear model, gbM genes did have significantly more conventional TSS compared to UM genes (P-value <2 × 10−16, Supplementary table S2, Supplementary Material online). We conclude that the Isoseq data do reflect some advantage of gbM in terms of avoiding internal transcription start in WT A. thaliana but not in maize (table 1).

33We explored these ideas further by turning to a different approach that relies on RNA-seq reads in gbM mutants. Because 5′ cap Isoseq data are not available for methylation mutants, we inferred internal transcription starts using the RNA-seq coverage ratio of exon3/exon1, following the approach of (see Supplementary Materials and Methods, Supplementary Material online). If gbM prevents aberrant internal transcription start, this ratio should increase in gbM mutants relative to WT. We therefore measured exon3/exon1 RNA-seq coverage in WT and gbM mutants. We performed this comparison for two data sets based on two different gbM mutants. In data set 1, RNA-seq data were generated on 13-day-old seedlings for three replicates of WT controls that were compared to three replicates of met1–3 mutants. In data set 2, RNA-seq data were generated on leaf tissue for three replicates of WT controls (these differed from the controls within data set 1) that were then compared to five replicates of met1,sdg7–8 triple mutants. In WT plants, gbM genes had a lower exon3/exon1 coverage ratio compared to UM genes (fig. 3b), suggesting superficially that gbM could prevent internal transcription start. This was observed consistently for WT plants from data set 1 (mean exon3/exon1 coverage 0.662 in gbM genes, 1.919 in UM genes, one-sided Wilcoxon test P-value <2.2 × 10−16) and also from data set 2 (mean exon3/exon1 coverage 1.013 in gbM genes, 3.684 in UM genes, one-sided Wilcoxon test P-value = 1.3 × 10−11). However, this same difference between gbM genes (as defined in WT plants) and UM genes was also apparent in mutant plants that lacked gbM (fig. 3b). That is, gbM genes had a lower exon3/exon1 ratio compared to UM genes in met1–3 mutants (0.512 vs. 1.876, one-sided Wilcoxon test P-value <2.2 × 10−16) and in met1,sdg7–8 triple mutants (1.232 vs. 1.566, one-sided Wilcoxon text P-value <2.2 × 10−16). These observations suggest that the fact that gbM genes have lower exon3/exon1 coverage in RNA-seq data is not due to their methylation state alone because the same pattern is observed in mutants without gbM. While one must always be careful that mutants can be complex and may reflect other (unknown) effects, the data again provide little support for the notion that gbM alone prevents aberrant internal transcription initiation in A. thaliana genes (table 1). Interestingly, showed that gbM and H1 play a redundant role in inhibiting aberrant reverse transcription, and hypothesized that another epigenetic mark, CHH islands, may have some redundant function with gbM. These redundancies could explain why we did not detect any change in internal transcription start in A. thaliana gbM mutants because these redundant epigenetic marks may have been functioning in gbM mutants, thus complicating inferences about gbM effects.

34Another aspect of aberrant transcription, aside from internal transcription initiation, is reverse transcription. We again used the Isoseq data from A. thaliana to estimate the proportion of full-length antisense reads, which correspond to reverse transcription events (see Supplementary Materials and Methods, Supplementary Material online). gbM genes had a significantly lower mean proportion of antisense reads (average 0.0046) compared to UM genes (average 0.016, Wilcoxon test P-value = 3.16 × 10−12). This result held after accounting for gene length and gene expression in a generalized linear model (Supplementary table S3, Supplementary Material online). However, in maize Isoseq data, the gbM genes had significantly more antisense transcription compared to UM genes (Supplementary table S4, Supplementary Material online). We conclude that there is evidence that gbM prevents antisense transcription in A. thaliana based on Isoseq data from WT plants (table 1), but we have uncovered no such evidence in maize. Altogether, however, we believe there are enough compelling observations—both from animals and from A. thaliana plants—to suggest that further dissection of this potential function may be worthwhile.

35Another hypothesis is that gbM improves splicing fidelity and prevents aberrant intron retention (fig. 2e), but this raises the question of how splicing fidelity may drive a relationship between gbM and gene expression. One possibility is that poor splicing in the absence of gbM leads to the retention of introns (fig. 2e) that contain premature stop codons. Aberrant transcripts containing premature stop codons are typically sent to the nonsense-mediated mRNA decay for destruction, which might in turn lower gene expression. This suggests a potential relationship among gbM, splicing fidelity, and gene expression.

36The effect of gbM on splicing fidelity has been tested across various taxa, and the results have been—like studies of aberrant transcription—somewhat inconsistent. In honeybee and mouse embryonic stem cells, for example, it is clear that alteration of DNA methylation impacts alternative splicing. In honeybees, DNA methylation is predominantly on gene bodies in the CG context. A knockdown of the expression of dnmt3, which is required for de novo DNA methylation, decreased global genomic methylation level and caused widespread changes in alternative splicing in fat tissue. In mouse embryonic stem cells, constructed an experimental system in which differential DNA methylation could be limited to a single gene while all other cellular factors remained identical. Using this system, they demonstrated a direct causal relationship between DNA methylation and the recruitment of splicing factors. Patterns of methylation near splice sites have also been studied in maize, where CHG methylation of the splicing acceptor site is associated with a lower efficiency of splicing and CHH methylation does not correlate with splicing efficacy. Surprisingly, however, the effect of CG methylation was not tested explicitly. followed the maize work by measuring splicing fidelity in A. thaliana. They found that gbM was negatively correlated with the amount of RNA-seq reads that map to introns, suggesting that gbM genes tend to retain fewer introns in their mRNA compared to UM genes. Similarly, recently used ONT DRS to characterize splicing in A. thaliana. They found that retained introns had significantly lower CG methylation levels around their splicing sites (both donor and acceptor sites) compared to spliced introns in WT and some CHG and CHH methylation mutants. This suggests that gbM facilitates splicing. However, found no evidence for this splicing effect when they compared met1 epiRILs to WT plants. In fact, they found that WT gbM genes retained significantly fewer intron reads than UM genes after they lost gbM in the met1 background. This work suggests that this intron effect is a property of the genes rather than gbM per se.

37Given contradictory results in the literature, we further tested the hypothesis that gbM prevents aberrant intron retention using A. thaliana Isoseq data (Supplementary Material and Methods, Supplementary Material online). The proportion of full-length Isoseq reads that retained at least one intron was higher in gbM genes (mean 0.149) compared to UM genes (mean 0.106). This result remained significant after taking gene length and gene expression into account in a generalized linear model (table 1, Supplementary table S5, Supplementary Material online). We also measured intron RNA-seq coverage in WT and mutant A. thaliana plants that lack gbM (Supplementary Material and Methods, Supplementary Material online). Similar to, we found that gbM genes had a lower intron read coverage compared to UM genes in WT plants (fig. 3c) both in data set 1 (mean gbM genes intron coverage 61.51 RPKM, vs. 133.42 for UM genes, one-sided Wilcoxon test P-value < 2.2 × 10−16) and in data set 2 (mean gbM genes intron coverage 53.13 RPKM, vs. 208.68 in UM genes, one-sided Wilcoxon test P-value < 2.2 × 10−16). However, the difference between gbM genes and UM genes was also found in mutant plants that lack gbM (fig. 3c). Indeed, gbM genes had a lower intron read coverage compared to UM genes in met1–3 mutants (59.04 vs. 131.65 RPKM, one-sided Wilcoxon test P-value < 2.2 × 10−16) and in met1, sdg7–8 triple mutants (51.27 vs. 132.97 RPKM, one-sided Wilcoxon test P-value < 2.2 × 10−16). This again suggests that the fact that gbM genes have lower intron coverage in RNA-seq data is not due to their methylation state alone (table 1). Another possibility is again that some other epigenetic mark plays a redundant role with gbM in preventing aberrant intron retention. In summary, there is not yet a clear consensus, or even a clear trend, as to whether gbM plays a role in splicing fidelity. However, we note again that most of the work in plants has focused on A. thaliana, and, as previously discussed, it may be helpful to extend these analyses to other species. The question of the potential role of gbM in splicing fidelity, like its role in aberrant transcription, remains unresolved and will benefit from the broader and more comprehensive investigation.

38Another hypothesized function of gbM is that methylation protects against the insertion of some TEs (fig. 2f). This hypothesis primarily stems from two studies in maize that focused on Robertson’s Mutator (Mu) transposons, which typically insert within or near genes and can be highly deleterious by disrupting gene function. In the first study, found that Mu transposons insert preferentially within unmethylated regions of the B73 genome. However, the methylation context could not be determined, which motivated to repeat the analyses with context-specific DNA methylation data. They found that Mu transposon insertion sites within genes were strongly depleted in CG-methylated regions (table 1), but these regions were not depleted in CHG nor CHH methylation relative to average gene methylation. This raises the possibility that gbM is beneficial because it deters transposon insertions. This hypothesis is also difficult to disentangle from covariates, particularly the observation that gbM genes tend to be under stronger selective constraint than unmethylated genes, as measured by nonsynonymous divergence, although there is conflicting evidence that gbM genes do or do not evolve more rapidly in the Poaceae. Nonetheless, it is possible that the apparent difference in TE insertion reflects different strengths of selection on gbM versus unmethylated genes, rather than a direct effect of gbM on TE insertion rate.

39It would be insightful to repeat analyses of the effect of DNA methylation on TE insertion in other plant species and with different TE types to test for the potential broader relevance of this idea. If this phenomenon occurs across diverse species and TE types, it could partially explain the link between gbM and expression stabilization within and between species. Indeed, a genic TE insertion might prevent proper gene expression because silencing of the TE by methylation in the three contexts might spread to the gene and affect expression. Also, if the TE is inserted within an exon, it might lead to the appearance of premature stop codons and the destruction of mRNA by nonsense-mediated mRNA decay. This would lead to more expression variation among accessions of a species and also among species (fig. 2c).

40The molecular mechanisms leading to gbM establishment and maintenance in plants have been remarkably well elucidated (fig. 1). However, the function and potential importance of gbM remain debated, as illustrated in this review by the numerous studies that are inconsistent or, in some cases, contradictory (table 1). The main sources for this persistent uncertainty come first from the difficulty in disentangling epigenetic from genetic effects and second from the complex system of redundancies and overlapping functions among gbM, histone variants, and other epigenetic marks. These dependencies and redundancies undoubtedly complicate the interpretation of experimental mutants, as illustrated by contrasting results based on met1 plants with or without h1 mutation (table 1;; ).

41Despite these difficulties, there has been substantive progress toward understanding the effect, dynamics, and potential adaptive impact of gbM. In the last few years, several studies have concluded that gbM is associated with both the level and the variance of gene expression (table {PI}1). An interesting corollary of these observations is that many experimental efforts to measure this association based on A. thaliana mutants have yielded negative results (table 1). While these experiments may reflect reality, our view is that experimental approaches have nonetheless suffered from a few common shortcomings. First, we suspect (but certainly do not know) that the reliance on A. thaliana is limiting; it is illogical to expect a strong experimental effect in a system that is studied in part because methylation mutants are viable and thus may not have a strong effect relative to other plant systems. Second, as noted above, it is not clear that all mutants are equivalent because some mutants may be unsuitable for detecting specific effects due to dependencies and redundancies. Finally, it can be exceedingly difficult to identify subtle effects using short-term experimental approaches. Unfortunately, however, the inability to detect an effect is often incorrectly interpreted as the absence of an effect.

42In contrast, evolutionary and comparative approaches based on genetic diversity or species comparison have often, but not always, found an association between gbM and gene expression, particularly expression homeostasis (table 1). These analyses also suffer from a number of potential drawbacks, including reliance on simplified models, discrete definitions of which genes are (or are not) gbM, and an inability to disentangle causation from correlation. One potential reason for detecting an effect using these approaches is that even very subtle effects can accrue over time and thus be detected by evolutionary contrasts. It is difficult to establish whether these associations are causal due to all the complex reasons cited above—that is, functional redundancies among epigenetic marks and difficulties in discriminating genetic from epigenetic effects. Nonetheless, the apparent association between gbM and expression is important, because it provides a potential phenotype on which selection can act. Although more investigation is needed to test whether gbM is shaped by selection, both phylogenetic and population genetic studies suggest that selection acts to maintain gbM status in some genes. Moreover, population variation in gbM has been shown to associate with fitness under water stress and selection for flowering time. These fitness effects appear to stem from a correlation between gbM and gene expression, as demonstrated by the fact that experimental modification of candidate gene expression affects the trait under study. Altogether, these results suggest that gbM may affect fitness and phenotype through an effect on gene expression, thereby potentially affecting species adaptation independently of genetic variation. We emphasize, however, that these effects are subtle, at best, and so, there is still much to learn, even about the simple question as to whether and when gbM associates with gene expression (table 1).

43If there is selection on gbM itself—or on another epigenetic feature that correlates with gbM—then this fact may be the basis for an “important conceptual change in evolutionary biology” because selection on epigenetic modifications has the potential to affect the timeframe of mutation and thus, potentially, of adaptative change. As an example, figure 4 illustrates the tempo of epigenetic versus genetic change. As we have noted, epigenetic change can be incredibly rapid. For example, CHH methylation is reset every generation, making it unlikely to be under direct selection because it is not transmitted transgenerationally. CHH methylation (and to some extent CHG methylation) is perhaps better described as tracking genetic (i.e., TE and CMT3) activity. In contrast, CG methylation is heritable and mutates approximately three and six orders of magnitude faster than gene duplication and nucleotides, respectively (fig. 4). It is worth noting, however, that the rate of change of an entire gene or allele from gbM to UM is substantially slower because it probably requires numerous changes of individual sites. The rate of change for an entire allele, based on our previous SFS analysis in A. thaliana populations, is ∼3 × 10−7 per gene per generation, an estimate based on models that require assumptions about the effective population size but, if accurate, is still faster than the rate of nucleotide change. In theory, then, methylation variation may provide a rapid source of phenotypic novelty that could be subjected to natural selection on rapid—and perhaps even ecological—time scales.