Skip to main content
Advertisement
Browse Subject Areas
?

Click through the PLOS taxonomy to find articles in your field.

For more information about PLOS Subject Areas, click here.

  • Loading metrics

Differential DNA methylation and transcription profiles in date palm roots exposed to salinity

  • Ibtisam Al-Harrasi,

    Roles Conceptualization, Data curation, Formal analysis, Investigation, Methodology, Validation, Writing – original draft

    Affiliation Department of Biology, College of Science, Sultan Qaboos University, Muscat, Oman

  • Rashid Al-Yahyai,

    Roles Conceptualization, Resources, Supervision, Writing – review & editing

    Affiliation Department of Crop Sciences, College of Agricultural and Marine Sciences, Sultan Qaboos University, Muscat, Oman

  • Mahmoud W. Yaish

    Roles Conceptualization, Funding acquisition, Investigation, Methodology, Project administration, Resources, Supervision, Writing – review & editing

    myaish@squ.edu.om

    Affiliation Department of Biology, College of Science, Sultan Qaboos University, Muscat, Oman

Abstract

As a salt-adaptive plant, the date palm (Phoenix dactylifera L.) requires a suitable mechanism to adapt to the stress of saline soils. There is growing evidence that DNA methylation plays an important role in regulating gene expression in response to abiotic stresses, including salinity. Thus, the present study sought to examine the differential methylation status that occurs in the date palm genome when plants are exposed to salinity, and to identify salinity responsive genes that are regulated by DNA methylation. To achieve these, whole-genome bisulfite sequencing (WGBS) was employed and mRNA was sequenced from salinity-treated and untreated roots. The WGBS analysis included 324,987,795 and 317,056,091 total reads of the control and the salinity-treated samples, respectively. The analysis covered about 81% of the total genomic DNA with about 40% of mapping efficiency of the sequenced reads and an average read depth of 17-fold coverage per DNA strand, and with a bisulfite conversion rate of around 99%. The level of methylation within the differentially methylated regions (DMRs) was significantly (p < 0.05, FDR ≤ 0.05) increased in response to salinity specifically at the mCHG and mCHH sequence contexts. Consistently, the mass spectrometry and the enzyme-linked immunosorbent assay (ELISA) showed that there was a significant (p < 0.05) increase in the global DNA methylation in response to salinity. mRNA sequencing revealed the presence of 6,405 differentially regulated genes with a significant value (p < 0.001, FDR ≤ 0.05) in response to salinity. Integration of high-resolution methylome and transcriptome analyses revealed a negative correlation between mCG methylation located within the promoters and the gene expression, while a positive correlation was noticed between mCHG/mCHH methylation rations and gene expression specifically when plants grew under control conditions. Therefore, the methylome and transcriptome relationships vary based on the methylated sequence context, the methylated region within the gene, the protein-coding ability of the gene, and the salinity treatment. These results provide insights into interplay among DNA methylation and gene expression, and highlight the effect of salinity on the nature of this relationship, which may involve other genetic and epigenetic players under salt stress conditions. The results obtained from this project provide the first draft map of the differential methylome and transcriptome of date palm when exposed to an abiotic stress.

Introduction

Soil salinity is one of the most serious environmental obstacles affecting the growth and the productivity of different plant species [1]. By the year 2050, it is expected that more than 50% of arable land will become unsuitable for agricultural uses due to salinization [2]. Therefore, production of salt-tolerant plants is an urgent need necessary to minimize this problem. However, an understanding of the mechanisms underlying salinity adaption in plants is a crucial step toward approaching this target [3].

Plants have evolved complex mechanisms to cope with salinity stress as a multigenic trait, including selective ion uptake and exclusion, synthesis of compatible solutes, photosynthesis and energy metabolism control, toxic ion compartmentalization, anti-oxidative enzyme accumulation, cell structure modification and hormonal regulation [4, 5]. Salinity adaptation involves various physiological processes; however, the efficiency of these processes in salinity tolerance is not only based on the presence of the required genes, but also on the presence of an appropriate pattern of gene expression, a procedure which is controlled in the cell by various epigenetic factors [6].

Epigenetic changes involve three main mechanisms: DNA methylation, histone modifications, and RNA-mediated gene silencing [7]. Cytosine methylation results from the addition of a methyl group (CH3) to the carbon 5-position of cytosine bases (5mdC) [810]. This process plays a crucial role in genomic imprinting, as well as in gene expression regulation, transposon silencing and transgenerational epigenetic memory [1114]. In mammalians and plants, DNA methylation occurs in mCG dinucleotide sequence context, as well as in mCHG and mCHH sequence contexts, where H stands for all nucleotides except guanine [1517]. The cytosine methylation within symmetrical mCG is maintained by the MET1 methyltransferase (DNMT1) [18] and methylation of symmetrical mCHG is maintained by CHROMOMETHYLASE 3 (CMT3) methyltransferase, while asymmetrical mCHH methylation is maintained by persistent de novo methylation [19].

There is an important correlation between DNA methylation status and gene expression when plants are exposed to salinity [2023]. However, the impact of DNA methylation on gene expression varies based on the location of methylated cytosines within the gene body and the methylated sequence contexts [24]. For example, while mCG methylation within the promoter inhibits gene expression [25], the same methylation pattern within the exons enhances gene expression in some species [8, 18, 25, 26].

The date palm (Phoenix dactylifera L.) is a relatively salt-tolerant plant. However, few reports have studied the genes and the mechanisms that code salinity tolerance [2730]. Despite their potential importance in gene regulation under salinity, epigenetic factors including DNA methylation have not yet been investigated in date palms.

DNA methylation changes in response to salinity have previously been studied in different plant species using various molecular techniques, including the methyl-sensitive amplification polymorphism (MSAP) [3133]. The outcomes of these studies provided evidence about the importance of this epigenetic process in controlling a plant’s response to salinity. However, the detection of differential methylation patterns in response to salinity based on the construction of a single-base resolution methylome has not yet been reported in date palms.

In this study, we hypothesized that cytosine methylation is differentially occurring in the genome when a date palm is prolonged exposed and then adapted to salinity stress. These changes can distinctively affect the expression level of some genes, including those that are associated with salinity tolerance mechanisms. Therefore, the main aims of this study were to explore the patterns of the DNA methylation across the genome of plants grown under control conditions, and to investigate the alterations in the DNA methylation that would ensure the response to salinity, based on a single nucleotide resolution. This study also was aimed to correlate the cytosine methylation with the transcriptome abundance at global levels. This includes those transcripts that have potential association with salinity tolerance mechanisms.

The results showed that the date palm methylome is heavily modulated by salinity treatment and that these modulations were statistically correlated with gene expression patterns.

Materials and methods

Plant materials, growth conditions and photosynthetic measurements

Date palm (Khalas, a moderately salinity-tolerant cultivar) seeds, derived from a single source of pollen grains, were thoroughly washed with tap water then surface-sterilized in 75% ethanol for 5 minutes. The ethanol traces were washed out by rinsing the seeds four times with sterile distilled water. Thereafter, the seeds were soaked overnight in sterile water at 37°C. The next day, the seeds were mixed with sterilized moist vermiculite and incubated at 37°C in the dark for one week. The germinated seeds were transferred to 2L pots containing peat moss and sand mixture (1:1, v/v). Seedlings were grown under controlled conditions of a 16/8-hour light/dark cycle, 350μE m−2 s−1 light intensity, 35/30°C day/night temperature and 60% humidity. The control treatment group was irrigated to field capacity with distilled water twice a week. The salt stress was gradually applied at an increment of 50 mM every 10 days, to a final concentration of 300 mM NaCl. The 300 mM NaCl treatment lasted for 20 days. The electrical conductivity (EC) of the soil was used as an indication of salinity. The level of soil salinity was measured three days after each salt application using the Em50 Digital Data Logger (Decagon Devices, WA, USA). FluorPen FP100 (Photon System Instruments, Czech Republic) was used to measure chlorophyll fluorescence-related parameters. Net photosynthetic rate (A), internal CO2 concentration (ci), stomatal conductance (gs) and transpiration (E) were measured using the LCpro-SD Advanced Photosynthesis Measurement System (ADC Bioscientific Ltd., UK).

Root analysis

Roots of the NaCl-treated and untreated seedlings were scanned, and some phenotypic parameters such as root length, diameter, total surface area and the number of root tips were analyzed using the WinRHIZO (RH-R XLR STD) software (version 5.0, Regent Instruments, Inc., Quebec, QC, Canada).

DNA and RNA extraction

