Water Stress Responsive Differential Methylation of Organellar Genomes of Zea mays Z59

DNA methylation is an important epigenetic change 
affecting gene expression in plants in both normal and stress conditions. The 
organelles, mitochondria and chloroplast play a significant role in sensing and 
initiating stress response. In this study, we report the methylation pattern in 
chloroplast and mitochondrial genomes in irrigated and water stressed 
conditions and its relationship with gene expression of a drought tolerant Zea mays cultivar, Z59. Whole genome bisulfite sequencing was done to analyze the 
pattern of methylation in both the conditions. Mapping of bisulfite reads to 
B73 reference of mitochondrial and chloroplast genomes showed hypomethylation 
in water stressed plants when compared to irrigated plants. Sliding window 
approach to the methylation count data showed highest peak at 419,800 to 
420,800 bp region in mitochondria and at 36,900 to 37,900 bp region in 
chloroplast genomes in both samples. Annotation of the methylated genomes 
showed that, genes related to photosystem I & II in chloroplast and nad4 
gene in mitochondria were hypo methylated in the water stressed sample. RNA-seq 
analysis of transcriptomics reads mapped to the same 
reference showed regulation of rps3, rps2A, ccmFC, atp1 and many 
uncharacterized genes in mitochondria and psbA, psbD, psbc, psaA, and atpA, 
genes in chloroplast.


Introduction
Epigenetics is a study on genetic modifications, which occur in a cell without How to cite this paper: Bhanu, B.D., Ulaganathan, K. and Shanker, A.K. (2020) Water Stress Responsive Differential Methylation of Organellar Genomes of Zea mays Z59. American Journal of Plant Sciences, any change in the DNA sequence. Epigenetic regulation influences gene expression by turning its switch on or off. They contribute mainly to DNA methylation, histone modifications and RNA processing. The role of epigenetic regulation is widely studied in various disciplines. Its role in abiotic stress tolerance is a challenge for plant researchers.
Drought is one of the major abiotic stresses and an important limiting factor of rain-fed agriculture. Evidences on genome wide methylation patterns suggest hypo [1] or hyper-methylation [2] as response to stress and thus regulating gene expression. Recently, with the advent of next generation sequencing, DNA methylation is detected using sodium bisulfite conversion of cytosines to uracils, where 5'-methylcytosines remain unaltered. Bisulfite-treated DNA is then, amplified and sequenced [3]. Genome wide methylation studies on various model plants are headway to address drought tolerance but there is no emphasis on the organellar genomes.
Chloroplast and mitochondria are crucial support to energy metabolism occurring in a plant. They respond in an effective, energy-efficient way to combat continuous changes in the environment. They play a significant role in cellular metabolism and therefore, are expected to be first in stress recognition [4]. In the year 1998, Simkova [5] had studied the methylation of Mitochondrial DNA in carrot using Southern hybridization and in 2012. Huang et al. [6] reported their study on mitochondrial and chloroplast DNA in the woody plant, Sequoia sempervirens using HPLC. In the recent years, tremendous data have been generated with respect to mitochondrial DNA methylation in humans. Mawlood et al. (2016) [7] suggest mtDNA methylation as biological marker for forensic age-prediction and health status measurement. A Study by Sun et al., (2018) [8] on tumorigenesis highlight the influences that the nuclear and mitochondrial genomes have in setting mtDNA methylation patterns to regulate mtDNA copy number. Another study on distribution and dynamics of MtDNA methylation suggests its role during gametogenesis [9]. Similarly, methylation patterns in mitochondrial DNA are related to cardiovascular disease [10], Down's syndrome [11], obesity [12] etc. However, limited work has been done on DNA methylation in plants as compared to Humans in general and mtDNA methylation in particular. Feng et al. (2016) [13] have studied cadmium exposed DNA methylation pattern in rice. Another study by Kawakatsu et al. (2017) [14] involves comparison of time-series methylomes of dry and germinating seeds to publicly available seed development methylomes and has concluded the role of DNA methylation in seed dormancy. Role of DNA methylation is seen in fruit ripening [15] and abiotic stresses [16].
However, mitochondrial and chloroplast DNA methylation has been faintly studied in plants. Very [19]. In this study, we have worked on the methylation pattern of mitochondria and chloroplast and their relation with the gene expression (transcriptome data of the same cultivar) in response to abiotic stress, drought.

