Epigenetic Variation Induced by Gamma Rays, DNA Methyltransferase Inhibitors, and Their Combination in Rice

DNA methylation plays important roles in the regulation of gene expression and maintenance of genome stability in many organisms, including plants. In this study, we treated rice with gamma rays (GRs) and DNA methyltransferase inhibitors (DNMTis) to induce variations in DNA methylation and evaluated epigenetic diversity using methylation-sensitive amplified polymorphism (MSAP) and transposon methylation display (TMD) marker systems. Comparative and integrated analyses of the data revealed that both GRs and DNMTis alone have epimutagenic effects and that combined treatment enhanced these effects. Calculation of methylation rates based on band scoring suggested that both GRs and DNMTis induce epigenetic diversity by demethylation in a dose-dependent manner, and combined treatment can induce variations more synergistically. The difference in the changes in full and hemi-methylation rates between MSAP and TMD is presumed to be caused by the different genomic contexts of the loci amplified in the two marker systems. Principal coordinate, phylogenic, and population structure analyses commonly yielded two clusters of individuals divided by DNMTi treatment. The clustering pattern was more apparent in TMD, indicating that DNMTis have a stronger effect on hypermethylated repetitive regions. These findings provide a foundation for understanding epigenetic variations induced by GRs and DNMTis and for epigenetic mutation breeding.