Root samples of eight seedlings grown either under control or salinity treatments were harvested, pooled and used for total DNA and RNA extraction. Root tissues were crushed in a pre-chilled mortar to fine powder in liquid nitrogen and the genomic DNA was extracted from the grinded samples as previously described [34]. Total RNA was extracted from the same plant tissues, following a protocol described by Xiao et al. [35].

Whole genome bisulfite sequencing (WGBS)

DNA and RNA samples were analyzed using Zymo Research CORP next generation sequencing facilities (Irvine, CA, USA), following the previously described strategy [36]. Genomic DNA (500 ng) was fragmented with dsDNA Shearase Plus enzyme (2 units). The resulted fragments were further processed in order to construct Methyl-MaxiSeq libraries. The DNA fragments were blunted to sticky ends using the End-It DNA End-Repair Kit. Subsequently, terminal A(s) tails were added to the 3’ end of the sticky-ended fragments while methylated adaptors, which contain 5’-methylcytosine instead of cytosine, were ligated into the 5’ end. The constructed libraries were underwent a bisulfite treatment, using the EZ DNA Methylation Kit (Zymo Research, Irvine, CA, USA), and the DNA fragments were amplified using specific amplification conditions as previously described [37]. Finally, the bisulfite-treated libraries were sequenced using the Illumina Hiseq platform to generate 50 bp paired-end reads on the Illumina HiSeq 1500.

Processing and alignment of the bisulfite-sequenced libraries

Bisulfite-sequenced reads were identified using standard Illumina base-calling software. After removal of the duplicate reads, unique best alignments were retained using Bismark version 0.14.3 software [38] and the genome coverage was calculate using Bedtools software version 2.26.0 [39]. A non-directional parameter was applied while running the alignment. In order to calculate the conversion efficiency of cytosines after treatment of DNA with sodium bisulfite, filtered reads were first aligned to the plastid genome of Phoenix dactylifera L. (Genbank accession number GU811709) [40]. Subsequently, the reads were aligned to the available reference date palm genome in the database (Genbank accession number PRJNA249070) [41].

Identification of differentially methylated sites (DMSs)

The methylation ratio of each site was calculated by dividing number of methylated cytosine(s) by the total number of cytosine(s) found in the reads and covered at the specific site. Only sites covered by at least five sequence reads were considered for further analysis. The differential methylation ratio due to salinity treatment was measured by subtracting the methylation ratio obtained from the control sample, from the methylation ratio obtained from the salt-treated sample at each 5mdC site. Significant differences in methylation ratios between treated and control samples were statistically determined by Fisher’s exact test (p < 0.05, FDR ≤ 0.05). The p value was then adjusted to q value using Benjamini-Hochberg correction for multiple testing. Accordingly, differentially methylated sites (DMSs) of the three sequence contexts (mCG, mCHG and mCHH) were classified based on the value of the difference in the methylation ratios into the following categories: DMSs were classified as “strongly hypermethylated” if the ratio difference was between +33% and 100%, and “hypermethylated” if the ratio difference was less than +33% but above 0%; however, they were considered “strongly hypomethylated” if the ratio difference was between ‒33% and ‒100%, and “hypomethylated” if the ratio difference was less than ‒33% but below 0%.

Identification of differentially methylated regions (DMRs)

The chromosomes that passed the p -value < 0.05 filter, from the methylation site files for each methylation context was selected for further analysis. Subsequently, 50-bp window sliding was used to obtain the DMRs in only those chromosomes that passed filter. For the methylation ratio for a region, the methylation ratio of each site (5X minimum read coverage) within the region of each sample was first calculated, then the values was averaged for the region and subsequently the difference of the average methylation ratio between the control and salinity-treated samples was considered as the methylation difference for the region. For p -value calculation, a fisher test was used. Unmethylated count was obtained from the methylation ratio and total count for a region for both samples. Based on this information, a 2 by 2 table was built to calculate a two-tailed p -value for each region.

The overall change (increase) in methylation percentage upon exposure to salinity was calculated based on the following formula: ((Hypermethylated DMRs ‒ Hypomethylated DMRs)/ Total Significant DMRs))*100.

RNA Sequencing and transcriptome analysis

After rRNA depletion, Illumina HiSeq 50 bp single-end reads from the two-pooled RNA samples were constructed and the resulted data were processed following the manufacturer’s instructions (Illumina). The re-sequencing strategy was used to determine the expression level where the short reads were aligned to the P. dactylifera L. reference genome (genbank accession number PRJNA249070) using TopHat software (version 2.0.13). Cufflinks software (version 2.2.1) was used for transcript assembly and differential expression, and cummeRbund software (version 2.0.0) was used for visualization of differential analysis. The expression level was normalized by considering the fragments per kilobase of exon model per million mapped fragments (FPKM) value. Differentially expressed genes (DEGs) with a significant value (≥ 2-fold changes, p-value ≤ 0.001, FDR ≤ 0.05) were selected for further analysis. To account for the computational necessity of applying a logarithmic transformation to the expression data, any no expression value (0.0) was replaced with a pseudocount (0.158882), which was calculated based on a tenth of the value of the lowest value among the expression data of the two samples. Subsequently, differential gene expression value was obtained using the formula: log2 (FPKM of treated/FPKM of control).

Identification of differentially methylated genes (DMGs)

In this analysis, the region comparison was carried out against the known genes in the annotation file, resulted from the RNA sequencing. The methylation ratios across the entire gene for each sample were calculated based on the methylation ratios of sites that are covered within the coordinates for each gene and have a minimum coverage of 5X. For each sample, the average of the methylation ratios for each site was assigned as the average methylation ratio of the gene for that sample. The methylation differences for the gene are the average of this methylation ratio. Since the total count for the entire gene for each sample and the average methylation ratio were obtained, a fisher test was run and a p-value for each gene was generated. The fisher calculation is similar to the one for the DMRs.

To investigate the correlation between the expression level and the methylation status, the expression data were first grouped based on the mapping location within the different gene features (promoter, exon, and intron regions). The DMGs with no methylation value were replaced by a pseudocount (2.75 X 10‒7), which was calculated as the tenth of the lowest value among the methylation ratios obtained from the samples. Pearson’s correlation was calculated between the log2 of mRNA expression (FPKM) values and the log2 of the methylation values using IBM SPSS Statistics 21. The correlation test applied for both the control and the treated samples.

Functional Annotation and analysis of DMGs

The DMGs containing mCG, mCHG or mCHH sequence contexts were functionally annotated using the Blast2GO PRO software package [42]. Only DMGs of differential gene expression values were included in this analysis. The annotation of the transcript was based on similarity with protein-coding mRNA sequences available in the protein databases. The protein sequences were classified into the categories of biological processes (BP), cellular components (CC) and molecular functions (MF) using gene ontology (GO) annotation. Differential functional enrichment analysis between the upregulated test group and the downregulated reference group were identified using Fisher’s exact test, based on p ≤ 0.05. The coded enzymes were mapped to the metabolic pathways using the Kyoto Encyclopedia of Genes and Genomes (KEGG) tools [43], which are implemented within the Blast2GO PRO software.

Quantification of global DNA methylation (5mdC) using mass spectrometer

Global genomic quantification of 5-hydroxymethyl-2′-deoxycytidine (5hmdC) and 5mdC was carried out using the mass spectrometer available at the Zymo Research, Irvine, CA, USA. DNA samples (400 ng) were digested to single nucleoside using DNA Degradase Plus (Zymo Research, USA). 5hmdC and 5mdC were then quantified by an SRM-based mass spectrometry assay. The assay quantified the concentration of 5mdC and 5hmdC as a percentage of 2′-deoxyguanosine (dG), the internal standard [44]. The calibration ranges (using 40 pmol of dG) were between 0 and 2.5% and between 0 and 25% for 5hmdC and 5mdC, respectively. Each DNA sample was run in triplicate.

Determination the global 5mdC level using enzyme-linked immunosorbent assay (ELISA) based colorimetric assay

The global (5mdC) level was quantified using DNA samples extracted from roots grown under control and salt stress conditions. Four biological replicates per treatment were used in this analysis. The quantification was carried out using the MethylFlash Methylated DNA 5mdC Quantification Kit (Epigentek Group Inc., Farmingdale, NY, USA), according to the manufacturer’s instructions.

The sequence DNA and RNA read data generated in this report were deposited in the read data archive of the NCBI under the sequence read project (SRP) numbers SRP108055 and SRP108259, respectively.

Results

Effect of salinity on the growth, photosynthetic gas exchange and root development