Plant Material & Growth Conditions
The seeds of Zea mays Z59 were sown in pots with uniform soil under greenhouse conditions. The growth conditions, viz., temperature, rainfall, relative humidity and other parameters were recorded for both the dry and wet seasons (Supplementary File 2). Pots were well-watered until the reproductive stage and drought was imposed from this stage till maturity by withholding irrigation for half of the plants (16 out of 32). After the onset of drought, the leaves were dark adapted for 15 mins and the readings were measured as OJIP parameters and Performance Index, PI (abs) using portable fluorometer, Handy PEA+ (Hansatech instruments ltd., Norfolk) at regular intervals. Performance Index (PI) includes three independent parameters: 1) density of fully active reaction centers (RCs); 2) efficiency of electron movement by trapped exciton into the electron transport chain beyond the QA; and 3) the probability that an absorbed photon will be trapped by RCs. PI (abs) is calculated by the formula where: F 0 means fluorescence intensity at 50 μs, F J is fluorescence intensity at the J step (at 2 ms), F M represents maximal fluorescence intensity, V J is relative variable fluorescence at 2 ms calculated as , M 0 represents initial slope of fluorescence kinetics, which can be derived from the equation: The extent of drought stress was observed by measuring relative water content and soil volumetric water content. The relative water content was measured from the third leaf at regular intervals and was calculated by the formula [20] (fresh weight − dry weight)/(turgid weight − dry weight) × 100. The soil volumetric water content was measured with Field Scout TDR 350 probe (Specturm Technologies Inc.). DNA was isolated from the leaf tissue (lamina, topmost leaf) from both the sets of irrigated and water-stressed plants and was subjected to bisulfite sequencing to study the methylation pattern during water stress.

Library Preparation and Bisulfite Sequencing
Libraries were prepared using 200 ng of fragmented DNA (~200 bp) using NEBNext Ultra DNA Library Prep Kit for Illumina (NEB), according to the manufacturer's protocol. DNA is then, treated with sodium bisulfite followed by PCR amplification. Unmethylatedcytosines are deaminated to uracil, and read as thymine after PCR amplification and sequencing, whereas methylated cytosines are not converted and therefore, are read as cytosine after PCR amplification and sequencing [21] [22].

Assembly of the Bisulfite Treated Sequencing Reads
The raw reads were quality checked using FASTQC [23] and the adapters were trimmed using cutadapt [24]. The clean reads were aligned to the B73 reference genome (AGPv3) using the Bismark aligner [25] using default parameters. Methylated cytosines were extracted from aligned reads using the Bismark methylation extractor available in the same package. The percent methylation in each CpG was calculated by the formula (number of reads with methylated C/total reads) × 100.

Sliding Window Approach
The data was sorted chromosome-wise and the methylation data of the organellar genomes were used in further downstream analysis. A sliding window of window size 1000 bp and sliding size 100 bp were used to generate methylation count for each window. The peaks showed the hyper-methylated and hypo-methylated regions.

Methylation Data Analysis
The methylation data in the CpG context of the organellar chromosomes were taken from the total genome wide methylation data. The regions of methylation were manually mapped to the genome feature file of the reference genome for both Mitochondria and Chloroplast for annotation of the methylated regions. Gene-wise methylation count was thus generated giving information of the genes that are hyper and hypo methylated due to stress.

Co-Relation with RNA-Seq Data
RNA-seq reads are available for Zea mays Z59 [18]. The clean RNA-seq reads were aligned to the organellar chromosomes of B73 reference genome (AGPv3) using Bowtie2 aligner [26]. The mapped reads for both Mitochondrion and Chloroplast were de novo assembled separately using Trinity [27] with K-mer = 25 for both irrigated and water stressed plants. The abundance of the trinity generated transcripts was calculated using RSEM [28]. Normalization of the expression values of irrigated and water stressed transcripts of both mitochondria and chloroplast was generated separately as TPM values considering the length of the transcript and the number of reads that mapped to the transcript. The transcripts were annotated using blastx [29] and GeSeq [30]. Circular maps for both mitochondria and Chloroplast Zea mays Z59 were generated using OGDRAW [31] after ordered contigs were obtained from Progressive Mauve alignment [32]. For each gene, the expression values were plotted against the methylation count to study the co-relation of gene expression and methylation. Statistical analysis of the expression and methylation data was performed using SAS [33].

Drought Stress and Photosynthesis
The plants were monitored at regular intervals. The sliding down of relative water content and soil water volumetric content readings shows the onset of drought in water-stressed plants. The OJIP parameters and performance index, PI (abs) did not show significant difference between irrigated and water-stressed samples initially, later which there was a slight decline in water-stressed plants compared to irrigated as the drought period prolonged ( Figure 1). This indicated that photosynthesis may be affected after a prolonged period of drought in this drought tolerant cultivar.

