Rice seeds (Oryza sativa L. Dongjin) were sterilized in 0.05% Spotak solution (Bayer crop science, South Korea) overnight at 4°C, and washed with distilled water 5 times. Sterilized seeds were then placed onto moist tissue paper and incubated at 28°C in the dark for germination. Germinated seeds were transferred to sterilized soil in a growth chamber (70% humidity, 25°C, and a light and dark cycle of 14 and 10 h) for four weeks.

For DNA microarray analysis, four-week-old rice leaves were sprayed with 5 ml of solution containing PBS (control) or 1 μM of purified recombinant MSP1-His protein dissolved in PBS, using a foliar spray. Recombinant MSP1-His protein was expressed in E. coli, and purified as previously described (Wang et al., 2016). Each treatment was performed using six plants in two replicates (6×2 rice plants/treatment), resulting in 12 plants per treatment. Treated leaves were harvested and frozen in liquid nitrogen at 6 h post-treatment (hpt) and 24 hpt, then stored at −70°C until further analysis. A total of six RNA samples were subjected to DNA microarray analysis.

Total RNA was isolated using an RNeasy Plant Mini Kit (QIAGEN, Hilden, Germany) according to the manufacturer’s instructions. Genomics DNA was removed with DNase I (QIAGEN, Hilden, Germany). Synthesis of target cRNA was performed using Agilent’s Low RNA Input Linear Amplification Kit PLUS (Agilent Technology, CA, USA) according to the manufacturer’s instructions. Briefly, each 0.2 μg of total RNA was mixed with the diluted Spike mix and T7 promoter primer mix, then incubated at 65°C for 10 min. Next, cDNA master mix (5× First strand buffer, 0.1 M DTT, 10 mM dNTP mix, RNase-Out, and MMLV-RT) was prepared and added to the reaction mixer. The procedures for DNA microarray were identical to those described in our previous study (Wang et al., 2014). All data normalization and selection of fold-changed genes was performed using GeneSpring GX 7.3 (Agilent Technology). Genes were filtered by removing flag-out genes in each experiment. The averages of normalized ratios were calculated by dividing the test channel intensity by the control channel intensity.

The MapMan program, version 3.6.ORC1 (Thimm et al., 2004), was used for pathway analysis. Fold-change values were log₂ transformed, after which their means were calculated. These genes were classified into MapMan BINs and their annotated functions were visualized using the MapMan program by searching against Oryzae sativa Osa_MSU_v7 mapping. The Page-Man visualization tool was used for GSEA with the Wilcoxon test (Benjamini_Hochberg correction and 1.0 as the ORA cutoff).

The pathway enrichment study was performed using the KEGG database (Kanehisa et al., 2014) with the DAVID Functional Annotation Tool against the background of Oryza sativa (Huang et al., 2008). Genes up-regulated or down-regulated in MSP1 treated samples at each sampling time point were loaded as gene lists in DAVID. A corrected p-value < 0.05 denoted the significance of enriched pathways. The top 10 significantly enriched pathways of each list were selected for comparison.

Two ug (this should be greek symbol micro) of total RNA from each treated sample (the same samples as for the microarray) was used to synthesize cDNA with RNA to cDNA EcoDry Premix kit (Clontech, Mountain View, CA, USA). qPCR was performed using a Rotor-Gene Q instrument (QIAGEN, Hilden, Germany) with QuantiNova SYBR Green RT-PCR kit (QIAGEN, Hilden, Germany). The following PCR reaction program was used: 95°C for 10 min followed by 40 cycles of 95°C for 10 s and 60°C for 20 s. Fold change was calculated relative to plants treated with PBS. The transcript level represented the mean of three replicates after normalization using Ubiquitin 5 (LOC_Os01g22490.1), which varied little in our microarray analysis and was previously used as a reference gene (Jain et al., 2006). Primers used are listed in Supplementary Table 1.