Photosynthesis performance was used in this project as an indication of plant health status upon salinity treatment. At the end of the salinity treatment, the average electrical conductivity of the soil of the control and NaCl treatment pots was 1.12 ± 0.06 (mean ± SD) and 18.03 ± 0.351 (mean ± SD) dSm-1, respectively. A morphological difference between two treatment-groups’ seedlings was noticed with respect to leaf length (Fig 1A). Salt treatment caused significant (p < 0.05) reduction of the photosynthetic rate (A) (Fig 1B), stomatal conductance (gs) (Fig 1C) and transpiration rate (E) (Fig 1D) in date palm seedlings. However, there was an insignificant increase in the internal carbon dioxide concentration (ci) under salt stress (Fig 1E).

thumbnail
Fig 1. The effect of salinity treatment on photosynthesis performance (A), photosynthesis rate (A) (B), stomatal conductance (gs) (C), transpiration rate (E) (D), and internal CO2 concentration (ci).

The values are means (n = 4) and error bars represent SE (p < 0.05). Significant differences are marked by an asterisk.

https://doi.org/10.1371/journal.pone.0191492.g001

Salinity showed a negative effect on the root system (Fig 2). This effect was demonstrated by a significant (p < 0.05) reduction in the root length (Fig 2A), surface area and number of root tips (Fig 2C and 2D). The salinity however, had a positive impact on the root diameter (Fig 2B).

thumbnail
Fig 2. The effect of salinity on root growth, which is measured as root length (A), diameter (B), surface area (C), and number of tips (D).

The values are means (n = 4) and error bars represent SE (p < 0.05). Significant differences are marked by an asterisk.

https://doi.org/10.1371/journal.pone.0191492.g002

WGBS revealed alterations in the methylation levels upon exposure to salinity stress

After bisulfite conversion, the whole genomes were sequenced using next generation sequencing tools. This type of DNA sequencing has the ability to detect the sites of the methylated cytosine as well as their relative abundance. The bisulfite conversion rates of the genome isolated from untreated (control) and treated plants were 99% and 98%, respectively. The plastid genome of date palm was used to determine the conversion efficiency because it is widely accepted that plastid genomes included methylation free cytosines [45], and this approach is commonly used in plant methylome studies [46].

The WGBS yielded 324,987,795 and 317,056,091 reads from the sequencing of the genomes isolated from the control and the salinity-treated plants, respectively (Table 1). Among these, 19,209,086 and 19,185,920 were aligned to unique positions of mCGs on the genomes extracted from the control and salinity-treated plants, respectively. The WGBS covered around 40% of the total genomic cytosines with an average read depth of 17-fold coverage per each DNA strand (Table 1). The percentage of genome coverage was 80.6% and 80.8% for the control and the salinity-treated samples, respectively. The analysis also revealed the presence of 25,155,548 mCHG and 120,635,111 mCHH sites in the genome isolated from the control plants. However, the same analysis showed the presence of 25,134,635 mCHG and 120,591,563 mCHH sites in the genome isolated from salinity-treated plants (Table 1).

thumbnail
Table 1. WGBS analysis of DNA samples extracted from roots grown under control (0 mM NaCl) and the saline (300 mM NaCl) conditions.

https://doi.org/10.1371/journal.pone.0191492.t001

WGBS revealed an increase in the significantly methylated 5mdC in response to salinity treatment

In order to determine the effect of salinity treatment on the DNA methylation status, a pairwise comparison of the DMRs was carried out. The analysis showed that 2,130,787 DMRs were significantly (p < 0.05) altered in response to salinity stress (Table 2). Statistical analysis indicated that the DMRs resulted from CHH methylation were more frequent in date palm genome (46.93%) comparing with CG (27.84%) and CHG methylation (25.23%). Using the methylation status in the DNA extracted from the control plants as a reference, the WGBS analysis showed that there were around 55.84% and 44.16% significantly (p < 0.05) hypermethylated and hypomethylated regions, respectively, in DNA extracted from plants exposed to salinity stress (Table 2). Collectively, statistical analysis of the DMRs resulted from three sequence contexts methylation showed a significant (p < 0.05) increase of about 24.5% in the net methylation level of the DMRs in the genome upon exposure to salinity. It was noticed that the percentages of hypermethylated DMRs were higher than hypomethylated counterparts in the case of mCHG and mCHH, but this situation was the opposite in the case of the mCG context, where the percentages of the hypermethylated DMRs were less than their hypomethylated counterparts. This implies that DNA methylation at the mCG region decreased when plants were exposed to salinity, while it increased at the mCHG and mCHH regions due to the same stress (Table 2) (Fig 3).

thumbnail
Table 2. DMRs obtained from the WGBS analysis.

The DMRs with a significant abundance value (p < 0.05) in response to salinity are showed.

https://doi.org/10.1371/journal.pone.0191492.t002

thumbnail
Fig 3. Box plots display the overall methylation level of DMRs (p < 0.05) under control and salinity stress conditions (p < 0.05 considered to indicate statistical significance.

*p < 0.05, **p < 0.01, ***p < 0.001, independent sample t test (two-tailed) was conducted to compare the means of two samples). Bars represent the median.

https://doi.org/10.1371/journal.pone.0191492.g003

Cluster analysis of the top 100 methylated sites

In order to clearly determine the effect of salinity on the status of DNA methylation, the DMSs that showed significant (p < 0.05) changes in methylation levels and had the highest methylation difference between the DNA samples (top 100 DMSs) in the mCG, mCHG and mCHH sequence contexts, regardless of their location within the genome, were selected for hierarchically clustering analysis. The resulted heat maps showed some variations in the level of DNA methylation between the two genomes at the mCG DMSs. However, these variations were enormous in the mCHG and mCHH sequence contexts (S1 Fig). This profound variation indicates the important effect of salinity treatment in modulating DNA methylation.

The level mCG methylation in the promoter and gene body regions

Fortunately, functionally annotated genes and their promoters were significantly covered in this analysis. DNA methylation percentage within the promoter, exon and intron regions of the annotated genes was estimated in the DNA extracted from treated and untreated roots. The DMRs that showed a significant methylation variation (p < 0.01) in the studied regions in response to salinity stress was used in the analysis. The results showed that the exonic DMRs of mCG methylation had higher methylation level than the promoter’ DMRs, regardless of the salt stress treatment (Fig 4A and 4B). Generally, there was an increase in the methylation level of DMRs located within promoters after exposure to salinity stress (Fig 4A). This trend was also noticed with exonic DMRs except for these related to mCHH sequence context (p < 0.01) (Fig 4B). However, the intronic DMRs of mCG and mCHH methylation showed a reduced (p < 0.01) DNA methylation level due to the salt treatment. Contradictory, the methylation level of the intronic DMRs of mCHG was increased (p < 0.01) due to the same treatment (Fig 4C).

thumbnail
Fig 4. Box plots display the overall methylation level of DMRs (p < 0.01,) among promoter (A), exon (B) and intron (C) regions, under control and salinity stress conditions (p < 0.01 considered to indicate statistical significance.

*p < 0.05, **p < 0.01, ***p < 0.001, independent sample t test (two-tailed) was conducted to compare the means of two samples). Bars represent the median.

https://doi.org/10.1371/journal.pone.0191492.g004

DNA methylation landscape across the genes and the flanking regions

In order to determine the landscape of the DMSs within the annotated genes and their flanking regions, and also to study the dynamics of DNA methylation in response to salinity stress in these regions, an overview of the global DNA methylation pattern was constructed based on the distribution of the mCG, mCHG and mCHH DMSs across the defined gene features that are available in the genome database. This includes the protein-coding region, the transcription start site (TSS), the transcription end site (TES) and their 3Kb flanking regions (Fig 5).

thumbnail
Fig 5. DNA methylation overview across the gene and the flanking regions based on the distribution of mCG, mCHG and mCHH differentially methylated sites (p < 0.001, FDR ≤ 0.05) under control (red lines) and salinity stress (black lines) conditions.

These sites were plotted against gene structural features, including the promoter, gene body and flanking regions, which include transcription start site (TSS), transcription end site (TES), and 3Kb up- and down-stream regions.

https://doi.org/10.1371/journal.pone.0191492.g005

The results revealed that the DNA methylation patterns among the mCG level were relatively low (~ 35%) in the ‒3Kb upstream the TSS regions; however, this level was dramatically dropped to the lowest level in the TSS regions, then gradually increased to nearly a constant level (~ 50%) in the coding region. This methylation level decreased again toward the TES regions before it increased again in the downstream +3Kb regions (Fig 4A). The methylation level in the mCHG sequence context was the highest (~ 20%) among the other regions in the promoter sequences and decreased close to the TSS, but this level slightly increased in the coding regions and decreased again near the TES region. Subsequently, this DNA methylation level gradually increased in the +3Kb flanking regions (Fig 5B). A similar trend of DNA methylation was also observed at the mCHH DMSs (Fig 5C); however, the DNA methylation level was too low (~ 3% at most) in the different gene regions. Collectively, salinity treatment did not show a clear effect on DNA methylation levels at the 5mdCs located within the gene regions.

