Genome-wide identification of mitogen-activated protein kinase (MAPK) cascade and expression profiling of CmMAPKs in melon (Cucumis melo L.).

Mitogen-activated protein kinase (MAPK) is a form of serine/threonine protein kinase that activated by extracellular stimulation acting through the MAPK cascade (MAPKKK-MAPKK-MAPK). The MAPK cascade gene family, an important family of protein kinases, plays a vital role in responding to various stresses and hormone signal transduction processes in plants. In this study, we identified 14 CmMAPKs, 6 CmMAPKKs and 64 CmMAPKKKs in melon genome. Based on structural characteristics and a comparison of phylogenetic relationships of MAPK gene families from Arabidopsis, cucumber and watermelon, CmMAPKs and CmMAPKKs were categorized into 4 groups, and CmMAPKKKs were categorized into 3 groups. Furthermore, chromosome location revealed an unevenly distribution on chromosomes of MAPK cascade genes in melon, respectively. Eventually, qRT-PCR analysis showed that all 14 CmMAPKs had different expression patterns under drought, salt, salicylic acid (SA), methyl jasmonate (MeJA), red light (RL), and Podosphaera xanthii (P. xanthii) treatments. Overall, the expression levels of CmMAPK3 and CmMAPK7 under different treatments were higher than those in control. Our study provides an important basis for future functional verification of MAPK genes in regulating responses to stress and signal substance in melon.


Introduction
Plants often suffer from various stresses including abiotic stresses (low temperature, drought, salt, etc.) and biotic stresses (insects, weeds, diseases, etc.) which attenuate plant growth, development, and adversely affect the quality and productivity. It is necessary to consider the damage and response mechanisms of plants under stresses for improving resistance [1]. Signal