DNA microarray analysis was used to investigate differentially expressed genes (DEGs) (p value < 0.01 and log₂ fold change > 1 or < −1) in response to MSP1 treatment in rice (Fig. 1, Supplementary Table 2). A total of 2268 DEGs were identified at both time points, whereas 475 and 1643 DEGs were up-regulated only at 6 hpt and 24 hpt, respectively (Fig. 1A), indicating a high similarity of DEGs at different sampling time points. These results indicate that MSP1 induces a dramatic change in rice at the transcriptional level. Among all DEGs, 1188 were up-regulated and 1555 were down-regulated at 6 hpt, while 1961 genes were up-regulated and 1950 genes were down-regulated at 24 hpt (Fig. 1B). Interestingly, almost 60% of DEGs showed less than 4-fold changes (Log₂ fold change < 2) at 6 hpt and the percentage of DEGs with higher fold changes (Log₂ fold change > 2) were increased in 24 hpt relative to the 6 hpt samples (Fig. 1C), indicating that at 24 hpt rice exhibits a more intense response to MSP1 treatments. Among all up-regulated genes, 1050 were present in both samples, with 138 and 911 genes specific for 6 hpt and 24 hpt, respectively (Fig. 1D). Among down-regulated genes, 1138 were present in both samples, with 417 and 812 genes specific for 6 hpt samples and 24 hpt samples, respectively (Fig. 1D). These results indicate that the transcriptomic changes in response to MSP1 occurred more dramatically at 24 hpt.

The Page-Man web tool of MapMan was used for functional annotation of the DEGs. Both 6 hpt and 24 hpt samples exhibited an elevation in glutathione S transferases, abiotic stress, regulation of transcription and protein degradation, while down-regulated genes were related to photosynthesis, protein synthesis and transport (Fig. 1E, Fig. 2, Supplementary Table 3). At 6 hpt, genes related to lipid degradation, amino acid degradation, abscisic acid metabolism, and development were up-regulated, while genes associated with secondary metabolism, ribosomal protein in protein synthesis, aspartate protease and protein folding were down-regulated (Fig. 1E, Fig. 2, Supplementary Table 3). At 24 hpt, up-regulated genes were mainly involved in minor CHO metabolism, ethylene metabolism, WRKY domain containing transcription factor and calcium signaling, while down-regulated genes were involved in processes of redox signaling (thioredoxin) and peroxidases (Fig. 1E, Fig. 2, Supplementary Table 3). These results show the processes responsive to MSP1 induced signaling at different sampling time points.

KEGG pathway enrichment analysis using DAVID bioinformatics tool helped to determine the metabolic pathways that were commonly or specifically affected by MSP1. Alpha-linolenic acid metabolism, glutathione metabolism and branched chain amino acid degradation (valine, leucine, and isoleucine) were commonly up-regulated, while photosynthesis, photosynthesis-antenna proteins, carbon fixation in photosynthetic organisms, glyoxylate and dicarboxylate metabolism and biosynthesis of antibiotics were commonly down-regulated at both 6 hpt and 24 hpt in response to MSP1. Furthermore, at 6 hpt, protein processing in the endoplasmic reticulum, nitrogen metabolism, starch and sucrose metabolism and tryptophan metabolism were up-regulated, while the pentose phosphate pathway and ribosome were specifically inhibited (Fig. 3A, B). At 24 hpt, genes related to plant hormone signal transduction, plant-pathogen interaction, galactose metabolism, amino acid metabolism (alanine, aspartate, and glutamate) and inositol phosphate metabolism were enhanced, while porphyrin and chlorophyll metabolism and amino acid metabolism (glycine, serine, and threonine) were repressed (Fig. 3C, D).

The expression pattern and functions of differentially regulated genes were visualized using MapMan to get insights into their biotic stress-related functions. Of the total 4836 DEGs, 1106 (22.8%) were classified as biotic stress response genes (Fig. 4, Supplementary Table 4). Functional classification of biotic stress response genes indicated a general elevation of ABA metabolism, ethylene metabolism, JA metabolism, signaling, genes encoding glutathione-S-transferases, heat shock proteins, β-glucanase, transcription factors (ERF, bZIP, WRKY, MYB, DOF) and genes involved in proteolysis. In cases of hormone signaling, 14 of 15 genes (e.g., NCED2, NCED3, BZIP23, aldehyde oxidase AO and carotenoid cleavage dioxygenase 8a) involved in ABA, 30 of 46 genes (e.g., ethylene response factor ERF3, ERF57, ERF77, ERF83, ERF91, ERF95, ERF96 and ERF98, OsACS2 and OsbHLH066) involved in ethylene and 11 of 16 genes (e.g., allene oxide synthase AOS2, AOS3, lipoxygenase LOX1, LOX6, LOX8, 12-oxo-phytodienoic acid reductase OPR1, OPR5, OPR6, OPR9 and OsJAZ5) involved in JA signaling were induced by MSP1 treatment. Thirty of 43 genes encoding heat shock proteins (including HSP90, HSP82, HSP70, HSP26, HSP24.1, HSP18.0, HSP17.9, HSP17.7, HSP17.3 and HSP16.9) were up-regulated in response to MSP1 treatment. Moreover, 26 of 33 genes encoding glutathione-S-transferases (e.g., GSTF10, GSTF11, GSTF5, GSTL2, GSTL3, GSTT1, GSTU1, GSTU12, GSTU19, GSTU21, GSTU24, GSTU28, GSTU31, GSTU37, GSTU40, GSTU47, GSTU50, GSTU7, GSTU8 and GSTU9) were induced by MSP1 treatment. These results indicate that MSP1 treatments induce a dramatic defense response in rice leaves.