The methylation landscape was also described for the annotated promoters and their flanking regions (Fig 6).

thumbnail
Fig 6. DNA methylation overview across the promoter and the flanking regions based on distribution of mCG, mCHG and mCHH differentially methylated sites (p < 0.001, FDR ≤ 0.05) under control (red lines) and salinity stress (black lines) conditions.

These sites were plotted against promoter sequences and flanking regions, which include 1Kb up- and down-stream promoter region.

https://doi.org/10.1371/journal.pone.0191492.g006

mCG showed the highest level of DNA methylation in the start regions of the promoters as compared to the other sequence contexts; however, it was observed that the DNA methylation of mCG and mCHG was dramatically decreased toward the end of the promoters (Fig 6). It was not possible to perceive obvious differences in DNA methylation profiles between the DNA extracted from the treated and untreated plants using this method.

Mass spectrometry and ELISA revealed an increase in the global genomic 5mdcs in response to salinity

The mass spectrometry analysis revealed undetectable amounts of 5-hydroxymethyl-2’-deoxycytidine (5hmdC). However, quantification of global 5mdC revealed the presence of 5.23 ± 0.08% (mean ± SD) and 7.88 ± 0.90% (mean ± SD) of 5mdC of the total cytosines of DNA samples extracted from the control and NaCl-treated plants, respectively. Statistical analysis revealed that these amounts account for a significant (p = 0.007) increase of 50.7% in global DNA methylation occurring in the genomes of the salinity-treated plants. Similarly, ELISA showed a relatively significant (p = 0.05) increase of 37.2% in the global 5-mCs of the DNA samples extracted from roots of plants grown in saline conditions.

DMSs profiles across the genes and their flanking regions showed distinct influence on gene expression

Sequencing of the mRNA purified from roots revealed a total number of 59,330,361 and 62,543,394 reads sequenced from the salinity-treated and untreated tissues, respectively (Table 3). The genome coverage was about 88% and 89% for the control and the salt treated samples, respectively. These reads were combined into 370,538 contigs, with a sequence length varied from 21 to 168,337 bp. The analysis indicated the presence of around 6,405 DEGs in response to salinity (S1 Table), such that 64% and 36% of these genes were significantly (p < 0.001, FDR ≤ 0.05) upregulated and downregulated, respectively.

thumbnail
Table 3. Re-sequencing results of total mRNA isolated from NaCl-treated and control root tissues.

https://doi.org/10.1371/journal.pone.0191492.t003

In order to investigate the relationship between DNA methylation and gene expression, the DEGs under salinity treatment were classified into five different groups based on their expression value (first high, second high, medium, lower, and no expression). Subsequently, DNA methylation profiles at the mCG, mCHG and mCHH sites within the functional (protein coding) genes (Fig 7) and pseudogenes (Fig 8) were correlated with the corresponding gene expression value of each group. Under the control conditions, the results showed that there was a negative correlation between the level of mCG DNA methylation within the ‒3Kb region and gene expression levels. This relationship was also found in the +3kb region (Fig 7A).

thumbnail
Fig 7. DNA methylation (mCG, mCHG and mCHH) profile among DEGs (p < 0.001, FDR ≤ 0.05), grouped based on their expression level (first high (blue line), second high (red line), medium (green line), low (beige line) and no expression (black line)).

The Y-axis represents the DNA methylation ratio of each group, while the X-axis represents studied gene regions of DMGs including protein-coding region, the transcription start site (TSS), the transcription end site (TES) and 3Kb up and down-steam the TSS and the TES, respectively.

https://doi.org/10.1371/journal.pone.0191492.g007

thumbnail
Fig 8. DNA methylation (mCG, mCHG and mCHH) profile among pseudogenes, grouped based on their expression level (first high (blue line), second high (red line), medium (green line), low (beige line) and no expression (black line)).

The Y-axis represents the DNA methylation ratio of each group, while the X-axis represents studied gene regions of DMGs including protein-coding region, the transcription start site (TSS), the transcription end site (TES) and 3Kb up and down-steam the TSS and the TES, respectively.

https://doi.org/10.1371/journal.pone.0191492.g008

On the other hand, under salinity conditions, high gene expression was associated with high mCG methylation levels, especially when it occurred in the gene body regions (Fig 7B). However, when plants grew under control conditions, DNA methylation at the mCHG DMSs located within all regions was associated with low gene expression (Fig 7C). When plants were exposed to salinity, a slight level of methylation within the coding regions was associated with gene expression (Fig 7D). The results also showed that there was unclear effect of the DNA methylation level at the mCHH sites on gene expression (Fig 7D). Irrespective of the environmental conditions, a clear negative relationship was observed between 5mdC methylation levels around TSS and gene repression levels (Fig 7). However, mCG methylation around TES was associated with gene activation only when plants were exposed to salinity (Fig 7B).

In order to compare the effect of DNA methylation on the functional gene and its effect on the pseudogenes in terms of transcriptome abundance, DNA methylation levels of the three sequence contexts were plotted against the different groups of pseudogenes classified based on their expression level, as described above (Fig 8). In general, the DNA methylation levels within the pseudogene regions were relatively higher than the levels found within the functional genes (Figs 7 and 8). For example, the highest level of mCG methylation observed in the pseudogenes was around 80%, while the highest level in the coding genes was around 50%. The results also showed that DNA methylation in the mCG and mCHG sequence contexts was positively correlated to gene expression over most of the studied regions when plants were grown under the control conditions (Fig 8A and 8C); the first high group of differentially expressed pseudogenes (blue lines) has the highest methylation ratio of mCHG sequence context. However, when plants were treated with salt, the level of DNA methylation was maintained, but it was mainly associated with a lower level of pseudogene expression (Fig 8B and 8D). As an example, the third group of expression level (green lines) showed the highest level of mCG methylation. Profiling of mCHH among pseudogene expression groups did not show a clear effect on gene expression, irrespective of the growth conditions of the plants (Fig 8E and 8F).

The effect of DMGs on gene expression

Among 6,405 DEGs (p < 0.001, FDR ≤ 0.05), 71.6% contained mCG, 79.7% contained mCHG, and 94.2% contained mCHH sites (S2 Fig). These genes were considered to be DMGs, and were recognized based on pairwise differential comparison using the t-test (p ≤ 0.05). In order to determine the effect of different methylated regions within the genes (promoter, exon, and intron) on gene expression under control and salinity conditions, a statistical correlation between DMGs and their gene expression values was studied using the measure of Pearson’s linear correlation. Specifically, the correlation was created between the log2 of mRNA expression (FPKM values) and the log2 of the methylation ratio. It was observed that the impact of DNA methylation on gene expression was modulated by the environmental conditions, in which there was generally an opposite effect of the DNA methylation on gene expression as a result of the changes in the environmental conditions (control versus salinity).

The results revealed that there was an insignificant negative correlation between mCG methylation levels located within the promoters and gene expression when plants were grown under both control and salinity conditions. However, there was a significant (p ≤ 0.05) positive correlation between the methylation level at the mCHG and mCHH sites located within the promoters and gene expression when plants were grown under control conditions. On the contrary, when plants were grown under salinity conditions, this relationship was significantly (p ≤ 0.01) negative (Table 4).

thumbnail
Table 4. Correlation between DMGs and gene expression based on the Pearson’s linear correlation.

Asterisk stands for statistical significance (* if p ≤ 0.05 and ** if p ≤ 0.01).

https://doi.org/10.1371/journal.pone.0191492.t004

The methylation profile of the mCG sequence contexts located within the exon regions of the DMGs showed a trend of a negative correlation with the gene expression when plants were grown under control conditions, but this correlation was reversed when the plants were grown under salinity conditions. The correlation between the methylation level of the mCHG and mCHH located within the exons and the gene expression was significantly (p ≤ 0.01) positive when plants were grown under control conditions; however, this relationship was significantly (p ≤ 0.01) negative when plants were exposed to salinity.

The methylation profile of the mCG sequence contexts located within the intron regions of the DMGs showed a significant (p ≤ 0.01) negative correlation with the gene expression when plants were grown under control conditions, but this relationship was insignificantly positive when the plants were grown under salinity conditions. The correlation between the methylation levels of the mCHG and mCHH encountered within the introns and the gene expression was significantly (p ≤ 0.01) positive when plants were grown under control conditions; however, this relationship was significantly (p ≤ 0.01) negative when plants were exposed to salinity.