Identification of MAPK cascade gene family members in melon
The method used to identify all the putative MAPK, MAPKK and MAPKKK family genes was similar to that used for cucumber and watermelon [36,37]. Firstly, reported amino acid sequences of MAPK, MAPKK and MAPKKK family members in Arabidopsis, cucumber and watermelon were downloaded from the database (Arabidopsis http://www.arabidopsis.org/, Cucumber http://cucurbitgenomics.org/organism/20, Watermelon http://cucurbitgenomics. org/organism/21). Then, these amino acid sequences were utilized to search against the melon database (http://melonomics.net/) using the BlastP program with an e-value of 1e -10 as the threshold to get candidate MAPK cascade gene family members. A further search was performed with Hidden Markov Model (HMM) analysis using HMMER 3.0 program. The HMMER hits containing the serine/threonine-protein kinase-like domain (PF00069) were compared with blast results and parsed by manual editing. All data were checked for redundancy by self-BLAST and no any alternative splice variants were considered. Finally, the online software Conserved Domain Database of NCBI (http://blast.ncbi.nlm.nih.gov) and SMART database (http://smart.embl-heidelberg.de/) were used to further confirm the predicted MAPK cascade genes.

Characteristics, conserved domain detection and subcellular location predication
ExPASy database (http://web.expasy.org/compute_pi/) was served to predict the molecular weight (MW) and isoelectric points (pI) of the melon MAPK, MAPKK and MAPKKK family proteins. MEME software (http://meme-suite.org/tools/meme) was used to analyze the conserved protein motif. The subcellular localization of MAPK cascade gene family members was analyzed using the TargetP software of the CELLO v2.5 server (http://cello.life.nctu.edu.tw/) [40].

Multiple sequence alignments and phylogenetic construction
Multiple sequence alignments were generated using DNAMAN 9.0 software. Phylogenetic analysis was carried based on the full-length protein sequences using the MEGA 6 program by the neighbor-joining (NJ) method and a bootstrap test was carried out with 1000 interactions [41]. Built on the nomenclature, the MAPK, MAPKK and MAPKKK genes in melon were named according to the homology of the corresponding genes in Arabidopsis, cucumber and watermelon.

Chromosomal location and gene structure analysis
The nucleotide sequences and chromosomal locations of all these genes were obtained from the Cucurbit Genomics Database (http://cucurbitgenomics.org/) and mapped to the chromosomes with MapInspect software, as well as the gene structures of the melon MAPK, MAPKK and MAPKKK family proteins were generated with the GSDS 2.0 (http://gsds.cbi.pku.edu.cn/).
For drought and salt stress treatments, seedlings at two-leaf stage grown in pots were transferred to Yamasaki melon nutrient solution and refreshed every 2 days. Then 8% (w/v) PEG-6000 and 100 mM NaCl were added to simulate drought and salt stress when seedlings were four fully expanded leaves, respectively. For SA and MeJA treatments, seedlings at four-leaf stage were sprayed with 5 mM SA and 100 μM MeJA (10 mL per plant) or with 0.1% ethanol as a control, respectively. Leaves were collected at 0 h, 1 h, 3 h, 24 h, 120 h after treatments for qRT-PCR analysis. For RL treatment, leaves were sampled after exposed to RL (160 μmol m -2 s -1 ) and white light (WL, 160 μmol m -2 s -1 ) for 6h. For biotic stress, P. xanthii were collected from the infected plants in the field and evenly distributed on the third and fourth leaves of the melon seedlings. In the meantime, RL were used as the inducer of the defense response in the P. xanthii infection experiment and untreated seedlings were taken as control. Leaves were collected after treatments at 0 d, 1 d, 3 d, 5 d, 7 d and 9 d for qRT-PCR analysis.

RNA isolation, cDNA synthesis and qRT-PCR
The ultrapure RNA kit (Cat.CW0581, Beijing, China) was used to extract total RNA from leaves according to the manufacturer's instructions. cDNA was synthesized using PrimeScript™ RT Master Mix (Perfect Real Time) (TAKARA, Japan) and the relative expression of the CmMAPK genes was quantified by qRT-PCR according to the SYBR 1 Premix Ex Taq TM (Tli RNaseH Plus) (TAKARA, Japan) manufacturer's protocol on an ABI PRISM 7500 real-time fluorescent quantitative PCR (Applied Biosystems, Thermo Fisher Scientific, USA). The primers for the qRT-PCR analysis were designed to avoid the conserved regions and listed in the S2 Table [43]. The 2 -ΔΔCt method was used to measure relative expressions of genes, and 18s rRNA from melon was used as the internal control gene. All experiments were performed with three biological and three technical replicates.

Statistical analysis
The data were analyzed using the SPSS 22.0 statistics program, and significant differences were compared using a one-way ANOVA following Duncan's multiple range test at p<0.05 level except expression profiles of CmMAPKs under RL induction were conducted by t test at p<0.05 level.

Identification of MAPK cascade genes in melon
To identify MAPK, MAPKK and MAPKKK family genes, we conducted respective BlastP searches against the melon protein databases using query protein sequences from Arabidopsis, cucumber and watermelon. Meanwhile, HMM search was used to identify the serine/threonine-protein kinase-like domains (PF00069) in all potential MAPK cascade sequences from melon. After removing redundancies, alternative splices, multiple steps of screening and validation of the conserved domains, we finally identified 14 CmMAPK genes, 6 CmMAPKK genes, and 64 CmMAPKKK genes. The sequence data of all above MAPK cascade genes was downloaded from the Cucurbit Genomics Database (http://cucurbitgenomics.org/ search/ genome/18). Each gene was named according to its homology with Arabidopsis MAPK, MAPKK, or MAPKKK proteins [44].
As shown in Table 1, 14 CmMAPKs were randomly distributed on 8 chromosomes. These CmMAPKs were predicted to encoding 370 (CmMAPK3, CmMAPK4-2, CmMAPK13) to 645 (CmMAPK9-1) amino acids (aa) in length, with putative MW ranging from 42.6 (CmMAPK13) to 73.1 kDa (CmMAPK9-1) and theoretical pIs ranging from 5.01 (CmMAPK13) to 9.28 (CmMAPK20-2). The subcellular localization was predicated and found that CmMAPKs were located in nuclear and cytoplasm, except CmMAPK7 and CmMAPK20-2 may be present in mitochondria. Meanwhile, 6 CmMAPKKs were predicted to encoding 340-518 aa, with MW ranging from 38-57.8 kDa and pIs ranging from 5.11-9.04. They were predicted to be located in cytoplasm and nuclear. All 64 CmMAPKKKs were randomly distributed on 13 chromosomes. The polypeptide lengths of these MAPKKKs ranged from 221 to 1208 aa with MW ranging from 24.6 to 133.5 kDa and pIs ranging from 4.63 to 9.41. Besides, 47 of the CmMAPKKKs were located in nuclear and cytoplasm, while the others may be present in plasma membrane, mitochondria or chloroplast (S1 Table).

Phylogenetic relationship and structures of MAPK cascade genes
To further characterize the MAPK cascade genes, unrooted phylogenetic trees were produced by aligning the full-length protein sequences of all MAPK cascade genes in Arabidopsis, watermelon, cucumber and melon using NJ method.

PLOS ONE
C and D. Except group A contained three members (CmMAPKK2-1, CmMAPKK2-2 and CmMAPKK6), each of the other groups contained one member. Notably, CmMAPKK5 and CmMAPKK9 were closely related in phylogenetic development, but belonged to different groups (Fig 1). Analysis of the gene structures showed that all three members of group A contained 8 introns, CmMAPKK3 contained 9 introns in group B, and CmMAPKK5 and CmMAPKK9 contained 1 intron in group C and D respectively (Fig 2). The MAPKKK gene family in plants is usually classified into 3 categories, namely MEKK, Raf, and ZIK [45][46]. Refer to the classification method of Arabidopsis, MEKK group contained 20 members, Raf 34 members and ZIK 10 members (S1 Fig). Genetic structures analysis showed that the number of introns of MAPKKK gene family members in melon varied greatly from 2 to 36. MAPKKK members in the same group also had different genetic structures. In MEKK group, 9 genes contained 2 introns, 3 genes contained 4 introns, and the remaining 8 genes contained 16-36 introns. In ZIK group, 2 genes contained 4 introns, 1 gene contained 16 introns, and the remaining 7 genes contained 14 introns. RAF group members had the largest differences in gene structure among 3 groups, with the number of introns ranging from 4 to 32. Through careful analysis, it was found that MAPKKK genes had similar or identical gene structure, these genes owned high homology and close phylogenetic relationship.

Conserved domains and motifs in MAPK cascade proteins
To survey the structural domain features of MAPK cascade gene families, a multiple sequence alignment of the predicted amino acid residues was carried out using the DNAMAN 9.0 software. The results showed that all CmMAPK proteins had conserved motif TxY and CD domain (LH) DxxDE (P) xC existing in group A and B (Fig 3). The conserved motifs of the CmMAPK proteins were analyzed using the online tool MEME. Results displayed that all CmMAPK proteins contained motif 2-8, motif 1 only existed in group A and B, and motif 9-10 only existed in group C and D (Fig 4). Multiple alignment of the CmMAPKK proteins showed that all proteins contained conserved D (L/I/V) K, the phosphorylation target site S/Txxxxx-S/T and signature VGTxxYMSPER, which were conserved in the catalytic domain ( Fig  3) [30]. Motif analysis found that CmMAPKK6, CmMAPKK5, and CmMAPKK9 contained motif 1-8 and the rest had motif deletion (Fig 4). These motifs may lead to different functions of each gene, which further upheld the accuracy of classification.
Multiple sequence alignment of CmMAPKKKs from melon showed that melon and other objects were similar to MEKK group containing conserved G Through the analysis of protein conserved motifs and protein functional domains, it was found that members of the MAPK cascade pathway gene family in melon were very conserved in structure, therefore we can predict their functions in melon based on the reported homologous gene functions in Arabidopsis [44], cucumber [36] and watermelon [37].

Chromosomal localization of MAPK cascade genes
BLASTN searches were performed against the melon genome database using the DNA sequence of each MAPK cascade gene to determine the chromosomal distribution and

PLOS ONE
transcriptional direction of the MAPK cascade genes in melon. By chromosome localization analysis, it was found that 84 genes in MAPK cascade pathway were distributed on 13 chromosomes (Fig 5). It was notable that some genes showed high homology, but they were distributed on different chromosomes such as CmMAPK4-1/CmMAPK4-2, CmMAPK9-1/ CmMAPK9-2/CmMAPK9-4, CmMAPK20-1/CmMAPK20-2, CmMKK2-1/CmMKK2-2 and CmZIK4-1/CmZIK4-2/CmZIK4-3. From the above results, we could infer that there were a large number of gene duplications and gene transfer events in the long-term evolution of members of the MAPK cascade pathway gene family in melon. They participated in both chromosome fragment replication events and tandem replication events, which were of profound significance for the expansion of MAPK cascade pathway gene family members in melon.

Expression profiles of CmMAPKs in melon seedlings under drought and salt stress
In open field cultivated of melon plants, it is often affected by drought and soil salinization, which further affects the yield and quality of melon. In order to study CmMAPKs expression in melon under drought and salt stress, we used nutrient solution cultivation way, by adding 8% PEG into the nutrient solution (simulation of drought stress) and 100 mM NaCl (simulation of salt stress) to treat melon seedlings. It was found that the expression pattern of CmMAPK genes in melon changed, and the expression trend of multiple genes was consistent under the same treatment (Fig 6). Under drought stress, except that the transcription levels of CmMAPK19, CmMAPK20-1 and CmMAPK20-2 decreased, the other members increased at 6 h, 12 h and 24 h after treatment, and then decreased to the original level. These results

PLOS ONE
indicated that these up-regulated genes might be critical genes respond to the drought stress. Under salt stress, only CmMAPK3, CmMAPK7, CmMAPK20-1 and CmMAPK20-2

Expression profiles of CmMAPKs in response to SA and MeJA treatment
In plant cells, hormones can affect physiological and biochemical reactions of plants through a variety of signal transduction pathways. SA and JA are important hormones in plant immune system. To study whether CmMAPKs in melon were regulated by hormone signaling substances, exogenous spraying of 5 mM SA and 100 μM MeJA on leaves of melon seedlings was carried out in this study. Cluster analysis was performed on the data results. After SA induction, except the down-regulated expression of CmMAPK20-1 and CmMAPK20-2, other genes induced to up-regulate expression, indicating that different members of CmMAPKs played different roles in SA-induced defense response. At the later stage of MeJA induction, all 14 CmMAPKs strongly induced, revealing that members of CmMAPKs positively responded to all the induced defense responses of MeJA and played an essential role in this process (Fig 7). In summary, members of CmMAPKs may be widely involved in signaling hormone-induced defense responses and play different roles.

Expression profiles of CmMAPKs in melon seedlings under RL induction
Studies have shown that moderate red light induction can improve the resistance of plants to Pst-DC3000 [47]. Our preliminary study of the research group showed that the disease resistance significantly increased by pre-irradiation of red light (160 μmol m -2 s -1 ) for 6h at melon seedling stage, when melon seedlings were inoculated with P. xanthii. In order to explore whether CmMAPKs were involved in the defense reaction induced by red light, this experiment conducted red light irradiation on melon seedlings and white light irradiation as the control. The results showed after red light induction, the transcription levels of CmMAPK3, CmMAPK6-1 and CmMAPK7 in melon leaves remarkably increased compared to the control. It was worth noting that the transcription levels of CmMAPK3 were up to 5 folds higher than that of the control, suggesting that they might play a positive regulatory role in the red light induced defense response. However, the transcription levels of CmMAPK20-1 and CmMAPK20-2 were not significantly different indicating that they were not sensitive to red light induction and played an auxiliary role in the red light induced defense response. The transcription levels of other gene family members significantly down-regulated after red light induction, which may act as a negative regulatory role in the red light induced defense response (Fig 8).

Expression profiles of CmMAPKs in response to pathogen infection
Gene expression patterns are typically used as an indicator of gene function. In the protected cultivation of melon, high temperature and high humidity environment make it vulnerable to the pathogen of P. xanthii, which restricts the growth of melon and seriously affects the yield and quality. Therefore, we selected P. xanthii as the biotic stress factor to infect the leaves of melon in this study. Cluster analysis was performed on the data results (Fig 9). In addition to CmMAPK 1, CmMAPK6-1, CmMAPK7 and CmMAPK 19, other CmMAPKs greatly up-regulated on the third day after infection by P. xanthii. Moreover, the expression level of these CmMAPKs gradually decreased with the increase of the invasion time of P. xanthii, except CmMAPK3, CmMAPK4-2 and CmMAPK13. At the same time, the transcription levels of most CmMAPKs with pre-irradiation of RL were significantly lower than those of P. xanthii inoculated alone at 3 and 5 days, and were significantly upregulated at 7 and 9 days after inoculation. It is speculated that CmMAPKs may play a negative regulatory role in the mechanism of disease resistance after irradiation of red light.

Discussion
During plant growth and development, it is inevitable to be subjected to stresses caused by external environments affecting its growth and development, while plants will form a series of defense reaction mechanisms. Reversible phosphorylation of proteins plays an important part in the recognition, gradual transduction and amplification of extracellular signals. Mitogenactivated protein kinase (MAPK), a type of serine/threonine protein kinase activated by the MAPK cascade (MAPKKK-MAPKK-MAPK), transmits the upstream signal cascade to the downstream response molecule via protein phosphorylation, thereby activating the expression of the stress-resistant genes, making the plant adapt to the adversity [48]. However, little is known about the identity and functions of melon MAPK cascade genes. Thus, we mapped out the major structural characteristics of all MAPK cascade genes and expression profiles of CmMAPKs in melon.
In the present study, a comprehensive set of 84 MAPK cascade genes were identified including 14 CmMAPKs, 6 CmMAPKKs and 64 CmMAPKKKs through BlastP homologous alignment combined with HMMER searches. Analyzing the amino acid sequences of MAPK gene members from Arabidopsis and other species, MAPK family members are divided into 4 groups: A, B, C, and D. Among them, MAPK proteins in A, B, and C contain TEY phosphorylation sites, while proteins in D contains TDY phosphorylation sites [44]. In our investigation of MAPK cascade, phylogenetic trees suggested that melon CmMAPK and CmMAPKK were mainly categorized into 4 groups, and CmMAPKKK were categorized into 3 groups, revealing the evolutionary homology between melon and Arabidopsis. Moreover, we found that the melon MAPK cascade genes are structurally conservative by analyzing the protein conserved motifs and protein functional domains. For example, CmMAPK and CmMAPKK proteins contained the conserved motif TxY and S/T-xxxxx-S/T, respectively. Previous research has suggested that gene duplication events play an important role in the expansion of gene families [49]. In our study, chromosome localization showed that 84 members of the MAPK cascade gene family of melon were distributed on 13 chromosomes indicating that these genes were involved in chromosome fragment replication events and tandem replication, which are indispensable for the expansion of family members. These features in melon have also been observed in Arabidopsis [44], cucumber [36] and watermelon [37], which supports our identification of MAPK cascade genes.
Stress resistance and disease resistance are the result of mutual adaptation and selection between plants and environmental stresses and pathogens. The expression of gene changes with environmental stress and pathogen infection. It is previously shown that Arabidopsis MAPK1/2/7/14 are induced by ABA [30,31]. Some of watermelon ClMAPKs (ClMAPK1, ClMAPK3, ClMAPK4-1, ClMAPK4-2, ClMAPK6, ClMAPK7, ClMAPK9-1, ClMAPK9-2, ClMAPK13, ClMAPK16, ClMAPK19, ClMAPK20-1, ClMAPK20-2) can respond to drought stress [37]. However, there is no information available about the response to CmMAPKs in stress-responsiveness. In our study, the transcripts of the 14 genes were firstly detected in leaves under drought, salt, SA and MeJA, RL and P. xanthii treatments using qRT-PCR. As shown from the heat map in Fig 6, some genes like CmMAPK19, CmMAPK20-1 and CmMAPK20-2 were down-regulated, and the remaining CmMAPKs decreased firstly after treatment and then returned to normal levels under drought stress. Besides, CmMAPK3, CmMAPK7, CmMAPK20-1 and CmMAPK20-2 significantly up-regulated under salt stress. After SA treatment, the expression levels of CmMAPK20-1 and CmMAPK20-2 were notably lower than that of the control, and other members were up-regulated. Nevertheless, CmMAPKs were strongly induced to express after spraying with MeJA for 120 h (Fig 7). After treatment with RL, CmMAPK3, CmMAPK6-1 and CmMAPK7 were up-regulated, even the expression of CmMAPK3 was up to 5-fold compared with control, but the expression levels of CmMAPK20-1 and CmMAPK20-2 changed little and other members were down-regulated (Fig 8). In addition, CmMAPK1, CmMAPK6-1 and CmMAPK19 were down-regulated, and the remaining members were expressed in different degrees after inoculation with P. xanthii. Interestingly, after pretreatment with RL, the expression of CmMAPKs was lower than that of P. xanthii inoculated alone, while some CmMAPKs were up-regulated, such as CmMAPK3, CmMAPK4, CmMAPK7 and CmMAPK9-2, suggesting that CmMAPKs may play a negative regulatory role in the red light-induced resistance to P. xanthii (Fig 9). The results suggested that CmMAPK3 and CmMAPK7 were up-regulated remarkably after different treatments indicating that these two genes may participate in multiple signal transduction pathways at the same time, playing important roles in the melon resistance network. It is found that Arabidopsis MAPK3 belongs to group A. The protoplasts are treated with bacterial flagellin (flg22) to detect MAPK activity: MAPK3 is activated, and the downstream transcription factor WRKY22/29/33 are activated to increase plant resistance [12,28,29]. Previous studies have proven that AtMAPK3 is involved in the synthesis of ethylene and phytoalexin, the biosynthesis of indole glucosinolates and stomatal immunity [32][33][34][35]. Besides, MAPK3/6 negatively regulates the tolerance of Arabidopsis to low temperature by phosphorylating ICE1 [50]. Further, on the 9th day after infection of Fusarium oxysporum, the expression of ClMAPK3 is significantly higher than that of the control [37]. Thus, MAPK3 is widely involved in resistance to stresses in plants. The MAPK members of different group are involved in many signal transductions while there are less reports and gene function studies of group C [51,52]. We found that CmMAPK1 and CmMAPK7 identified in melon belonged to group C of the MAPK gene family. In Arabidopsis, AtMAPK1, AtMAPK2, AtMAPK7, and AtMAPK14 are classified as group C, of which AtMAPK7 may participate in pathogen signaling pathways that regulate PR gene expression [53]. GhMAPK7 has been documented to take part in plant disease resistance, which expression can be activated by pathogen infection and multiple defense response signaling molecules [54]. In addition, ZmMAPK7 can improve the tolerance to oxidative stress and found that it mediates ABA and H 2 O 2 signaling processes in transgenic tobacco [55]. The expression patterns of watermelon MAPK family members in abiotic stress, hormone induction and Fusarium oxysporum infection show that the overlapping expression characteristics of ClMAPK7 gene to various stress factors predicted the pleiotropic effect of ClMAPK7 [37]. Our findings suggested that CmMAPK3 and CmMAPK7 homologous to AtMAPK3 and AtMAPK7 were involved in various abiotic stresses, hormone signal induction, RL induction and P. xanthii infection. However, whether the function is consistent with that of the identified and reported genes in each of the above plants should be verified by transient expression and gene editing techniques. Besides, in-depth research of CmMAPK7 can complement and improve the function of group C.

Conclusions
In our study, 14 CmMAPKs, 6 CmMAPKKs and 64 CmMAPKKKs were identified in the genome of melon. CmMAPKs and CmMAPKKs were grouped into 4 groups (A, B, C and D) while CmMAPKKKs were grouped into 3 groups (MEKK, ZIK and RAF). Besides, the study provided the first comprehensive evaluation of the phylogenetic relationships, gene structures, conserved domains and motifs, and chromosomal localizations of MAPK cascade genes. Furthermore, we analyzed expression profiles of CmMAPKs under drought, salt, SA, MeJA, RL and P. xanthii treatments. Most of the CmMAPK genes were induced by biotic and abiotic stress, as well as signal substance (RL, SA and MeJA). Our results also suggested the expression levels of CmMAPK3 and CmMAPK7 in different treatments were higher than those in control, suggesting that they were widely involved in various stress resistance mechanisms, which could be used as future research priorities.