Bisulfite Illumina Sequencing and Analysis
Whole genome bisulfite sequencing of Zea mays Z59 resulted in 260485668 and 204847228 number of paired-end reads in irrigated and water-stressed samples.
Assembly of the reads with the reference genome B73 (AGPv3) with Bismark

Co-Relation of Methylation Data with RNA-Seq Data
The RNA-seq reads of Zea mays Z59 were aligned to the organellar B73 reference genome using bowtie 2 [26]. De novo assembly of the mapped reads gener-

Data Availability
The data is submitted to NCBI and can be accessed using the following accession numbers.

Discussion
Plant vitality during water stress is an integrative study, involving network of genes and their relationships. The organellar genomes of the plants also play an essential role like the nuclear genome during abiotic stress. Our study on methylation and gene expression patterns of the organellar genes during water stress in Maize has given wealth of information for plant breeders about this drought tolerant genotype Z59. The relationship between methylation, gene expression and PI (abs) & OJIP parameters is like a triangular hypothesis. DNA methylation affects gene expression which in turn monitors PI (abs) & OJIP readings. Hypo methylation of genes in water stressed samples showed up regulation of the same genes, thereby American Journal of Plant Sciences showing tolerance to water stress which is reflected in the PI (abs) and the OJIP readings.
When leaves are dark-adapted for a minimum of 15 minutes all photochemistry in the chloroplasts stops, on subsequent exposure to light, photochemistry starts and chlorophyll fluorescence rises from a low to a high. This fast fluorescence transient rise is directly related to the Photosystem II activity. OJIP transient rise reflect the reduction reactions taking place in the PSII which involves Q A , one electron accepting bound PQ, Q B and two electron accepting mobile PQ. When all the reaction centres are open upon light interception it is termed as level O where there is no reduction of Q A , here the fluorescence intensity lasts for 10ms. After this the fluorescence transient rises from O to J which is when there is the net reduction of QA which is the stable primary electron acceptor of PS II to Q − A , this rise lasts for 2 ms. The J to I phase is due to the resultant reduced state of closed Reaction Centres (RCs) QA − QB − , QA QB2 − and QA − QB H2 which lasts for 2 -30 ms. At 300 ms the plastoquinol pool is maximally reduced and this along with maximum concentration of QA − QB2 is represented by the level P in the OJIP curve [34].
The quantitative information on the plant performance under stress was studied using PI (abs) . It reveals the functionality of both photosystem I & II [35]. Both PI (abs) and OJIP transient values are observed akin in water stressed samples as in irrigated samples probably because of the upregulation of photosystem II genes in water stressed samples. Up regulation of functional genes during stress can be attributed as one of the mechanisms to combat stress in drought tolerant genotypes like Zea mays Z59.
Joel (2013) [36] speculated hyper-methylation as an indicator of drought susceptibility and hypo-methylation for drought tolerance. We have also observed hypomethylation, overall, in water stressed sample for both Mitochondria and Chloroplast indicating up regulation of genes like rps3, rps2A, ccmFC, atp1 and many uncharacterized genes in mitochondria and psbA, psbD, psbC, psaA, at-pA, ndhk, ndhB1 & B2 genes in chloroplast. Ribosomal protein S3 and S2 (rpS3 & rps2) are known to play critical roles in ribosome biogenesis and DNA repair during increase of cellular ROS levels as response to stress [37]. Cytochrome c biogenesis (ccmFC) gene is known to form a complex with CCMH, CCMFN1 and CCMFN2 that performs the assembly of heme with c-type apocytochromes in mitochondria. According to the study by Evans et al., (2019) [38], thisgene is essential for mitochondrial function. Disruption of this might lead the mitochondrion to be non-functional. Mitochondrial ATP1 and atpA of chloroplast produces ATP from ADP and ATP is an immediately available energy form. The up regulation of this gene indicates enhancement of energy requirements as response to stress [39].
In chloroplast, the upregulated genes are psbA, psbD, psbC of photosystem II and psaA of photosystem I. The genes ndhk, ndhB1 & B2 are subunits of the chloroplast NAD (P) H dehydrogenase (NDH) complex involved in photosystem I cyclic electron transport and chlororespiration [40]. This indicates that the American Journal of Plant Sciences photosystem I and II activity remains unaltered during water stress.

Conclusion
Whole genome bisulfite sequencing of irrigated and water stressed Zea mays Z59 plants and assembly of the reads on to the reference organellar genomes showed hypomethylation in water stressed plants when compared to irrigated plants.
Annotation of the methylated genomes showed that genes related to photosystem I & II in chloroplast and nad 4 gene in mitochondria were hypo methylated in the water stressed sample; RNA-seq analysis of transcriptomic reads mapped to the same reference showed regulation of hypomethylated genes in water stressed plants.