When the methylation profiles of the whole gene regions were correlated with their gene expression levels, the statistical analysis showed that there was a significant (p ≤ 0.01) negative relationship between the DNA methylation levels of mCG and gene expression when plants were grown under control conditions. The same negative relationship was found in mCHG and mCHH when plants were exposed to salt stress (Table 4). On the other hand, when plants were exposed salt stress, there was a significant (p ≤ 0.01) positive correlation between the levels of DNA methylation at the mCG sites and gene expression. The same positive relationships were observed between the methylation levels of mCHG and mCHH and gene expression when plants were grown under control conditions.

Functional annotation of DMGs

In order to obtain a comprehensive idea of the biological functions of the DMGs in the date palm, the protein-coding DNA sequence of the DMGs from the three sequence contexts were separately annotated and designated gene ontologies (GOs).

Each sequence context was separately analyzed based on three GOs: biological process (BP), molecular function, and cellular component (S3S5 Figs). The functional annotation of the DMGs included one of the three sequence contexts: mCG (S3 Fig, S2 Table), mCHG (S4 Fig, S3 Table) and mCHH (S5 Fig, S4 Table). The results showed that the DMGs were involved in the biological processes of oxidation-reduction, regulation of transcription, protein phosphorylation, proteolysis, transport and methylation. The DMGs were found to be also involved in the integral components of the membrane, nucleus, mitochondrion, cytoplasm, chloroplast, plastids, plasma membrane and cytosolic cellular components. These genes were also encoded for molecular functions, including the binding of proteins, ATP, DNA, ions, cations and nucleic acids.

Differential enrichment analysis based on Fisher’s exact test (p < 0.05) was carried out in order to identify the differential presented GOs. The analysis showed that downregulated genes containing mCG were enriched in salinity-responsive categories, such as cellular ion and cation homeostasis, carbohydrates transmembrane transport, cation and anion binding, protein kinase activity, protein dephosphorylation and ion binding. These categories were also enriched with GOs acting on plasma membranes and extracellular regions (Fig 9, S5 Table).

thumbnail
Fig 9. Functional enrichment analysis of significant DMGs (p ≤ 0.05) including mCG based on the categories biological process, cellular component and molecular function.

Fisher’s exact test analysis considered the upregulated DMGs as a test group.

https://doi.org/10.1371/journal.pone.0191492.g009

Upregulated genes containing mCHG were enriched with DNA and protein methylation, vesicle-mediated transport, intracellular transport, hydrogen and proton transport, vacuole organization, protein transporter, N-methyltransferase, 3'-5' exonuclease. These GOs were found to be involved in the whole membrane, the bounding membrane of organelles, endomembrane systems and chloroplast parts (S4 and S6 Figs, S6 Table). The upregulated genes containing mCHH were abundantly enriched with cellular potassium ion transport, vacuolar transport, DNA methylation, amino acid activation, protein, proton, and hydrogen ion transmembrane transport, solute:cation and solute: proton antiporter activity, hydrogen ion transmembrane transporter activity, protein transporter activity, 3'-5'-exodeoxyribonuclease activity and monovalent inorganic cation transmembrane transporter activity (S5 and S7 Figs, S7 Table).

Functional annotation analysis of the DMGs also revealed that the DMGs encoded 620 enzymes. The list included 199 DMGs that embraced mCG sites, 206 that embraced mCHG sites, and 215 that embraced mCHH sites. These enzymes were mapped within 318 unique pathways (S8S10 Tables). This includes those pathways that are involved in cysteine and methionine metabolic pathway (S8 Fig).

Discussion

DNA methylation is an epigenetic mechanism that regulates key cellular processes in eukaryotes [1113, 47]. Recent studies that focused on methylome alterations in response to abiotic stresses have shown variable patterns of DNA methylation among different plant species. However, these patterns are inconsistent, especially in perennial species [23]. Changes in DNA methylation are usually associated with alterations in plant growth and development [22]. In the current study, salinity showed a negative effect on the growth of the leaves and roots of the treated seedling. An abnormal development of the leaves was also associated with a reduction in photosynthetic performance. This result is consistent with previously published research on date palms in which transcriptome analysis showed a negative effect of salinity on the expression of some key genes associated with photosynthesis [48, 49]. These changes in gene expression could be due to alterations in DNA methylation in the leaf tissues.

Root systems represent the very first line of defense against salinity stress damage. Therefore, epigenetic alterations associated with physiological changes are important sources of information needed to decode the salinity tolerance mechanisms in plants. Impaired root growth and development in date palms due to salinity could be the result of alterations in DNA methylation, because it was previously documented that changes in DNA methylation caused abnormal root development in Arabidopsis [50].

In this project, a high throughput analysis of single-base resolution DNA methylome as well as transcriptome analysis was used to investigate the regulation of salinity-responsive genes by DNA methylation in date palm roots. WGBS is the most versatile method used to obtain a comprehensive idea about DNA methylation, especially for studies aiming at constructing the first methylome map of an organism such as the date palm [51].

In this project, DNA of each sample was extracted from a pool of root tissues collected from eight different plants. Using of this pool of plant tissues in the WGBS analysis can minimize the epigenetic variations, which may occur between the individual seedlings of the same treatment. Each pool of these DNA was considered as a biological sample. Obviously, including more experimental and experimental replicates in the WGBS analysis may improve the quality of the results however, WGBS sequencing is expensive and few laboratories can afford it. Although we consider the results obtained from this project as the first draft methylome map for the date palm, these results are novel and largely consistent with the known facts about methylome and transcriptome relationships.

The unique reads obtained by bisulfite-converted sequencing of the control and the salinity-treated samples had a mapping efficiency of around 41%. This value is consistent with the range of the typical mapping efficiency of bisulfite-converted DNA, which is 30–50% [52]. In addition, our methylome sequencing analysis results have an average read depth of 17X per each DNA strand (34X total coverage), a coverage which was previously recommended by the National Institutes of Health (NIH) Roadmap Epigenomics Project as sufficient [53]. The results presented in this study were obtained despite the fact that the date palm genome is not yet completely sequenced or assembled. It is known that the performance of the WGBS technique mainly depends on the quality of the reference genome [53, 54] and the alignment efficiency of the bisulfite-converted sequencing reads to the reference genome [55].

Global quantification of 5-mCs using mass spectrometry showed a similar amount to those previously found in other plant species, which ranged from 5% to 30% of total cytosine [56, 57]. Mass spectrometry also showed a substantial increase in the DNA methylation level after exposure to salt stress. Frequently, overmethylation is associated with the salt stress and salinity tolerance phenotype in plants. For example, when Arabidopsis is exposed to salinity, the global DNA methylation level increases [58]. Such an increase was also detected in Medicago species [31] and Mesembryanthemum crystallinum [59]. Moreover, salt-tolerant wheat cultivar are more methylated than their susceptible counterparts [60]. On the other hand, previously published reports showed that global DNA demethylation is associated with salinity tolerance in plants [61]. Accordingly, the results obtained in this study are not surprising, because the effect of salinity on DNA methylation varies based on plant species, genotype and the organ within the same plant [6264].

Mass spectrometer analysis could not detect 5HmdC in the date palm genome. This is consistent with the fact that only trace amounts of 5HmdC were detected in Arabidopsis and these amounts could be the result of random oxygen species reactivity with 5-mCs [65]. In mammals, the ten eleven translocation 1 (TET1) family of enzymes oxidizes 5mdC to 5HmdC [9, 66, 67]; however, none of this family of related enzymes were previously reported in plant species. Consequently, 5HmdC was not detected in a remarkable amount [68].

The methylated mCHH was the most abundant type among the total significant DMRs, followed by mCG methylation (Table 2). In fact, the methylation in the mCG sequence context is the most persistent and stable in the genome [8, 46]. The differential cytosine methylation analysis showed that the date palm genome was hypermethylated in response to salt stress at these DMRs of a significant (p < 0.05) number of reads (Table 2).

DMRs may include more than one type of methylated sequence contexts thus, the presence of overlapped methylated sequences will affect the calculations regarding the total levels of methylation within DMRs. However, since the number of hypermethylated mCHH within the DMRs is significantly high (26% higher than their hypomethylated counterparts), the chance that the overlapped methylated sequence contexts could change the conclusion that salinity increases global DNA methylation is unlikely especially because this conclusion is supported by mass spectrometry and ELISA analyses.