In plants, receptor-like kinases (RLKs) and receptor-like proteins (RLPs) are employed as PRRs to monitor the environment and detect potential dangers (Boutrot and Zipfel, 2017). The RLKs are a large protein family with over 600 members in the model plant Arabidopsis thaliana, which have 1131 members in rice (Shiu et al., 2004). Of 89 RLK genes differentially regulated in MSP1 treated samples, 31 were up-regulated in at least in 24 hpt samples. Most of these 31 RLKs were leucine-rich repeat RLKs and S-locus glycoprotein like RLKs. When compared with transcriptomics changes in M. oryzae (R line, KJ401; S line, KJ301) infected rice leaves we reported previously (Wang et al., 2014), these RLKs showed a high similarity of expression patterns in MSP1 treated samples (Fig. 5B), demonstrating their importance in rice immune responses.

MAPK (mitogen-activated protein kinase) cascades, a set of kinases for signal transduction through phosphorylation, were activated as one of the earliest signaling events after recognition of PAMPs by plant PRRs, which consist of a MAPKKK-MEK-MPK module (Meng and Zhang, 2013). To identify important MAPK signaling cascades responsive to MSP1, we focused on expression changes in all MAPKKK, MEK and MPK genes. In total, 59 MAPKKK, four MEK, and 12 MPK genes were differentially regulated in MSP1 treated samples. Of these MAPK genes, 14 MAPKKK, one MEK, and two MPK genes were induced at 24 hpt (Fig. 5A). Among the induced genes, 13 MAPKKK genes, one MKK gene and one MPK gene were also up-regulated at 6 hpt. In 24 hpt samples, the most strongly induced MAPK gene was MAPKKK62, which was also highly up-regulated in 6 hpt samples. Comparison of MAPK genes between MSP1 treated samples and M. oryzae infected samples revealed a high similarity of expression patterns.

WRKY transcription factors were reported to regulate the transcriptional reprogramming related to plant immune responses (Eulgem and Somssich, 2007). The rice genome contains a total of 113 WRKY genes, which is greater than the number of WRKY genes in Arabidopsis (Ryu et al., 2006). In total, 55 differentially regulated WRKY genes were identified in MSP1 treated samples. Of these, 30 were up-regulated in the 24 hpt samples and selected for comparison between MSP1 treated sample and M. oryzae infected samples, revealing a high similarity of expression patterns (Fig. 5A, Supplementary Table 5).

Pathogenesis-related proteins have been reported to possess antimicrobial activities via hydrolytic activities on cell walls, contact toxicity, and involvement in defense signaling (van Loon et al., 2006). In total, 19 genes encoding PR proteins were detected in response to MSP1 treatment. Among these, three (JIOSPR10, CHT1, and PR1-71) were commonly induced in samples at both time points, while five genes (PR10B, PR4C, CHT3, PR5 and PR1b) and one genes (PR1a) were specifically up-regulated in the 24 and 6 hpt samples, respectively. These results indicate that MSP1 treatment induces an intensive defense response in rice leaves at both time points.

To validate the expression patterns of the identified DEGs, eight were further selected and confirmed by qPCR in 24 hpt samples. Relative expression levels obtained by qPCR were converted to log₂ fold change (MSP1/PBS treatment) for comparison with microarray data. Analysis of these eight genes using both techniques showed a high correlation (R = 0.826) in 24 hpt samples, indicating the high fidelity of the obtained data (Fig. 6, Supplementary Table 6).

