Identification and characterization of the cytosine-5 DNA methyltransferase gene family in Salvia miltiorrhiza

SmC5-MTase gene identification

Sequences of 11 Arabidopsis AtC5-MTase proteins were downloaded from the Arabidopsis Information Resource (TAIR, http://www.arabidopsis.org). Arabidopsis AtC5-MTase proteins were blast-analyzed against the S. miltiorrhiza 99-3 whole genome sequence using tBLASTn algorithm (Altschul et al., 1997; Xu et al., 2016). S. miltiorrhiza C5-MTase gene models were predicted from retrieved genomic DNA sequences through alignment with C5-MTase genes from other plants and S. miltiorrhiza transcriptome data (http://www.ncbi.nlm.nih.gov/sra) using the BLASTx and BLASTn program, respectively (http://www.ncbi.nlm.nih.gov/blast/). Obtained raw genomic sequences, open reading frame (ORF) sequences and deduced protein sequences were listed in Data S1.

Gene structure determination, protein sequence analysis and phylogenetic tree construction

Gene structures were determined on the Gene Structure Display Server (GSDS 2.0, http://gsds.cbi.pku.edu.cn/index.php) using the coding sequences and the corresponding genomic sequences as inputs. Protein amino acid number, molecular weight (Mw) and theoretical isoelectric point (pI) were analyzed using the EXPASY PROTOPARAM tool (http://www.expasy.org/tools/protparam.html). Multiple sequence alignment was carried out using ClustalW (Larkin et al., 2007). Conserved motifs were determined using the MEME suite (Bailey et al., 2009). The neighbor-joining (NJ) tree was constructed using MEGA version 7.0 (Kumar, Stecher & Tamura, 2016). Node robustness was detected using the bootstrap method (1,000 replicates).

RNA extraction and quantitative real-time reverse transcription polymerase chain reaction (qRT-PCR)

Roots, stems, leaves and flowers of 2-year-old S. miltiorrhiza (line 99-3) plants grown in a field nursery at the Institute of Medicinal Plant Development, Chinese Academy of Medical Sciences, were collected in August and stored in liquid nitrogen until use. Total RNA was extracted from three biological replicates for each tissue using the Quick RNA isolation kit (Huayueyang Biotechnology, Beijing, China). Each biological replicate represents independent single plant. RNA integrity and quantity were analyzed using agarose gel and NanoDrop 2000C Spectrophotometer (Thermo Scientific, Waltham, MA, USA), respectively. Total RNA was reverse-transcribed using the PrimeScript™  RT reagent kit (TaKaRa, Japan). qRT-PCR primers listed in Table S1 were designed using Primer3 (http://frodo.wi.mit.edu/primer3/). qRT-PCR analysis was conducted in triplicate on the CFX96™ real-time PCR detection system (Bio-Rad, Hercules, CA, USA) for each biological replicate using SYBR premix Ex Taq™ kit (TaKaRa, Shiga, Japan) as described (Ma et al., 2012). The 2^−ΔΔCt method was used to evaluate relative transcript levels (Livak & Schmittgen, 2001). Differential transcript abundance among tissues was determined by the one-way ANOVA (analysis of variance) method using IBM SPSS 20 software (IBM Corporation, Armonk, NY, USA).

RNA-seq data and bioinformatic analysis

Transcriptome data of periderm, phloem and xylem from S. miltiorrhiza roots was downloaded from SRA database of NCBI (SRA accession number SRR1640458) (Xu et al., 2015). RNA-seq reads from S. miltiorrhiza hairy roots treated with yeast extract (YE) and methyl jasmonate (MeJA) were downloaded from SRA database under the accession number SRP111399 (Zhou et al., 2017). RNA-seq data for S. miltiorrhiza cell cultures treated with salicylic acid (SA) were downloaded from SRA database under the accession number SRX1423774 (Zhang et al., 2016). Differential transcript abundance was analyzed using TopHat2.0.12 and Cufflinks2.2.1 (Trapnell et al., 2012). Heat maps were created using R statistical software (Le Meur & Gentleman, 2012).

Genome-wide identification and sequence feature analysis of SmC5-MTase genes

Blast analysis of eleven Arabidopsis cytosine DNA methyltransferase proteins against the whole S. miltiorrhiza 99-3 genome sequence resulted in the identification of eight SmC5-MTase genes. The number of SmC5-MTase genes in S. miltiorrhiza is comparable with that in other plants, such as Arabidopsis (Ashapkin, Kutueva & Vanyushin, 2016), rice (Ahmad et al., 2014) and maize (Qian et al., 2014), which contain eleven, ten and eight C5-MTases, respectively. The presence and distribution of conserved domains and motifs was confirmed by PFAM analysis of protein sequences (Fig. 1). Although the conserved catalytic C-terminal domains with conserved motifs are ubiquitous in all identified SmC5-MTase proteins, the N-terminals have diverse combinations of conserved domains among subfamilies (Fig. 1). It allows us to easily name newly identified genes. C5-MTases with two bromo adjacent homology (BAH) domains in the N-terminal belong to the MET subfamily, whereas C5-MTases having only one BAH domain belong to the CMT subfamily. C5-MTases containing ubiquitin associated (UBA) domains in the N-terminal were classified as the DRM subfamily. In addition, C5-MTases without N-terminal conserved domain belong to the DNMT2 subfamily. Accordingly, the eight SmC5-MTase genes were named SmMET1, SmCMT1, SmCMT2a, SmCMT2b, SmCMT3, SmDRM1, SmDRM2, and SmDNMT2, respectively (Table 1).

The length of ORFs of SmC5-MTase genes varies from 1,152 bp to 5,034 bp (Table 1). The deduced proteins vary between 383 amino acids and 1,677 amino acids, respectively. The theoretical isoelectric point (pI) varies from 4.89 (SmCMT3) to 8.03 (SmCMT2a) and the molecular weight (Mw) varies from 43.7 kDa (SmDNMT2) to 189.4 kDa (SmMET1) (Table 1). These sequence features were similar to their counterparts in Arabidopsis. On the other hand, the genomic DNA length of SmC5-MTase genes was obviously greater than their Arabidopsis counterparts (Table 1). Gene structure analysis showed that introns were responsible for the enlargement of SmC5-MTase genes (Fig. 2). The underlying mechanism of intron size expansion is unknown.

Sequence alignment and conserved motif analysis in C5-MTases

C5-MTases contain N- and C-terminal domains, which are necessary for establishment and maintenance of cytosine DNA methylation at distinct sequence contexts. In order to understand the conservation and divergence of these conserved domains in SmC5-MTases, multiple sequence alignment was performed for C5-MTases from S. miltiorrhiza and Arabidopsis. Four highly conserved motifs, including I, IV, VI and X (Malone, Blumenthal & Cheng, 1995), were found in the C-terminal methyltransferase domain of all C5-MTases analyzed (Figs. S1–S5). Among them, motifs I and X are involved in S-adenosyl-L-methionine (AdoMet) binding, motif IV has the proline-cysteine doublet acting as the functional catalytic site, and motif VI has a glutamic acid involved in target cytosine binding (Malone, Blumenthal & Cheng, 1995). In addition to the highly conserved motifs, various less conserved motifs were found in the C-terminal domain (Figs. S1–S5). These motifs could be associated with gene-specific functions. Compared with the C-terminal domain, the N-terminal domain is more diverse. The N-terminal of SmMET1 contains two replication foci domains (RFD) and two BAH domains (Fig. S1). Except for SmCMT2a, which lacks the Chr domain, the N-terminal of SmCMTs contains a BAH domain and a chromo (Chr) domain (Figs. S2–S3). The N-terminal of SmDRMs has two conserved UBA domains (Fig. S4). The divergence of the N-terminal domains could be important for distinct roles of different SmC5-MTases.

To further explore sequence conservation and divergence of SmC5-MTases, MEME motif search tool was used to identify conserved motifs in SmC5-MTases and AtC5-MTases. A total of 15 conserved motifs were identified. The length of motifs ranges from 28 to 100 amino acids (Table 2). The number of motifs in each C5-MTase varies between 2 and 10. Motifs 1, 4, 5, 6, 7, 8, 12 and 13 are located in the C-terminal DNA methyltransferase domain (Fig. 3). Among them, motifs 1, 5, 7 and 8 are highly conserved in 12–14 C5-MTases, whereas motifs 4, 6, 12 and 13 are conserved in less than ten C5-MTases. The other seven, including motifs 2, 3, 9, 10, 11, 14 and 15, are located in the N-terminal regions of C5-MTases. Among them, motifs 2, 3, 9, 14 and 15 only exist in the MET subfamily. Motifs 2 and 14 correspond to two BAH domains, whereas motifs 3 and 9 are located in the two RFD domains. Motifs 10 and 11 only exist in the CMT and the DRM subfamily, respectively. Motif 10 is specific to the BAH domain of the CMT subfamily, however motif 11 is not the UBA domain of the DRM subfamily (Fig. 3). Motifs commonly existing in C5-MTases are probably associated with the conserved functions of C5-MTases, whereas those specific to a few C5-MTases seem to be related to gene-specific functions.

Phylogenetic analysis of C5-MTases in S. miltiorrhiza and other plant species

To understand the phylogenetic relationships and evolutionary history of plant C5-MTases, an unrooted neighbor-joining phylogenetic tree was constructed for full-length protein sequences of C5-MTases from S. miltiorrhiza and other 10 plant species (Table  S2). Four subfamilies, including MET, DNMT2, CMT and DRM, are exhibited in the tree (Fig. 4). It is consistent with previous studies on other plants (Cao et al., 2014; Qian et al., 2014), confirming the domain based classification and nomenclature. The MET subfamily contains SmMET1 and can be further divided into dicot and monocot groups. The CMT subfamily can be grouped into three clades, including CMT1, CMT2 and CMT3. Among them, CMT1 and CMT3 can be further divided into dicot and monocot groups, whereas CMT2 exists in dicots only. The DRM subfamily can be divided into two dicot and two monocot groups. Similarly, The DNMT2 subfamily can be divided into dicot and monocot groups. The results imply that each subfamily is phylogenetically monophyletic and includes dicot and monocot evolutionary lineages. Taken together, each of the MET and the DNMT2 subfamilies contains a SmC5-MTase, whereas the CMT and the DRM subfamilies include four and two SmC5-MTases, respectively. SmC5-MTases have close phylogenetic relationships with the counterparts from Erythranthe guttatus (Fig. 4). It is consistent with the fact that both S. miltiorrhiza and E. guttatus are members of the Lamiales order.

Transcript abundance analysis of SmC5-MTase genes in S. miltiorrhiza

C5-MTases are important for plant growth and development (Cao & Jacobsen, 2002a; Finnegan, Peacock & Dennis, 1996; Hu et al., 2014). To preliminarily understand the function of SmC5-MTases, transcript abundance of genes in roots, stems, leaves and flowers of 2-year-old, field nursery-grown S. miltiorrhiza plants were analyzed. Transcripts of all eight SmC5-MTase genes could be detected in organs analyzed, although differential transcript levels were observed. SmMET1, which may be responsible for CG methylation maintenance, shows relatively uniform transcript levels in roots, stems, leaves and flowers of S. miltiorrhiza (Fig. 5). It is consistent with the transcript abundance patterns of AtMET1 in vegetative organs (Ashapkin, Kutueva & Vanyushin, 2016). Among the four CMT subfamily members, SmCMT1 and SmCMT3 exhibited similar transcript abundance patterns, showing more stem-specific pattern, whereas no significant difference was observed for the transcript levels of SmCMT2a and SmCMT2b in the organs analyzed (Fig. 5). The transcript abundance patterns of SmCMTs are consistent with their phylogenetic relationships (Fig. 4). The transcript level of SmDRM1 was significantly higher in flower and stem than root and leaf (Fig. 5). SmDRM2 was higher in stem, followed by flower, and less in root and leaf (Fig. 5). The transcript abundance patterns of SmDRMs are consistent with the role of DRM subfamily members in de novo methylation (Cao & Jacobsen, 2002a; Cao & Jacobsen, 2002b; Law & Jacobsen, 2010). SmDNMT2 exhibited similar transcript levels in all tissues analyzed, suggesting its importance in plant development.

S. miltiorrhiza roots are major medicinal materials of various TCMs. To further investigate possible roles of SmC5-MTases in S. miltiorrhiza, RNA-seq data from periderm, phloem and xylem of roots were analyzed. Differential transcript abundance was observed for SmCMT2a and SmDRM1 (Fig. 6, Table S3), indicating the role of SmCMT2a and SmDRM1 in the biosynthesis of bioactive compounds.

Responses of SmC5-MTases to yeast extract and methyl jasmonate treatments

Yeast extract (YE) and methyl jasmonate (MeJA) are effective elicitors for tanshinone and phenolic acid production in S. miltiorrhiza (Chen et al., 2010; Kai et al., 2012; Zhou et al., 2017). The responses of SmC5-MTases to YE and MeJA treatments were analyzed using RNA-seq data from S. miltiorrhiza hairy roots (Zhou et al., 2017). The transcripts of all SmC5-MTase genes were detected in hairy roots (Fig. 7, Table S3). Compared with non-treated control, five SmC5-MTases, including SmMET1, SmCMT1, SmCMT2a, SmCMT2b and SmCMT3, were slightly down-regulated at the time-point of 12 h of YE treatment and at the time-point of 6 h of MeJA treatment. SmDRM1 and SmDRM2 were slightly down-regulated at the time-point of 1 h after YE and MeJA treatments, whereas they were up-regulated with treatment time extension. SmDNMT2 was slightly up-regulated at the time-point of 1 h after two treatments, whereas it was down-regulated after 6 h MeJA treatment. Although all SmC5-MTases exhibited differential transcript abundance compared with control, the variance was not significant. It could be due to the short treatment time. Further analysis of SmC5-MTase genes in response to treatment with extended treatment time will help to elucidate the role of SmC5-MTases in secondary metabolite biosynthesis.

Salicylic acid-responsive SmC5-MTases

DNA methylation plays a vital role in plant defense by regulating the expression of a subset of stress responsive genes and is involved in priming of salicylic acid (SA)-dependent immunity (López Sánchez et al., 2016; Yu et al., 2013). In order to investigate the potential role of SmC5-MTases in plant defense, RNA-seq data from SA-treated S. miltiorrhiza cell cultures was analyzed (Zhang et al., 2016). The transcripts of all SmC5-MTase genes were determined in S. miltiorrhiza cell cultures and the majority of them were significantly suppressed at 8 h after SA treatment (Fig. 8, Table S3), suggesting the significance of SmC5-MTases in SA-dependent defenses.