Gene repression can take place when the transcription factor is unable to access specific motifs located within the promoter regions [69]. An increase in the methylation level within the promoter regions is usually associated with a low or even no gene expression [69]. However, a positive correlation was previously documented between gene body methylation and gene expression [8, 18, 26]. In the date palm, mCG methylation within the promoter was associated with gene repression; however, other types of cysteine methylation located within other gene features showed a salinity stress-dependent gene expression. 5-mC bisulfite sequencing and gene expression analysis in maize seeds revealed the opposite effect of mCG methylation located within the coding regions on transcription activation, as compared to mCHG and mCHH methylation [70]. This is consistent with our results. However, when plants were exposed to salt stress, the mCHG and mCHH methylation levels were inversely correlated with gene expression. Therefore, it was clear from the results of our study that salinity alters the commonly established relationship between DNA methylation and gene expression. This is likely because under stress conditions, there are additional epigenetic factors that may be involved in gene regulation networks, such as sRNA molecules [71, 72].

It was noticed from the results that the DNA methylation level of pseudogenes was higher than those found in the protein-coding genes. This is consistent with the previous findings of the Arabidopsis methylome [13, 73]. There were some differences between the effect of DNA methylation on gene expression of the functional genes and of that of the pseudogenes. DNA methylation in the pseudogenes had a lower impact on the expression level it did in the case of the coding genes. Recently, it was found that pseudogene expression plays a crucial role in regulating the expression of their parental genes and other unrelated genes; however, the underlying mechanism is unclear [74].

Under both control and salt stress conditions, our analysis indicated that mCG and mCHG methylation around the TSS has negatively affected gene expression. However, mCG methylation around the TES was positively correlated with gene expression when plants were exposed to salinity. It was previously demonstrated that mCG methylation around the TSS and TES have a negative effect on gene expression in rice [75] and maize [70].

GO enrichment analysis showed that the upregulated DMGs of the three sequence contexts were abundantly involved in DNA methylation. This is could be evidence of an increased methylation level under salinity stress. In addition, most of the transporter and antiporter coding proteins involved were in abundantly enriched categories in the downregulated DMGs containing mCG. However, these salinity response-associated ontologies were abundantly enriched in the upregulated DMGs containing mCHG and mCHH. This finding is in agreement with the transcriptome/methylome analysis based on Pearson’s correlation, which showed that under salinity stress there was a negative correlation between the mCG methylation ratio and expression level, while there was a positive correlation between mCHG and mCHH methylation ratios and expression levels.

Mapping DMGs on the KEGG maps revealed that the cysteine and methionine metabolic pathway was enriched with enzymes coded by methylated genes in three sequence contexts. Interestingly, it was previously observed that this pathway is actively involved in abiotic stress response [76]. Sulfur-containing amino acids (cysteine and methionine) are sensitive to ROS-mediated oxidation. For example, oxidation changes methionine to S- and R-methionine sulfoxide diastereoisomers. Two types of methionine sulfoxide reductases, MSRA and MSRB, can specifically reverse this reaction. Based on functional characterization analysis in rice, overexpression of the gene coding MSRA (OsMSRA4.1) showed an enhanced growth under salt stress[76].

In conclusion, this initial study identified widespread salinity-induced changes in 5mC methylation and transcriptome abundance in the roots of the date palm. The methylation changes that occurred within various regions of the gene showed an impact on gene expression. While the DNA methylation demonstrated the common regulation pattern known in eukaryotes on gene expression when plants grow under optimum environmental conditions, this study revealed that when plants grow under salinity stress, DNA methylation is not the primary agent that controls gene expression in date palm roots. The results obtained from this study could highlight the importance of DNA methylation in epigenetic remodeling under stress conditions, and the information encoded by the methylome could be beneficial for crop improvement strategies.

Supporting information

S1 Fig. Clustering heatmap analysis of the top 100 differentially methylated sites (DMSs) (p < 0.05).

Dendogram of hierarchical clustering was obtained based on DNA methylation ratio of mCG, mCHG and mCHH sequence contexts. Red and yellow color scale represents individual sites that are 0% and 100% methylated, respectively.

https://doi.org/10.1371/journal.pone.0191492.s001

(TIF)

S2 Fig. Classification of the differentially expressed genes (DEGs) (6,405 genes) which were recognized as differentially methylated genes (DMGs) in response to salinity stress based on their methylation status.

https://doi.org/10.1371/journal.pone.0191492.s002

(TIF)

S3 Fig. Gene ontology (GO) analysis of DMGs for each sequence context mCG of date palm epigenome in response to salinity stress, based on three sets of ontologies: (A), molecular function (B), and cellular component (C) against the GO databases.

https://doi.org/10.1371/journal.pone.0191492.s003

(TIF)

S4 Fig. Gene ontology analysis of DMGs for each sequence contexts mCHG of date palm epigenome in response to salinity stress, based on three sets of ontologies; (A), molecular function (B), and cellular component (C) against the GO databases.

https://doi.org/10.1371/journal.pone.0191492.s004

(TIF)

S5 Fig. Gene ontology analysis of DMGs for each sequence contexts mCHH of date palm epigenome in response to salinity stress, based on three sets of ontologies; (A), molecular function (B), and cellular component (C) against the GO databases.

https://doi.org/10.1371/journal.pone.0191492.s005

(TIF)

S6 Fig. Functional enrichment analysis of significant DMGs (p ≤ 0.05) including mCHG based on the categories biological process, cellular component and molecular function.

Fisher’s exact test analysis considered the upregulated DMGs as a test group.

https://doi.org/10.1371/journal.pone.0191492.s006

(TIF)

S7 Fig. Functional enrichment analysis of significant DMGs (p ≤ 0.05) including mCHH based on the categories biological process, cellular component and molecular function.

Fisher’s exact test analysis considered the upregulated DMGs as a test group.

https://doi.org/10.1371/journal.pone.0191492.s007

(TIF)

S8 Fig. Schematic representation of enzymes coded by mCHG DMGs.

These enzymes might function in the cysteine and methionine metabolic pathway.

https://doi.org/10.1371/journal.pone.0191492.s008

(TIF)

S1 Table. Differentially expressed genes (DEGs) (p < 0.001, FDR ≤ 0.05) under salinity stress.

https://doi.org/10.1371/journal.pone.0191492.s009

(XLSX)

S2 Table. Functional annotation and gene ontology of DMGs for the mCG sequence context.

https://doi.org/10.1371/journal.pone.0191492.s010

(XLSX)

S3 Table. Functional annotation and gene ontology of DMGs for the mCHG sequence context.

https://doi.org/10.1371/journal.pone.0191492.s011

(XLSX)

S4 Table. Functional annotation and gene ontology of DMGs for the mCHH sequence context.

https://doi.org/10.1371/journal.pone.0191492.s012

(XLSX)

S5 Table. Functional enrichment analysis of significant DMGs (p ≤ 0.05) for the mCG sequence context.

Fisher’s exact test analysis considered the upregulated DMGs as a test group.

https://doi.org/10.1371/journal.pone.0191492.s013

(XLSX)

S6 Table. Functional enrichment analysis of significant DMGs (p ≤ 0.05) for the mCHG sequence context.

Fisher’s exact test analysis considered the upregulated DMGs as a test group.

https://doi.org/10.1371/journal.pone.0191492.s014

(XLSX)

S7 Table. Functional enrichment analysis of significant DMGs (p ≤ 0.05) for the mCHH sequence context.

Fisher’s exact test analysis considered the upregulated DMGs as a test group.

https://doi.org/10.1371/journal.pone.0191492.s015

(XLSX)

S8 Table. Mapping of DMGs for the mCG sequence context on the metabolic pathways based on annotated coding enzymes.

https://doi.org/10.1371/journal.pone.0191492.s016

(XLSX)

S9 Table. Mapping of DMGs for the mCHG sequence context on the metabolic pathways based on annotated coding enzymes.

https://doi.org/10.1371/journal.pone.0191492.s017

(XLSX)

S10 Table. Mapping of DMGs for the mCHH sequence context on the metabolic pathways based on annotated coding enzymes.

https://doi.org/10.1371/journal.pone.0191492.s018

(XLSX)

Acknowledgments

This work was supported by a generous grant from the Research Council of Oman (TRC) number 151.