During pathogen attack, recognition of pathogens or elicitors leads to a massive reprogramming of plant cells to activate and deploy defense responses (Bolton, 2009; Rojas et al., 2014). Assignment of recourses for immune response and biosynthesis of protective compounds leads to an increased demand for energy and carbon skeleton, which are supported by primary metabolic pathways (Bolton, 2009; Rojas et al., 2014).

Numerous studies have reported that photosynthesis is impaired when plants shift towards non-assimilatory metabolism in response to pathogen attack (Kangasjärvi et al., 2012; Major et al., 2010; Roberts and Paul, 2006). In the present study, 59 of 63 genes involved in light reaction of photosynthesis, 18 of 19 Calvin cycle genes and 14 of 18 tetrapyrrole biosynthesis genes were down-regulated by MSP1 treatment (Fig. 2, 3, Supplementary Table 3), indicating broad down-regulation of photosynthesis in response to MSP1 treatment. These results are in line with those of previous reports that biotic stress causes uniform down regulation of the transcription of genes encoding the major components of photosynthesis (Bilgin et al., 2010). Thirty of 45 genes (e.g., fatty-acid desaturase OsSSI2, Acyl-CoA-binding protein OsACBP1 and OsLACS1) involved in lipid synthesis and 73 of 83 genes related to protein synthesis were down-regulated in MSP1 treated samples, suggesting the activation of non-assimilatory metabolism in response to MSP1 treatment. Previously, it was observed that OsSSI2-knock down plants exhibit markedly enhanced resistance to blast and leaf blight diseases in rice and Arabidopsis (Liu et al., 2014a; Shah et al., 2001), indicating it is a negative regulator of stress tolerance in plants. Here, we observed down-regulation of OsSSI2 in response to MSP1 treatments, indicating an enhanced resistance by MSP1.

Inhibition of photosynthesis activity and collapse of assimilatory metabolism results in a metabolic transition in which the demand for carbohydrates and energy is achieved through increased activities of respiratory metabolism and other types of metabolism (Smedegaard-Petersen and Stolen, 1981). In the present study, three main pathways of respiration were activated during the resistance response, including glycolysis, the TCA cycle, and the electron transport chain. Three of eight genes encoding glycolysis-related proteins, ten of 18 genes involved in the TCA cycle, and seven of 16 genes involved in mitochondrial electron transport were up-regulated by MSP1 treatment, revealing a transcriptional change in the three main pathways of respiration in response to MSP1. Totally, 23 of 28 genes involved in lipid degradation and 115 of 198 genes involved in protein degradation were up-regulated by MSP1 treatment, providing an alternative respiratory substrate for MSP1-treated rice. Moreover, the up-regulation of protein and lipid degradation indicates that senescence processes in rice leaves may be highly elevated and nutrients released are mobilized from the leaves to other parts of the plants to support the immune response.

PTI is a form of basal defense that confers resistance to broad-spectrum pathogens in plants (Boller and Felix, 2009; Niks et al., 2015). Pattern recognition receptors (PRRs) trigger PTI, a relatively weak immune response restricting pathogen growth by recognizing highly conserved PAMPs. In this study, most of the 31 up-regulated RLK genes were leucine-rich repeat (LRR) RLKs and S-locus glycoprotein like RLKs, indicating potential RLKs for MSP1 triggered immunity. In FLS2, EFR, and Xa21, there is a large ectodomain including 28 LRRs, 21 LRRs, and 23 LRRs, respectively (Gómez-Gómez and Boller, 2000; Song et al., 1995; Zipfel et al., 2006), which are directly involved in elicitor binding (Antolín-Llovera et al., 2012). The DUF26 RLKs, which are cysteine-rich RLKs, comprise a large subgroup of RLKs involved in plant defense and programmed cell death (Chen et al., 2004). In the present study, the gene encoding OsDUF26 (also named OsRLK) was highly up-regulated by MSP1 treatment, which is in agreement with the results of a previous study in which OsDUF26 was induced by M. oryzae infection and fungal elicitor at both the RNA and protein levels (Kim et al., 2004, 2009b).