Introduction
DNA methylation contributes to the epigenetic regulation of gene expression, genomic stability, and transposon silencing and plays important roles in multiple biological processes, including plant and mammalian development [1,2]. DNA methylation has been a focus of recent research that has sought to identify beneficial new traits, and it has been established that hypermethylation is generally correlated with the downregulation of gene expression, whereas hypomethylation is correlated with upregulation [3,4]. DNA methylation maintains the methylation state after a methyl group has been added via the action of DNA methyltransferases and may also occur passively when methylation does not occur in the newly synthesized DNA strand during DNA replication [5,6]. Additionally, 5-methylcytosine can be actively removed by the base excision repair pathway [2]. The chemical inhibition of DNA methyltransferases using DNA methyltransferase inhibitors (DNMTis), such as 5-azacytidine (AZA) and zebularine (ZEB), as demethylating agents is an approach used to study the effects of DNA methylation loss in plants [7,8]. Although AZA and ZEB have similar inhibitory effects on DNA methylation, ZEB has a longer half-life and better stability than AZA [8,9]. Many studies have investigated the effects of DNMTis on gene expression and the silencing, growth, development, and induction of DNA damage in plants [10][11][12][13].
Transposable elements (TEs) are ubiquitous in the genomes of eukaryotes, and their relocalization within genomes serves to generate genomic plasticity by inducing various types of genic and intergenic mutations that can increase allelic diversity [14][15][16]. After being discovered first by McClintock in the 1940s, TEs were initially considered junk or selfish DNA [17,18]. However, they are now recognized to play important roles as a driving force of evolution in eukaryotic taxa [19]. TEs can be classified into two types based on their mechanism of transposition: class I retrotransposable elements and class II DNA transposons. Class I retrotransposons make new copies in different chromosomal locations via a "copy-and-paste" mechanism of transposition, whereas class II DNA transposons move via a "cut-and-paste" mechanism [20]. Transposition activity, including the movement of DNA transposons or the insertion of new copies of retrotransposons in a genome, can affect genomic stability [21]. Accordingly, TEs tend to be hypermethylated relative to other genomic sequences in eukaryotic genomes [2,22,23]. In Arabidopsis, methylation occurs at high levels in heterochromatin and euchromatin, including in some transposons, and this high level of DNA methylation seems to inhibit transposon mobility [23]. A high methylation rate in TEs and adjacent genomic regions was also reported in rice [24].
Molecular marker systems are useful tools that can be used in various fields, such as taxonomy, plant breeding, and genetic engineering, to assess genetic diversity, identify and fingerprint genotypes, determine genetic distances between populations, strains, and breeding materials, and be the basis of marker-assisted selection [25][26][27]. Amplified fragment length polymorphism (AFLP), one of the most popular PCR-based marker systems, uses two sticky-end restriction enzymes and ligate adapters to detect mutations in restriction sites dispersed throughout genomes. Methylation-sensitive amplified polymorphism (MSAP) is a modified AFLP method used to detect DNA methylation polymorphisms using methylation-sensitive HpaII and its methylation-insensitive isoschizomer MspI, and both enzymes recognize 5 -CCGG-3 sequences [28,29]. Transposon methylation display (TMD) is a modified MSAP method that uses an additional primer specific to transposons to identify methylation patterns at CCGG sites [4].
Rice is an important cereal crop for human consumption, with over 100,000 cultivars worldwide [30]. It is a diploid species with a small genome and has become the most well-studied crop with respect to functional genomics, as it is used as a monocot model plant [31,32]. Various methods of mutation breeding have been developed to identify and improve rice characteristics, and extensive studies have been conducted on mutagenesis using physical and chemical agents [31]. In terms of physical mutagens, gamma rays (GRs) have been used in the development of more than 90% of rice varieties, in addition to X-rays, neutrons, and protons [33][34][35]. Regarding chemical mutagens, studies using ethyl methane sulfonate [36], sodium azide [37], and N-methyl-N-nitrosourea [38] have been reported. Furthermore, the efficiency, effectiveness, and induced chlorophyll mutation rate of combined treatment with GRs and ethyl methane sulfonate have been reported [36,[39][40][41]. To date, however, there have been no studies that have examined the effects of combined treatment with GRs and DNMTis on epigenetic mutations.
Here, we investigated methylation variation in rice plants generated from seeds subjected to treatment with GRs, DNMTis, and their combination by MSAP and TMD analyses. The enzyme combinations HpaII/MseI and MspI/MseI were used for MSAP analysis. A p-SINE1 retrotransposon-specific primer with HpaII/MspI adapters was adopted for TMD analysis.

Evaluation of Polymorphisms Detected by MSAP and TMD Marker Systems
Six and four primer combinations were used for MSAP and TMD analyses, respectively (see details in Materials and Methods). Each primer combination generated 10 to 40 discernable loci ranging in size from 150 to 500 bp. Ultimately, MSAP and TMD yielded 102 and 60 loci, 66 and 41 of which were polymorphic, respectively ( Table 1). All indices related to polymorphism level, including the percentage of polymorphic loci, Nei's gene diversity (H), and Shannon's information index (I), were higher in TMD than in MSAP (Table 1).

Methylation Changes Induced by Treatment with GRs and DNMTis
Methylation levels and patterns in rice treated with GRs and DMNTis were estimated based on the band scoring of individual plants. In the MSAP analysis, the total methylation rates in the control (58.69%) and GR-irradiated groups (100 Gy [G100], 58 54.35%) and combined (G100 + A80, 53.94%; G150 + A80, 54.04%; G100 + Z80, 54.66%; and G150 + Z80, 52.63%) treatment groups showed lower rates ( Table 2). The combined treatment group showed the most reduced methylation rate (average 4.87%) relative to the control ( Table 2). Both the full and hemi-methylation rates were lower in all treatment groups than in the control, although this was more apparent in the DMNTi and combined treatment groups ( Table 2).

Analysis of Molecular Variance
We performed analysis of molecular variance (AMOVA) in order to determine whether the epigenetic variations observed among rice groups treated with either GRs or DNMTis were differentiated from each other. The MSAP and TMD data showed that 74% and 67% of the variances could be attributed to differences within groups, respectively ( Table 4). The variance among groups shown by TMD (33%) was higher than that shown by MSAP (26%), indicating that higher epiallele diversity was induced in the loci amplified in TMD. The AMOVA result for the total data (MSAP and TMD data combined) showed that the variance within groups and among groups was 71% and 29%, respectively. The correlation coefficient between the epigenetic distance of the two marker systems was calculated using the Mantel test, which revealed that there was no correlation (r 2 = 0.1612).  4 Estimates of variance components. 5 Percentage of total variance contributed by each component. 6 Probability of a random value greater than or equal to the observed data value based on permutation test (n = 999).

Population Structure and Phylogenetic Analysis
Principal coordinate analyses (PCoA) were performed to investigate the epigenetic relationships among the different treatment groups (Figure 1). In the PCoA plot for MSAP, the control and GR-irradiated groups were distributed closely, and some individuals were located closer to those in the DNMTi and combined treatment groups, which showed a scattered distribution (Figure 1a). The PCoA plot for TMD showed two distinct clusters: a cluster of the control and GR-irradiated groups and another cluster of the DNMTi and combined treatment groups (Figure 1b). In other words, these two clusters were divided whether the individuals were treated with DNMTis or not. This trend was similar in the PCoA plot for merged MSAP and TMD data, although a few individuals were separated from the clusters (Figure 1c). Population structure and phylogenetic relationships were analyzed based on the MSAP and TMD data. The calculation of delta K suggested that the optimal number of clusters (K) was three for all MSAP, TMD, and merged data ( Figure S1). The three phylogenetic trees showed similar clustering trends: the control and GR-irradiated groups formed one cluster (depicted by red branches in Figure 2) and the DNMTi and combined treatment groups formed two intermingled clusters (depicted by yellow and blue branches in Figure 2). There were some exceptions: seven plants in the combined treatment group (A80 + G100) were found in the MSAP tree. However, these exceptions were not notable because they clustered outside of the main cluster with different epigenetic properties in the MSAP tree, and they were distributed far from the clusters of the control and GR-irradiated groups in both the TMD and merged data trees. The three clusters (K = 3) were therefore clearly divided according to data source.

Discussion
Methylation changes in plant DNA are molecular mechanisms that confer adaptations to various environmental factors, thus contributing to biodiversity [43,44]. In the present study, we found that epigenetic variations were induced by GRs, DNMTis, and their combination in rice plants based on MSAP and TMD analyses. Indices related to epigenetic diversity were increased in all treatment groups compared to the control ( Table 1), indicating that GRs as well as DNMTis can affect DNA methylation status. The two types of DNMTi used here, AZA and ZEB, induce epigenetic variation by reducing the methylation rate [13]. Ionizing irradiation is also known to have a demethylating effect, which was deduced from the activation of TEs by hypomethylation in irradiated plants [33,45]. Low-energy ion beam irradiation induced both methylation and demethylation in Arabidopsis [46]. Together, the results of the present study indicate that both GRs and DNMTis play important roles in increasing epigenetic diversity and that combined treatment can induce more methylation variations.
The total methylation rates, calculated based on banding patterns, were reduced in all treatment groups (Tables 2 and 3). As the treatment dose became higher, the MSAP data showed decreases in both the full and hemi-methylation rates (Table 2), whereas the TMD data showed a decrease in the full methylation rate and an increase in the hemi-methylation rate. This may be explained by the different genomic context of the loci amplified in each marker system. DNA methylation in plant genomes is distributed mainly in repetitive regions to stabilize genomes by inhibiting the activation of TEs [2,22,47]. The short interspersed retroelement p-SINE1 used in this study for TMD, which was first identified in rice, is dispersed across the genome with a high copy number (approximately 6500 copies per haploid genome) [48]. This retroelement is a repetitive sequence of 70-500 bp with an internal promoter for RNA polymerase III, which is involved in transcription, although it lacks open reading frames [48,49]. A high level methylation of p-SINE1 has been reported not only within the retroelement itself but also in adjacent regions in the rice genome [24]. Therefore, the increase in hemi-methylation shown by TMD may be a result of fully methylated repeat-rich regions being partially demethylated by treatment with GRs and DNMTis, which converted fully methylated loci into hemi-methylated loci.
AMOVA revealed that a large proportion of the epigenetic diversity observed in the present study was due to within-group variations (Table 4). This suggests that the methylation variations induced by GRs and DNMTis occur almost randomly. However, a higher percentage of variation among groups was observed in the TMD data than in the MSAP data (Table 4), indicating that the amplified loci in TMD, which may be highly methylated, were more dynamically affected by treatment with different doses of GRs and DNMTis. The differential influence on methylation variation between GRs and DNMTis could be inferred from the PCoA and population structure analysis (Figures 1 and 2). The two clusters of rice plants clearly divided by DNMTi treatment in the PCoA plot (Figure 1b) reflect the fact that demethylation induced by DNMTis can produce more epigenetic diversity in hypermethylated repetitive regions. This tendency is similarly observed in the TMD-based phylogenetic tree and its population structure (Figure 2b). The intermixed cluster containing the DNMTi and combined treatment groups branches from the cluster of the control and GR-irradiated groups, which reflects increased methylation variation due to demethylation by DNMTis.
In this study, we revealed methylation variation induced by GRs and DNMTis in rice, which can be synergistically increased when the two mutagens are administered in combination. Epigenetic diversity is potentially of value for crop breeding, including for widening the range of phenotypic variation, increasing stress tolerance with less gene erosion, and balancing major agronomic traits in a species [50]. Although the question of whether epigenetic changes can be transmitted to progeny has remained controversial, a great deal of experimental evidence of transgenerational epigenetic inheritance has been reported, particularly with regard to DNA methylation [51]. The present study was conducted only on M 1 plants, so further investigation into subsequent generations is needed. The findings of this study will lay the foundation for understanding the variations in DNA methylation artificially induced by GRs, DNMTis, and their combination.

Plant Materials and Genomic DNA Extraction
Solutions of two types of DNMTi, AZA (Sigma-Aldrich, St Louis, MO, USA) and ZEB (TCI, Tokyo, Japan), were prepared in 0.5% dimethyl sulfoxide (DMSO; Sigma-Aldrich, St Louis, MO, USA). Rice (Oryza sativa L. ssp. japonica cv. Ilpum) seeds were soaked in 80 µM DNMTi solution for 24 h at room temperature. Seeds in the control group (no treatment with either GRs or DNMTis) and the GR irradiation only group were pre-soaked in 0.5% DMSO solution for 24 h at room temperature. Seeds were irradiated using a 60 Co gamma irradiator at the Korea Atomic Energy Research Institute (Jeongeup, Republic of Korea) for 8 h at a dose of 100, 150, or 250 Gy (dose rates of 12.5, 18.75, and 31.25 Gy h −1 , respectively). The DNMTi solution treatment dose and time (80 µM for 24 h, respectively) were determined in our preliminary study, in which rice treated at this dose for this length of time had a survival rate of approximately 70%. In the combined treatment, irradiation doses of 100 and 150 Gy were applied because plants from seeds treated with DNMTi and ≥200 Gy GR irradiation did not survive after emergence.
Young leaves were harvested from plants during the three-leaf stage, and genomic DNA was extracted using a DNeasy Plant Mini Kit (Qiagen, Germantown, MD, USA). Quantification of the extracted DNA was performed using a NanoDrop spectrophotometer, and DNA integrity was assessed based on 1.0% agarose gel electrophoresis. Aliquots of the quantified total genomic DNA were subsequently diluted to a final concentration of 500 ng µL −1 .

Molecular Marker Analysis
The MSAP analysis was conducted according to Reyna-López et al. [52]. Genomic DNA (500 ng) was digested for 16 h at 37 • C with 5 U of each restriction enzyme combination (MseI/MspI and MseI/HpaII) in 10 µL of buffer solution (Invitrogen, Carlsbad, CA, USA). The digested DNA was ligated in a final reaction volume of 25 µL with 15 pmol enzyme adapters (listed in Table S1) using T4 DNA ligase (Invitrogen, Carlsbad, CA, USA) for 16 h at 20 • C and diluted 1:4 with Tris-EDTA (TE) buffer. The ligated DNA (3 µL) was amplified with 10 pmol pre-amplification primers in a final reaction volume of 25 µL. The PCR conditions for pre-amplification were as follows: 95 • C for 2 min, followed by 20 cycles of 95 • C for 30 s, 55 • C for 30 s, and 72 • C for 1 min, and a final extension at 72 • C for 5 min. The pre-amplified product was diluted 1:100 with TE buffer, and 3 µL were used for selective amplification with the MSAP primer combinations listed in Table 4 under the following PCR conditions: 95 • C for 2 min, followed by 25 cycles of 95 • C for 30 s, 55 • C for 30 s, and 72 • C for 1 min, and a final extension at 72 • C for 5 min. The PCR products were diluted 1:4 with TE buffer and segregated using a Fragment Analyzer 5300 (Agilent, Santa Clara, CA, USA).
TMD analyses with p-SINE1 were carried out according to the protocol described by Takata et al. [24], with minor modifications. DNA digestion, pre-and selective amplification, and segregation of PCR products were carried out under the same conditions as MSAP with the TMD primer combinations listed in Table S2.

Data Analysis
DNA bands between 150 and 600 bp were scored as either present (1) or absent (0), and indistinct bands were excluded to avoid potential genotyping errors. POPGENE version 3.2 [53] was used to calculate the polymorphic index, Nei's genetic diversity (H), Shannon's information index, and the effective number of alleles. A genetic diversity matrix was produced using GenAlEx version 6.5 [54] for Mantel's tests, PCoA, and AMOVA. Genetic distances and similarity coefficients were calculated using PowerMarker version 3.25 [55]. The neighbor-joining clustering method was used to construct an unrooted phylogenetic tree of rice individuals using MEGA X [56]. The population structure was analyzed based on Bayesian inference using STRUCTURE version 2.3.4 [57]. Calculations of K values from 1 to 10 were repeated 10 times to obtain log Pr (X|K) values. STRUCTURE HARVESTER web version 0.6.94 was used to estimate the optimal value of K [58].