References

  1. 1. Stepien P, Johnson GN. Contrasting responses of photosynthesis to salt stress in the glycophyte Arabidopsis and the halophyte Thellungiella: role of the plastid terminal oxidase as an alternative electron sink. Plant Physiol. 2009;149(2):1154–65. pmid:19052149
  2. 2. Wang MC, Peng ZY, Li CL, Li F, Liu C, Xia GM. Proteomic analysis on a high salt tolerance introgression strain of Triticum aestivum/Thinopyrum ponticum. Proteomics. 2008;8(7):1470–89. pmid:18383010
  3. 3. Assaha DVM, Ueda A, Saneoka H, Al-Yahyai R, Yaish MW. The Role of Na+ and K+ Transporters in Salt Stress Adaptation in Glycophytes. Front Physiol. 2017;8:509. pmid:28769821
  4. 4. Munns R, Tester M. Mechanisms of salinity tolerance. Annu Rev Plant Biol. 2008;59:651–81. pmid:18444910
  5. 5. Zhang H, Han B, Wang T, Chen S, Li H, Zhang Y, et al. Mechanisms of plant salt response: insights from proteomics. Journal of proteome research. 2011;11(1):49–67. pmid:22017755
  6. 6. Kumari A, Das P, Parida AK, Agarwal PK. Proteomics, metabolomics, and ionomics perspectives of salinity tolerance in halophytes. Front Plant Sci. 2015;6.
  7. 7. Jablonka E, Raz G. Transgenerational epigenetic inheritance: prevalence, mechanisms, and implications for the study of heredity and evolution. Quarterly Review of Biology. 2009;84(2):131–76. pmid:19606595
  8. 8. Cokus SJ, Feng S, Zhang X, Chen Z, Merriman B, Haudenschild CD, et al. Shotgun bisulphite sequencing of the Arabidopsis genome reveals DNA methylation patterning. Nature. 2008;452(7184):215–9. pmid:18278030
  9. 9. Saze H, Tsugane K, Kanno T, Nishimura T. DNA Methylation in Plants: Relationship to Small RNAs and Histone Modifications, and Functions in Transposon Inactivation. Plant and Cell Physiol. 2012;53(5):766–84. pmid:22302712
  10. 10. Hardcastle TJ. High-throughput sequencing of cytosine methylation in plant DNA. Plant Methods. 2013;9. pmid:23758782
  11. 11. Li E, Bestor TH, Jaenisch R. Targeted mutation of the DNA methyltransferase gene results in embryonic lethality. Cell. 1992;69(6):915–26. pmid:1606615
  12. 12. Lippman Z, Gendrel AV, Black M, Vaughn MW, Dedhia N, McCombie WR, et al. Role of transposable elements in heterochromatin and epigenetic control. Nature. 2004;430(6998):471–6. pmid:15269773
  13. 13. Zhang XY, Yazaki J, Sundaresan A, Cokus S, Chan SWL, Chen HM, et al. Genome-wide high-resolution mapping and functional analysis of DNA methylation in Arabidopsis. Cell. 2006;126(6):1189–201. pmid:16949657
  14. 14. Yaish M. Epigenetic Modifications Associated with Abiotic and Biotic Stresses in Plants: An Implication for Understanding Plant Evolution. Front Plant Sci. 2017;8:1983. pmid:29209347
  15. 15. Stroud H, Do T, Du J, Zhong X, Feng S, Johnson L, et al. Non-CG methylation patterns shape the epigenetic landscape in Arabidopsis. Nat Struct Mol Biol. 2014;21(1):64–72. pmid:24336224
  16. 16. Lindroth AM, Cao X, Jackson JP, Zilberman D, McCallum CM, Henikoff S, et al. Requirement of CHROMOMETHYLASE3 for maintenance of CpXpG methylation. Science. 2001;292(5524):2077–80. pmid:11349138
  17. 17. He Y, Ecker JR. Non-CG methylation in the human genome. Annu. Rev. Genomics Hum. Genet. 2015;16:55–77. pmid:26077819
  18. 18. Feng SH, Cokus SJ, Zhang XY, Chen PY, Bostick M, Goll MG, et al. Conservation and divergence of methylation patterning in plants and animals. Proc Natl Acad Sci, USA. 2010;107(19):8689–94. pmid:20395551
  19. 19. Law JA, Jacobsen SE. Establishing, maintaining and modifying DNA methylation patterns in plants and animals. Nature Rev Genet. 2010;11(3):204–20. pmid:20142834
  20. 20. Yaish MWF, Peng M, Rothstein SJ. AtMBD9 modulates Arabidopsis development through the dual epigenetic pathways of DNA methylation and histone acetylation. Plant J. 2009;59(1):123–35. pmid:19419532
  21. 21. Boyko A, Blevins T, Yao YL, Golubov A, Bilichak A, Ilnytskyy Y, et al. Transgenerational Adaptation of Arabidopsis to Stress Requires DNA Methylation and the Function of Dicer-Like Proteins. Plos One. 2010;5(3). pmid:20209086
  22. 22. Yaish MW, Colasanti J, Rothstein SJ. The role of epigenetic processes in controlling flowering time in plants exposed to stress. J Exp Bot. 2011;62(11):3727–35. pmid:21633082
  23. 23. Song Y, Ci D, Tian M, Zhang D. Stable methylation of a non-coding RNA gene regulates gene expression in response to abiotic stress in Populus simonii. J Exp Bot. 2015:erv543.
  24. 24. Wang J, Marowsky NC, Fan C. Divergence of gene body DNA methylation and evolution of plant duplicate genes. PLoS One. 2014;9(10):e110357. pmid:25310342
  25. 25. Bell AC, Felsenfeld G. Methylation of a CTCF-dependent boundary controls imprinted expression of the Igf2 gene. Nature. 2000;405(6785):482–5. pmid:10839546
  26. 26. Lister R, Pelizzola M, Dowen RH, Hawkins RD, Hon G, Tonti-Filippini J, et al. Human DNA methylomes at base resolution show widespread epigenomic differences. Nature. 2009;462(7271):315–22. pmid:19829295
  27. 27. Yaish MW. Proline accumulation is a general response to abiotic stress in the date palm tree (Phoenix dactylifera L.). Genet Mol Biol. 2015;14(3):9943–50. pmid:26345930
  28. 28. Yaish MW, Kumar PP. Salt tolerance research in date palm tree (Phoenix dactylifera L.), past, present, and future perspectives. Front Plant Sci. 2015;6. pmid:26042137
  29. 29. Yaish MW, Sunkar R, Zheng Y, Ji B, Al-Yahyai R, Farooq SA. A genome-wide identification of the miRNAome in response to salinity stress in date palm (Phoenix dactylifera L.). Front Plant Sci. 2015;6. pmid:26594218
  30. 30. Al Kharusi L, Assaha DV, Al-Yahyai R, Yaish MW. Screening of Date Palm (Phoenix dactylifera L.) Cultivars for Salinity Tolerance. Forests. 2017;8(4):136.
  31. 31. Al-Lawati A, Al-Bahry S, Victor R, Al-Lawati AH, Yaish MW. Salt stress alters DNA methylation levels in alfalfa (Medicago spp). Genet Mol Biol. 2016;15(1). pmid:26985924
  32. 32. Pandey G, Yadav CB, Sahu PP, Muthamilarasan M, Prasad M. Salinity induced differential methylation patterns in contrasting cultivars of foxtail millet (Setaria italica L.). Plant Cell Rep. 2016. pmid:27999979
  33. 33. Yaish MW, Peng M, Rothstein SJ. Global DNA methylation analysis using methyl-sensitive amplification polymorphism (MSAP). Methods Mol Biol. 2014;1062:285–98. pmid:24057373
  34. 34. Khan S, Qureshi MI, Kamaluddin , Alam T, Abdin MZ. Protocol for isolation of genomic DNA from dry and fresh roots of medicinal plants suitable for RAPD and restriction digestion. Afr J Biotechnol. 2007;6(3):175–8.
  35. 35. Xiao Y, Yang YD, Cao HX, Fan HK, Ma ZL, Lei XT, et al. Efficient isolation of high quality RNA from tropical palms for RNA-seq analysis. Plant Omics. 2012;5(6):584–9.
  36. 36. Al-Harrasi I, Al-Yahyai R, Yaish MW. Detection of Differential DNA MethylationDNA methylationBisulfite sequence Under Stress Conditions Using Bisulfite Sequence Analysis. In: Sunkar R, editor. Methods Molecular Biology. 1631. 2nd ed. Humana Press: Springer; 2017. 1631:121–137.
  37. 37. Urich MA, Nery JR, Lister R, Schmitz RJ, Ecker JR. MethylC-seq library preparation for base-resolution whole-genome bisulfite sequencing. Nat Protoc. 2015;10(3). pmid:25692984
  38. 38. Krueger F, Andrews SR. Bismark: a flexible aligner and methylation caller for Bisulfite-Seq applications. Bioinformatics. 2011;27(11):1571–2. pmid:21493656
  39. 39. Quinlan AR, Hall IM. BEDTools: a flexible suite of utilities for comparing genomic features. Bioinformatics. 2010;26(6):841–2. pmid:20110278
  40. 40. Yang M, Zhang X, Liu G, Yin Y, Chen K, Yun Q, et al. The complete chloroplast genome sequence of date palm (Phoenix dactylifera L.). PLoS One. 2010;5(9):e12762. pmid:20856810
  41. 41. Al-Mssallem IS, Hu SN, Zhang XW, Lin Q, Liu WF, Tan J, et al. Genome sequence of the date palm Phoenix dactylifera L. Nature Commun. 2013;4:9. pmid:23917264
  42. 42. Götz S, García-Gómez JM, Terol J, Williams TD, Nagaraj SH, Nueda MJ, et al. High-throughput functional annotation and data mining with the Blast2GO suite. Nucleic Acids Res. 2008;36(10):3420–35. pmid:18445632
  43. 43. Kanehisa M, Goto S, Kawashima S, Okuno Y, Hattori M. The KEGG resource for deciphering the genome. Nucleic Acids Res. 2004;32(suppl 1):D277–D80.
  44. 44. Song LG, James SR, Kazim L, Karpf AR. Specific method for the determination of genomic DNA methylation by liquid chromatography-electrospray ionization tandem mass spectrometry. Anal Chem. 2005;77(2):504–10. pmid:15649046
  45. 45. Fojtová M, Kovařıík A, Matyášek R. Cytosine methylation of plastid genome in higher plants. Fact or artefact? Plant Sci. 2001;160(4):585–93. pmid:11448733
  46. 46. Lister R, O'Malley RC, Tonti-Filippini J, Gregory BD, Berry CC, Millar AH, et al. Highly integrated single-base resolution maps of the epigenome in Arabidopsis. Cell. 2008;133(3):523–36. pmid:18423832
  47. 47. Bestor TH. The DNA methyltransferases of mammals. Hum Mol Gen. 2000;9(16):2395–402. pmid:11005794
  48. 48. Yaish MW, Patankar HV, Assaha DV, Zheng Y, Al-Yahyai R, Sunkar R. Genome-wide expression profiling in leaves and roots of date palm (Phoenix dactylifera L.) exposed to salinity. BMC Genomics. 2017;18(1):246. pmid:28330456
  49. 49. Sperling O, Lazarovitch N, Schwartz A, Shapira O. Effects of high salinity irrigation on growth, gas-exchange, and photoprotection in date palms (Phoenix dactylifera L., cv. Medjool). Environ Exper Bot. 2014;99:100–9.
  50. 50. Finnegan EJ, Peacock WJ, Dennis ES. Reduced DNA methylation in Arabidopsis thaliana results in abnormal plant development. Proc Natl Acad Sci, USA. 1996;93(16):8449–54. pmid:8710891
  51. 51. Yong WS, Hsu FM, Chen PY. Profiling genome-wide DNA methylation. Epigenetics Chromatin. 2016;9:16.
  52. 52. Tran H, Porter J, Sun M-a, Xie H, Zhang L. Objective and Comprehensive Evaluation of Bisulfite Short Read Mapping Tools. Advances in Bioinformatics. 2014;2014:11. pmid:24839440
  53. 53. Ziller MJ, Hansen KD, Meissner A, Aryee MJ. Coverage recommendations for methylation analysis by whole-genome bisulfite sequencing. Nature Methods. 2015;12(3):230–+. pmid:25362363
  54. 54. Stevens M, Cheng JB, Li D, Xie M, Hong C, Maire CL, et al. Estimating absolute methylation levels at single-CpG resolution from methylation enrichment and restriction enzyme sequencing methods. Genome Res. 2013;23(9):1541–53. pmid:23804401
  55. 55. Krueger F, Kreck B, Franke A, Andrews SR. DNA methylome analysis using short bisulfite sequencing data. Nature Methods. 2012;9(2):145–51. pmid:22290186
  56. 56. Wagner I, Capesius I. Determination of 5-methylcytosine from plant DNA by high-performance liquid chromatography. Biochimica et Biophysica Acta (BBA)-Nucleic Acids and Protein Synthesis. 1981;654(1):52–6.
  57. 57. Leutwiler LS, Hough-Evans BR, Meyerowitz EM. The DNA of Arabidopsis thaliana. Mol Gen Genet. 1984;194(1–2):15–23.
  58. 58. Boyko A, Blevins T, Yao Y, Golubov A, Bilichak A, Ilnytskyy Y, et al. Transgenerational Adaptation of Arabidopsis to Stress Requires DNA Methylation and the Function of Dicer-Like Proteins. PLoS One. 2010;5(3):e9514. pmid:20209086
  59. 59. Dyachenko OV, Zakharchenko NS, Shevchuk TV, Bohnert HJ, Cushman JC, Buryanov YI. Effect of hypermethylation of CCWGG sequences in DNA of Mesembryanthemum crystallinum plants on their adaptation to salt stress. Biochemistry (Moscow). 2006;71(4):461–5.
  60. 60. Zhong L, Wang J-B. The Role of DNA Hypermethylation in Salt Resistence of Triticum aestivum L. J Wuhan Botanical Research. 2007;1:019.
  61. 61. Viggiano L, de Pinto MC. Dynamic DNA Methylation Patterns in Stress Response. Plant Epigenetics: RNA Technologies. Springer, Cham 2017. p. 281–302. https://doi.org/10.1007/978-3-319-55520-1_15
  62. 62. Secco D, Wang C, Shou H, Schultz MD, Chiarenza S, Nussaume L, et al. Stress induced gene expression drives transient DNA methylation changes at adjacent repetitive elements. Elife. 2015;4:e09343.
  63. 63. Kumar S, Beena AS, Awana M, Singh A. Salt-Induced Tissue-Specific Cytosine Methylation Downregulates Expression of HKT Genes in Contrasting Wheat (Triticum aestivum L.) Genotypes. DNA and Cell Biol. 2017.
  64. 64. Karan R, DeLeon T, Biradar H, Subudhi PK. Salt stress induced variation in DNA methylation pattern and its influence on gene expression in contrasting rice genotypes. PloS One. 2012;7(6):e40203. pmid:22761959
  65. 65. Liu S, Dunwell TL, Pfeifer GP, Dunwell JM, Ullah I, Wang Y. Detection of oxidation products of 5-methyl-2′-deoxycytidine in Arabidopsis DNA. PLoS One. 2013;8(12):e84620. pmid:24391970
  66. 66. Tahiliani M, Koh KP, Shen Y, Pastor WA, Bandukwala H, Brudno Y, et al. Conversion of 5-methylcytosine to 5-hydroxymethylcytosine in mammalian DNA by MLL partner TET1. Science. 2009;324(5929):930–5. pmid:19372391
  67. 67. Ito S, Shen L, Dai Q, Wu SC, Collins LB, Swenberg JA, et al. Tet proteins can convert 5-methylcytosine to 5-formylcytosine and 5-carboxylcytosine. Science. 2011;333(6047):1300–3. pmid:21778364
  68. 68. Erdmann RM, Souza AL, Clish CB, Gehring M. 5-Hydroxymethylcytosine Is Not Present in Appreciable Quantities in Arabidopsis DNA. G3. 2014. pmid:25380728
  69. 69. Suzuki MM, Bird A. DNA methylation landscapes: provocative insights from epigenomics. Nature Rev Genet. 2008;9(6):465–76. pmid:18463664
  70. 70. Lu XD, Wang WX, Ren W, Chai ZG, Guo WZ, Chen RM, et al. Genome-Wide Epigenetic Regulation of Gene Transcription in Maize Seeds. Plos One. 2015;10(10). pmid:26469520
  71. 71. Ganie SA, Dey N, Mondal TK. Promoter methylation regulates the abundance of osa-miR393a in contrasting rice genotypes under salinity stress. Funct Integr Genomics. 2016;16(1):1–11. pmid:26319531
  72. 72. Sunkar R, Li Y-F, Jagadeeswaran G. Functions of microRNAs in plant stress responses. Trends Plant Sci. 2012;17(4):196–203. pmid:22365280
  73. 73. Zilberman D, Gehring M, Tran RK, Ballinger T, Henikoff S. Genome-wide analysis of Arabidopsis thaliana DNA methylation uncovers an interdependence between methylation and transcription. Nature Genet. 2007;39(1):61–9. pmid:17128275
  74. 74. Poliseno L. Pseudogenes: Newly Discovered Players in Human Cancer. Sci Signal. 2012;5(242). pmid:22990117
  75. 75. Zemach A, Kim MY, Silva P, Rodrigues JA, Dotson B, Brooks MD, et al. Local DNA hypomethylation activates genes in rice endosperm. Proc Natl Acad Sci, USA 2010;107(43):18729–34. pmid:20937895
  76. 76. Guo X, Wu Y, Wang Y, Chen Y, Chu C. OsMSRA4. 1 and OsMSRB1. 1, two rice plastidial methionine sulfoxide reductases, are involved in abiotic stress responses. Planta. 2009;230(1):227–38. pmid:19415325