Following RLKs recognition, MAPK cascades are activated. In the present study, 14 MAPKKK, one MEK, and two MPK genes were up-regulated by MSP1 (e.g., OsMPK5 gene and MAPKKK1 gene). OsMPK5, also known as OsBIMK1 and OsMPK3, has been reported to be involved in disease resistance responses (Song and Goodman, 2002) and to positively regulate the JA signaling pathway and plant resistance to a chewing herbivore in rice (Wang et al., 2013). Moreover, MAPKKK1 (also known as OsEDR1 and OsACDR1) is a potential positive regulator of fungal disease resistance (Kim et al., 2009a) and a negative regulator of rice bacterial resistance via activation of ethylene biosynthesis (Shen et al., 2011). These results indicate that MAPKs play a vital role in the transduction of MSP1 induced signal to disease resistance in rice.

As a downstream signal of MAPKs, numerous WRKY transcription factors are involved in the development and initiation of stress and defense responses. Among 30 up-regulated genes encoding WRKYs, OsWRKY23, OsWRKY24, OsWRKY28, OsWRKY45, OsWRKY51, OsWRKY53, OsWRKY62, OsWRKY70, OsWRKY71, OsWRKY72, and OsWRKY76 have been thoroughly investigated for their roles in rice innate immune responses. Several previous transcriptional analyses showed expression of WRKYs were increased in response to attack by M. oryzae (e.g., OsWRKY45, OsWRKY47, OsWRKY53, OsWRKY55, OsWRKY62 and OsWRKY71), and over-expression of some of these genes enhanced resistance to rice blast infection (Bagnaresi et al., 2012; Chujo et al., 2007; Ryu et al., 2006; Shimono et al., 2012; Wei et al., 2013). Moreover, four WRKYs (OsWRKY51, OsWRKY62, OsWRKY71, and OsWRKY76) were reported to be involved in resistance against X. oryzae pv. oryzae infection in rice (Hwang et al., 2016; Liu et al., 2016, 2007; Peng et al., 2008). For instance, OsWRKY45 plays important roles in rice resistance to both M. oryzae and X. oryzae pv. oryzae via SA hormone signaling (Shimono et al., 2012) and negatively modulates resistance against the brown planthopper, Nilaparvata lugens (Huangfu et al., 2016). Moreover, OsWRKY45 and Pb1 interaction contributes to blast resistance through protection of OsWRKY45 from ubiquitin-proteasome system degradation (Matsushita et al., 2013). These results indicate that complicated signaling of WRKY transcription factors could be initialized by the MSP1 treatment that can be transduced to defense responses.

Phytohormones are central regulators of plant defense. In the present study, the JA, ET, and ABA signaling pathways were significantly upregulated by MSP1 treatment, while only three genes involved in the SA signaling pathway were detected in MSP1 treated samples, which is in accordance with our previous finding that JA and ABA positively regulate MSP1 induced cell death (Wang et al., 2016). It is well known that JA regulates growth development and stress responses, especially defense responses to herbivores and necrotrophic pathogens (Browse, 2009). Previous studies in rice have revealed that JA signaling was involved in rice basal defense against bacterial and fungal pathogens (Yamada et al., 2012; Yang et al., 2012, 2013). Moreover, ABA was reported to be both a positive and negative regulator of rice disease resistance, depending on the pathogen type and infection stage (Yang et al., 2013). There are also several lines of evidence that suggest ET confers broad-spectrum resistance against the fungal infection (Helliwell et al., 2013). Additionally, numerous studies have shown that SA acts as a positive regulator of immunity against biotrophic pathogens and hemibiotrophic pathogens (Berens et al., 2017; Glazebrook, 2005). However, rice plants did not show any increase in SA after inoculation with M. oryzae (Silverman et al., 1995). Limitation of the microarray may be one of the reasons why only three SA signaling genes were identified. Moreover, a recent study revealed that overexpression of MSP1 in rice confers broad-spectrum resistance through modulation of the SA- and JA-mediated signaling pathways (Hong et al., 2017). In addition, a study of the newly identified CP protein, SsCP1, revealed that the SA signaling pathway is involved in CP-mediated plant defense (Yang et al., 2018). Together, these results indicate that intricate networks of phytohormone signaling are activated by MSP1 treatment, leading to increased rice immunity.

Overall, our study provides a broader picture of the major physiological reactions underlying MSP1-induced PTI signaling at the transcriptional level. Upon treatment with MSP1, elicitor perception and a series of signal transductions (including MAPKs, WRKY transcription factors, hormone signaling and PR proteins) were activated. Our results suggest MSP1 triggers the activation of plant defense responses against various pathogens.