Full-Length Transcriptome-Wide Characteristic and Functional Identiﬁcation of WRKY Family in Malus sieversii during the Valsa Canker Disease Response

: WRKY transcription factors are one of the largest families in plants, playing important roles in regulating plant immunity. Malus sievesii has abundant genetic diversity and can offer various and high-quality gene resources. In this study, 112 putative MsWRKY proteins were identiﬁed from a full-length transcriptome of M. sieversii during the Valsa canker disease (caused by Valsa mali ). The MsWRKY proteins were phylogenetically divided into three groups (I–III). Motif compositions of the MsWRKY proteins were clustered and ﬁfteen conserved motifs were observed. Expression pattern analysis showed that thirty-four MsWRKY transcripts strongly responded to the V. mali infection, demonstrating that MsWRKY transcripts might play different roles during the response. Functional identiﬁcations were subsequently conducted with transient expressions, demonstrating that MsWRKY16, MsWRKY21, MsWRKY70, MsWRKY74 and MsWRKY85 positively regulated the resistant response. Besides, the MsWRKY21 , MsWRKY70 and MsWRKY85 were dramatically induced by salicylic acid (SA), methyl-jasmonate acid (MeJA) and 1-aminocyclopropane-1-carboxylate (ACC), indicating that they play important roles in the regulatory resistance of V. mali infection. This work provides a comprehensive understanding of the WRKY family in M. sieversii and will build a foundation for future research of the potential disease resistances MsWRKY transcripts.


Introduction
WRKY transcription factor (TF) family is one of the largest families in plants and plays an important role, both positively and negatively, in transcription regulatory networks to respond to the biotic and abiotic stresses [1,2]. The first WRKY gene was identified in 1994, named SWEET POTATO FACTOR 1 (SPF1) [3]. The WRKY TF family contains 74, 109 and 113 members in Arabidopsis thaliana, Oryza sativa and Malus domestica (GDDH13) [4][5][6].
WRKY transcription factors contain the highly conserved amino acid sequence WRKYGQK domain in the N-terminal and zinc-finger motifs C2H2 (CX 4-5 CX [22][23] HXH) or C2HC (CX 7 CX 23 HXC) in the C-terminal, hence named WRKY TF family [7]. Recent research reported that the highly conserved heptapeptide WRKYGQK also had other forms, which are WRKYGKK, WKKYGQK, WRKYGQR and WRKYGEK [8,9]. WRKY transcription factors have sequence-specific DNA-binding activity, named W box (C/T) TGAC (T/C) element [10]. The WRKY transcription factors are classified into three main groups based on the number of WRKY domains and the type of zinc-finger [7]. WRKY proteins with two WRKY domains (in N-terminal and C-terminal, respectively) were classified into group I, while proteins with only one WRKY domain belong to group II or III. Generally, the greenhouse at 24 ± 2 • C in 16/8 h light/dark cycles. These apple seedlings were used for the transient transformation assay. The tissue-cultured apple seedlings of M. sieversii were cultured in a climate chamber with a day/night cycle of 16/8 h and at 24 ± 2 • C at a light intensity of 25 µmol m −2 s −1 . The fungal pathogen V. mali isolate EGI 1 was isolated from the M. sieversii in the Tianshan Wild Fruit Forest area of Xinjiang, China, and identified by Liu et al. [31]. The V. mali isolate EGI 1 was grown on the (potato dextrose agar) PDA plate at 25 • C for three days.

Identification of the WRKY Family Proteins in M. Sieversii
The transcriptome database of M. sieversii upon V. mali infection was derived from the laboratory of Professor Daoyuan Zhang, Xinjiang Institute of Ecology and Geography Chinese Academy of Sciences, which was published with the project accession number PRJNA687214 (https://dataview.ncbi.nlm.nih.gov/object/PRJNA687214?reviewer=dv e4c2uosci45ag7dnbunbrt44, accessed on 22 December 2020) [32]. For the sampling of the PacBio and Illumina sequencing, barks of twigs in M. sieverisii near the canker were separately harvested at the time points of 0, 1, 2, and 5 dpi and each sample contained three biological replicates. These uninfected bark samples at 0 dpi time point were collected for RNA extraction as controls. The hidden Markov models (PF03106) of WRKY protein was downloaded from pfam database (http://pfam.xfam.org/, accessed on 1 January 2019). Hmmscan program of HMMER package was used to search the potential WRKY protein form the transcriptome database of M. sieversii with parameters as "E cut-off 1e". The sequences were confirmed using the NCBI-CDD (https://www.ncbi.nlm.nih.gov/Struc ture/cdd/cdd.shtml, accessed on 9 May 2019) and SMART databases (http://smart.em bl-heidelberg.de/, accessed on 2 August 2019) [33]. After removing the redundant and incorrect predicted-protein sequences, the MsWRKY proteins were validated finally.

Classification and Phylogenetic Analysis of MsWRKYs
Multiple alignments were conducted by DNAMAN Version 9.0 with an identity protein weight matrix to obtain detailed information on the structure of MsWRKYs. Weblogo (https://weblogo.berkeley.edu/logo.cgi, accessed on 8 May 2019) was used to generate the conversed motif logo of MsWRKYs. MEME suite 5.1.1 (http://meme-suite.org, accessed on 8 May 2019) was employed to analyze the possible conversed motifs in MsWRKY protein sequences setting a maximum value of motifs to 15

RNA Extraction and qRT-PCR Analysis
The total RNAs used for the constructions of Illumina and PacBio sequencing libraries were used to implement the qRT-PCR assay, which was for validation of the consistency of RNA-seq analysis. Apple leaves in the transient transformation assay were collected and full ground with liquid nitrogen. Total RNAs were isolated using Plant RNA Kit (Omega, No. R6827, China). The first-strand cDNA was synthesized with 1 µg of total RNA from each sample using PrimeScript RT reagent Kit with gDNA Eraser (Takara, No. RR047Q, China). Gene-specific primers used for the validation of expression levels of MsWRKYs (Table S3), were designed using primer-BLAST (https://www.ncbi.nlm.nih.gov /tools/primer-blast/, accessed on 29 November 2018) and synthesized by Sangon Biotech (Shanghai, China). Real-time quantitative reverse transcription polymerase chain reaction (qRT-PCR) was performed with TB Green Premix Ex Taq II kit (Takara, No. RR820A, China) on CFX96 Real-Time PCR Detection System (Bio-Rad, CFX96, American). The reaction system contained 12.5 µL TB Green Premix Ex Taq (2×) buffer, 1 µL each of the primers (10 µM) and 2 µL template and 8.5 µL ddH 2 O. The thermal profile for qRT-PCR was as follows: preheating at 95 • C for 30s; cycling stage: 95 • C for 5 s, 60 • C for 30 s, 40 cycles; melt curve analysis was produced to confirm the specificity of amplification. The relative transcript abundances were analyzed with the 2 −∆∆CT method [34]. The EF-1α (Elongation factor 1-α) gene was used as an internal reference gene [35]. Each sample contained three biological replicates and each biological replicate contained three technical replicates. Statistical analysis was performed by one-way ANOVA method, Tukey's HSD test using SPSS v.18 software. The heatmap of each group of MsWRKYs was generated by TBtools (http://cj-chen.github.io/tbtools/, accessed on 5 May 2019) based on the FPKM values of the Illumina transcriptome database.

Transient Expression Assay
To transiently express the selected MsWRKY transcripts, pBI121-MsWRKYs-GFP overexpression vectors were constructed using an in-fusion HD cloning kit (Code No. 639648). The method of transient expression was performed followed as previously described with modifications [36]. The A. tumefaciens in 20% glycerol stock at −80 • C was recovered on a Luria-Bertani (LB) agar plate for 3-day incubation at 28 • C. A single colony from the plate was picked up to inoculate in 5 mL LB medium containing kanamycin (50 µg/mL) and rifampicin (50 µg/mL) for shaking (220 rpm) at 28 • C for 20-24 h. Moreover, 1 mL of the overnight grown bacterial solution was taken into the sterile 100 mL LB medium with antibiotics and continued to grow for 6-8 h until the OD 600 of 0.6. For pre-induction, the pellet of A. tumefaciens was suspended with AB-MES solution (17.2 mM K 2 HPO 4 , 8.3 mM NaH 2 PO 4 , 18.7 mM NH 4 Cl, 2 mM KCl, 1.25 mM MgSO 4 , 100 µM CaCl 2 , 10 µM FeSO 4 , 50 mM MES, 2% glucose (wt/vol), pH 5.5) and acetosyringone (AS, 200 µM) to OD 600 of 0.2 without antibiotics, then shaken (220 rpm) at 28 • C for 12-16 h. To perform the transient transfection procedure, the inducted A. tumefaciens were collected and resuspended with the solution AB-MES and MS medium (v/v: 1:1) containing the AS (100 µM). The bacterial solution used for transient transfection was adjustee to OD 600 of 0.8. The detached leaves were washed with ddH 2 O and disinfected the surface with alcohol cotton. The leaves were immersed in the transfection solution containing A. tumefaciens in the dark at 28 • C for 1 h. Subsequently, these washed leaves with ddH 2 O were placed in sterile dishes and preserved moisture. After a one-day cocultivation in the dark at 28 • C, the transiently expressed leaves were prepared for inoculations of V. mali.

Infection Assays with Transient Expression Leaves
The procedures for the infection assay used the method optimized by Liu et al. [31]. The mycelial plugs (diameter: 5 mm) were excised aseptically from the edge of the 3-dayscultured isolate EGI1 on PDA medium. After the 2-day transient expression, these leaves were wounded by a fabric pattern wheel and inoculated with these mycelial plugs. Leaves inoculated with ddH 2 O acted as controls. The mocks were that leaves were inoculated with A. tumefaciens for 2 days and then inoculated with V. mali. These inoculated leaves were placed in dishes tapped with parafilm for keeping moisture. After a 4-day inoculation of the V. mali, the lesion and leaf areas were measured by software ImageJ for lesion ratios (%) calculations. Each sample contained ten leaves with two lesions. Each sample contained three biological replicates and each biological sample contained three technical replicates. The mean and standard error (SE) were analyzed by the one-way ANOVA using Tukey's HSD method (p < 0.05) in intra-group comparisons compared with the Mock.

Hormone Treatment
For the hormone treatment, apple seedlings (M. domestics 'M26') were planted in the greenhouse condition, which was 25 ± 2 • C and a relative humidity of 50 ± 10%. Leaves were detached to use to inoculate with V. mali and spray in the solutions, which were SA solution (5 mM, dissolved in sterile ddH 2 O), MeJA solution (1 mM, dissolved in sterile ddH 2 O) and ACC solution (1 mM, dissolved in sterile ddH 2 O) [37,38]. After the inoculation of V. mali and hormone treatments, the leaf samples were immediately frozen in liquid nitrogen at 0 h, 1 h, 3 h, 6 h, 1 d, 2 d, 5 d timepoints, and stored at −80 • C for further use. The hormone treatment was implemented in a growth chamber at 24 • C in 16/8 h light/dark cycles and 60% relative humidity.

Identification of WRKY Gene Family in M. sieversii during Canker Disease Response
A total of 112 unique MsWRKY transcripts encoding the conserved WRKY domain were identified from the genome of M. domestica, named as MsWRKY1 to MsWRKY112. The MsWRKY106 had the largest amino acid sequence of 995 aa (MW: 11.32 kDa), while the MsWRKY55 had the smallest amino acid sequence of 80 aa (MW: 9.28 kDa). The isoelectric points (PIs) ranged from 4.81 (MsWRKY83)-9.99 (MsWRKY78). The highly conserved domain WRKYGQK was presented in 110 MsWRKY proteins, whereas the remaining two members contained WRKYGKK and WRKYGMK domains ( Table 1). The biochemical properties of MsWRKY proteins were predicted, including the grand average of hydropathicity (GRAVY), instability index, subcellular location, nuclear location signal (NLS) and phosphorylation site. One hundred and eight MsWRKY proteins were hydrophilic protein (GRAVY < −0.5); the remaining four MsWRKY proteins were hydrophobic (GRAVY > −0.5). The instability indexes of 105 MsWRKY proteins were more than 40, which showed that they were unstable proteins. The predicted subcellular location showed that MsWRKY proteins localized in the nucleus (75), chloroplast (27) and cytoplasm (10). It was inferred that the MsWRKY transcription factors outside of the nucleus might play transcriptional regulatory functions through the process of transmembrane transfer. NLS prediction analysis showed that 27 of 112 MsWRKY proteins were without NLSs. The phosphorylation site prediction analyzed by KinasePhos showed that MsWRKY proteins had abundant phosphorylation sites, including predominant serine (S), threonine (T) and tyrosine (Y) sites. There were 43 phosphorylation sites in the MsWRKY29 protein, while there were 2 phosphorylation sites in MsWRKY55 and MsWRKY60 separately (Table S1).

Phylogenetic Analysis and Classification of MsWRKY Proteins
Based on the number of WRKY domains and the zinc finger type, the putative WRKY proteins could be grouped I, II (a-e) and III [7]. Multiple alignment results showed that Group I (MsWRKY1-20) had two WRKYGQK domains and zinc-finger motifs of C2H2 type. Group II (MsWRKY21-95) had one WRKYGQK domain and two zinc-finger motifs of C2H2 type except for MsWRKY50/MsWRKY61 (WRKYGKK domain) and MsWRKY62 (WRKYGMK domain). Aligned with AtWRKY protein sequences, Group II of MsWRKY proteins were classified into the five subgroups (II-a, II-b, II-c, II-d II-e). Sixteen members of Group III (MsWRKY96-112) had one WRKYGQK domain and a zinc-finger motif of C2HC (Table 1). Additionally, the zinc-finger motifs of MsWRKY proteins were relatively conserved. Group I N, C and group II-c had the zinc-finger motif of CX 4 X 23 HXH, except for MsWRKY6 (CX 4 CX 24 HXH) and MsWRKY14 (CX 4 CX 24 HXH). Group II-a, II-b, II-d and II-e had the zinc-finger motif of CX 5 CX 23 HXH. Excluding MsWRKY100 and MsWRKY101 (CX 5 CX 23 HXC), the Group III MsWRKYs had zinc-finger of CX 7 CX 23 HXC (Table 1). To determine the phylogenetic relationship between all 112 MsWRKYs, the neighbor-joining (NJ) phylogenetic tree was constructed with 112 MsWRKYs, 65 AtWRKYs and 113 Md-WRKYs proteins (Figure 1). The NJ phylogenetic tree of MsWRKY full-length proteins was clustered into three groups (I, II and III). The Group II diverged into three clades, of which were clade IIa/IIb, clade IId/IIe and clade IIV. It indicated that close evolutionary relationships were between subgroup IIa and IIb; also, IId and IIe.

Conserved Motif Composition Analysis of MsWRKY Proteins
The conserved motifs of 112 MsWRKY proteins were predicted using MEME suite 5.

Transmembrane Domain (TMD) and Conserved Domain in MsWRKY Proteins
Through TMHMM and SMART prediction analysis, four MsWRKY proteins with transmembrane domains were identified (Table S1). The TMDs of MsWRKY7, MsWRKY97 and MsWRKY111 were in the N-terminus, while the TMD of MsWRKY44 protein was in the C-terminus. Additionally, six specific conserved domains were identified in certain MsWRKY proteins, including RPW8, rad50, UPF0242, nuc hydro, Plant zn clust and PLN03210 domains using NCBI CDD database. There were two RPW8 domains in the N-terminal of the WRKY domain of the MsWRKY7 protein.

Transmembrane Domain (TMD) and Conserved Domain in MsWRKY Proteins
Through TMHMM and SMART prediction analysis, four MsWRKY proteins with transmembrane domains were identified (Table S1). The TMDs of MsWRKY7, MsWRKY97 and MsWRKY111 were in the N-terminus, while the TMD of MsWRKY44 protein was in the C-terminus. Additionally, six specific conserved domains were identified in certain MsWRKY proteins, including RPW8, rad50, UPF0242, nuc hydro, Plant zn clust and PLN03210 domains using NCBI CDD database. There were two RPW8 domains in the N-terminal of the WRKY domain of the MsWRKY7 protein.

Expression Patterns of MsWRKY Transcripts upon the V. mali Infection
To confirm whether MsWRKYs are involved in response to the infection of V. mali, expression patterns of the 112 MsWRKY transcripts were analyzed upon the infection of V. mali. (Figure 4A, Table S2). The expression levels were represented using FPKM values at 0, 1, 2, 5 dpi. The transcripts were considered to be differentially expressed depending on the cut-off that was |log 2 (Fold change)| ≥ 1 with the Q-value < 0.05. Thirty-three MsWRKY transcripts were significantly differentially expressed during the response stage, indicating that these MsWRKY transcripts may play important roles in response to V. mali. The differentially expressed MsWRKY transcripts were clustered into three clusters based on expression patterns. Cluster 1 contained 10 MsWRKY transcripts that were all up-regulated at 1, 2 and 5 dpi, including MsWRKY6, MsWRKY27, MsWRKY29, MsWRKY30, MsWRKY34, MsWRKY35, MsWRKY51, MsWRKY52, MsWRKY67, MsWRKY85. Cluster 2 contained 10 MsWRKY transcripts, which were all significantly down-regulated at 1, 2 and 5 dpi, including MsWRKY7, MsWRKY12, MsWRKY32, MsWRKY70, MsWRKY71, MsWRKY94, MsWRKY101, MsWRKY102, MsWRKY103, MsWRKY105. In cluster 3, the expression levels of MsWRKY10, MsWRKY16, MsWRKY28, MsWRKY33, MsWRKY75, MsWRKY90 were significantly up-regulated at 1 dpi, and then down-regulated at 2 and 5 dpi, except for MsWRKY10 and MsWRKY33. MSWRKY13 and MsWRKY86 showed the same expression patterns, which were significantly increased at 2 and 5 dpi. The expression level of MsWRKY108 was significantly decreased at 1 and 2 dpi, and then were significantly increased at 5 dpi ( Figure 4B). Overall, the results show that MsWRKYs might play important roles in the response to the V. mali infection in M. sieversii.

Expression Patterns of MsWRKY Transcripts upon the V. mali Infection
To confirm whether MsWRKYs are involved in response to the infection of V. mali, expression patterns of the 112 MsWRKY transcripts were analyzed upon the infection of V. mali. (Figure 4A, Table S2). The expression levels were represented using FPKM values at 0, 1, 2, 5 dpi. The transcripts were considered to be differentially expressed depending on the cut-off that was |log2 (Fold change)| ≥ 1 with the Q-value < 0.05. Thirty-three MsWRKY transcripts were significantly differentially expressed during the response stage, indicating that these MsWRKY transcripts may play important roles in response to V. mali. The differentially expressed MsWRKY transcripts were clustered into three clusters based on expression patterns. Cluster 1 contained 10 MsWRKY transcripts that were all up-regulated at 1, 2 and 5 dpi, including MsWRKY6, MsWRKY27, MsWRKY29, MsWRKY30, MsWRKY34, MsWRKY35, MsWRKY51, MsWRKY52, MsWRKY67, MsWRKY105. In cluster 3, the expression levels of MsWRKY10, MsWRKY16, MsWRKY28, MsWRKY33, MsWRKY75, MsWRKY90 were significantly up-regulated at 1 dpi, and then down-regulated at 2 and 5 dpi, except for MsWRKY10 and MsWRKY33. MSWRKY13 and MsWRKY86 showed the same expression patterns, which were significantly increased at 2 and 5 dpi. The expression level of MsWRKY108 was significantly decreased at 1 and 2 dpi, and then were significantly increased at 5 dpi ( Figure 4B). Overall, the results show that MsWRKYs might play important roles in the response to the V. mali infection in M. sieversii.  To confirm the contributions of the MsWRKY transcripts on Valsa canker resistance, twelve significantly differential expression MsWRKY transcripts in different groups were selected to quantify by qRT-PCR at different time points after infection. At 1 dpi, the relative expression level of MsWRKY21, MsWRKY28, MsWRKY30 and MsWRKY90 were significantly increased. At 5 dpi, the relative expression levels of MsWRKY5, MsWRKY13, MsWRKY27, MsWRKY35, MsWRKY67, MsWRKY75, MsWRKY102 were significantly upregulated. On the contrary, at 1 dpi the relative expression levels of MsWRKY32 and MsWRKY102 were significantly down-regulated. At 5 dpi, the relative expression level of MsWRKY28 was significantly decreased. The qRT-PCR results of the twelve MsWRKY transcripts were consistent with the FPKM values from RNA-seq data ( Figure 5). regulated. On the contrary, at 1 dpi the relative expression levels of MsWRKY32 and MsWRKY102 were significantly down-regulated. At 5 dpi, the relative expression level of MsWRKY28 was significantly decreased. The qRT-PCR results of the twelve MsWRKY transcripts were consistent with the FPKM values from RNA-seq data ( Figure 5).  The relative expression levels were showed on the left y-axis. MsEF1-α was used as an internal reference gene. Each sample contained three biological replicates and each biological replicate contained three technical replicates. The relative expression data were analyzed by one-way ANOVA method, Tukey's HSD test. Lowercase letters represented significant differences (p < 0.05).

Functional Identification of WRKY Family Members in M. sieversii Responded to V. mali
Dependent on the different expression patterns of these differentially expressed MsWRKY transcripts, six MsWRKY transcripts were cloned to identify the resistance to V. mali, including MsWRKY16, MsWRKY21, MsWRKY70, MsWRKY74, MsWRKY85 and MsWRKY105. To identify the potential resistance of MsWRKY proteins that responded to the infection of V. mali, we performed disease resistance verifications using the transient expression method mediated by A. tumefaciens [36]. Six MsWRKYs transcripts were separately transiently transfected into leaves using A. tumefaciens containing pBI121-MsWRKYs-GFP vectors. These relative expressions of MsWRKYs were correspondingly increased in leaves, after 4 days of the transiently overexpression by A. tumefaciens ( Figure 6A), which indicated that these MsWRKY transcripts were successfully overexpressed in the leaves. In the dis-ease resistance test, after 4-day inoculations of V. mali, the lesion ratios of the MsWRKY16, MsWRKY21, MsWRKY70 and MsWRKY85 were significantly decreased compared with the mock (Figure 6B). These results suggested that MsWRKY16, MsWRKY21, MsWRKY70, MsWRKY74 and MsWRKY85 could be positive regulators for the resistance to the Valsa canker. Yet, the mean lesion ratio of MsWRKY105 transiently expressed leaves was similar to the mock ( Figure 6B), indicating that the MsWRKY105 was no resistant to the canker disease. Interestingly, MsWRKY70 and MsWRKY74 obviously induced severe necrosis at 6 days after transient expression compared with WT ( Figure 6B).  Plant hormones play important roles in signaling networks of immune defense involved in plant responses to pathogen infection [39]. To detect whether these six MsWRKY transcripts were induced by different phytohormones during the infection of V. mali, qRT-PCR was performed to determine the expression levels ( Figure 6C). As the positive regulator for the canker disease resistance, the expression of MsWRKY16 was significantly up-regulated under the ACC + V. mali, MeJA + V. mali and the SA + V. mali conditions. With the ACC + V. mali treatment, the relative expression of MsWRKY16 significantly decreased at 6 h and 2d. Under the MeJA + V. mali condition, the expression of the MsWRKY16 peaked at 6h, yet was down-regulated after 2 dpi, as well as up-regulated at 5 dpi. The expression of MsWRKY16 was continuously increased after the SA + V.

Discussion
Wild apple is an ancestor of the cultivated apple distributed in Tianshan Forest Mountain from Central Asia to West Europe along the Silk Road and holds underlying potential for the germplasm improvement of apple, dependent on genetic diversity [26,40]. However, M. sieversii was suffering from attacking by the canker disease caused by necrotrophic pathogen V. mali, resulting in a decrease in the area of this population [31]. The pathogen V. mali not only infects M. sieversii, but also causes extensive damage and even death to apples in eastern Asia with severe yield losses [41]. Thus, it is important to explore the disease resistance gene resources in M. sieversii upon the infection of V. mali, which could aid in the improvement of apple germplasm resources.
The WRKY family is among the ten largest families of transcription factors in higher plants and major positive or negative regulators in plant immunity [2]. WRKY families were identified in unicellular green alga Chlamydomonas rein (1), the moss Physcomitrella patens (37), Arabidopsis (74), Glycine max (197) and Populus triuncatula (103) [42]. In apple, based on the apple genome database of "Golden Delicious" [43], 116 WRKY transcription factor genes [44] and 119 MdWRKY genes in response to the fungal pathogen (Alternaria alternata) and phytohormone treatment [45] were separately identified. Recent research reported that a high-quality de novo assembly of the apple (GDDH13) was produced [4]. Based on this database, 114 WRKY transcription factor genes in M. domestica were identified [6]. In this study, a total of 112 MsWRKY transcripts were firstly identified from the full-length transcriptome data of M. sieversii upon the V. mali infection. The number of members of the WRKY family of M. sieversii identified from the full-length transcriptome responded to V. mali infection, is almost equal to the number in M. domestica (GDDH13) (Figure 1). Additionally, the number of members of the MsWRKY family expanded compared with V. rein, P. patens, Arabidopsis and Populus, suggesting that the number of WRKY genes has an expansion from lower plants to higher.
The phylogenetic analysis demonstrated that MsWRKYs could be classified into three groups, among which group II was further divided into five subgroups. The number of MsWRKYs group II accounted for the largest proportion of 67%, which was consistent with former reports regarding the WRKYs in Arabidopsis, Populus and M. domestica [2,4,6]. The NJ tree of MsWRKY proteins indicates that groups IIa and IIb, and groups IId and IIe have close phylogenetic relationships separately. The monophyletic clades of three groups seemed to be formed major five lineages regarding group I, group IIa + IIb, group IId +IIe, group IIc and group III (Figure 1). Fifteen conserved motifs were identified in the different groups of MsWRKYs ( Figure 2B), indicating that the conserved motif composition of the group of WRKY family could be the cause of MsWRKY transcription factor with various biological functions. WRKY proteins are defined by the conserved WRKY domain, the amino acid sequence is WRKYGQK [7]. In a few WRKY proteins, the WRKY domain has some variants, including WRKYGKK, WRKYGEK, WRKYGMK and WKKYCEDK [46,47]. In this study, two variants were mainly observed in the subgroup IIc of MsWRKY proteins, which were WRKYGKK and WRKYGMK. It is inferred that MsWRKY50/MsWRKY61 and MsWRKY62, which possessed variants in WRKY domains (Table 1), might have a variety of biological functions. In N. tabacum, NtWRKY12 specifically binds the promotor of PR1a genes contained in theWK-box (TTTTCCAC), not W-box TTGAC(C/T), thereby improving disease resistance in the plant [48]. It could be inferred that MsWRKY50 and MsWRKY61 might have similar cis-acting elements, of which one is NtWRKY12. Furthermore, the biochemical properties of 112 MsWRKY proteins were characterized, including the grand average of hydropathicity (GRAVY), instability index, subcellular location, and phosphorylation site. The molecular weight of MsWRKY55 is the smallest of all, and consists of 80 amino acids. MsWRKY106 is the largest one of all, and consists of 995 amino acids. Twenty-seven MsWRKY proteins are without NLSs, may need other corresponding ligands to transport into the nucleus for biological functions (Table S1). MAP kinase pathways are involved in regulating WRKY transcription activation activity [16,49,50]. OsWRKY53 is a downstream target gene of OsMKK4-OsMPK3/OsMPK6, of which multiple clustered serine-proline residues (SP cluster) can be phosphorylated by OsMPK3/OsMPK6 to regulate the resistance to rice blast [51]. The abundant phosphorylation sites in MsWRKY proteins were observed (Table S1), speculating that they may play important roles in biological functions. The membrane-bound transcription factor is a kind of transcription factor with the C-terminal transmembrane domain (TMD). Under normal conditions, it localizes in the cytoplasm, while under stress conditions, the cleavage of the TMD allows the transcription factor to translocate to the nucleus [52,53]. In this study, four putative membrane-bound MsWRKY transcription factors were identified (Table S1). The conserve domain identification result showed that MsWRKY7, MsWRKY28, MsWRKY29, MsWRKY70-82 and MsWRKY106 had different specific conserved domains (Figure 3), inferring that these MsWRKY might not only be involved in transcriptional regulation as a transcription factor, but also play other roles in biological functions. Among these, the expression level of the MsWRKY7 transcript was significantly down-regulated at 5 dpi suggested that it might be involved in the response to V. mali infection ( Figure 4B). The RPW8 domain can mediate defense resistance and cell death in the cell nucleus and cytoplasm, respectively, and finally drive broad-spectrum resistance to powdery mildew [54,55]. There are also two disease resistance RPW8 domains in the N-terminus of MsWRKY7 protein ( Figure 3). It speculated that MsWRKY7 containing two RPW8 domains might also have the resistant function to powdery mildew.
The full-length RNA-seq transcriptome database was performed to profile MsWRKY transcripts expression patterns to respond to the Valsa canker disease. Thirty-four MsWRKY transcripts exhibited differential expression patterns upon V. mali infection ( Figure 4B). Based on the expression patterns of the RNA-seq and qRT-PCR data, it is indicated that the differentially expressed MsWRKY transcripts may perform diversified functions in response to Valsa canker disease ( Figure 5). Based on the different expression patterns, the selected six MsWRKY transcripts subsequently conducted the functional identification with a transient expressions assay. Among these, MsWRKY16, MsWRKY21, MsWRKY70, MsWRKY74 and MsWRKY85 were positive regulators in canker disease resistance ( Figure 6). In the previous study, numerous WRKY genes have been identified as functions of disease resistances and clarified the mechanisms. AtWRKY18, AtWRKY40 and AtWRKY60 exhibit a response to the infection of P. syringae and B. cinerea [56]. In Arabidopsis, AtWRKY18 and AtWRKY40 play major roles in PAMP-triggered basal defense in defense to the powdery mildew (caused by fungus: Golovinomyces orontii) during the early stages of infection [57]. The elaborate functions and mechanisms of the canker disease resistance MsWRKY transcripts should be implemented for further research.
The WRKY family plays an important in the phytohormones synthesis and signaling pathway to immune defense. The co-regulation of AtWRKY33 and AtMYB46 mediates the coordination of SA and JA/ET pathways to orchestrate the defense against the phytopathogen [58]. The AtWRKY18, AtWRKY40 and AtWRKY33 act as repressors for PAMP-induced secreted peptide 3 (PIP3) to turn on the immunity via regulating SA and JA biosynthesis and signaling [59]. EF-Tu receptor (EFR) activates the SA and JA synthesis to defense against the infection, through the direct phosphorylation of a receptorlike cytoplasmic kinase (BIK1) and WRKY33, WRKY50 and WRKY57 [60]. Thus, the phytohormone-induced expression pattern of these disease resistance-related MsWRKY transcripts was verified. Upon the hormone treatments of ACC + V. mali, MeJA + V. mali and SA + V. mali, MsWRKY16, MsWRKY21, MsWRKY70 and MsWRKY85 showed different expression patterns ( Figure 6C).
The expressions of MsWRKY16 and MsWRKY105 were induced through the treatments of MeJA + V. mali and SA + V. mali. The expressions of MsWRKY21 and MsWRKY85 were induced in the treatment of MeJA + V. mali. It speculated that these MsWRKYs could be key regulatory genes in the ET, JA and SA signaling pathways for defense of the V. mali infection. The potential functions and disease resistance regulatory mechanism of these MsWRKY TFs will be experimentally illuminated in our further research.

Conclusions
One hundred and twelve MsWRKY transcription factors were identified and divided into three groups (I, II, III) based on the phylogenetic analysis. Fifteen conserved motifs were identified in the different groups of MsWRKYs, of which compositions of similar motifs were usually clustered into the same clade. Based on the Illumina sequencing database, 34 differentially expressed MsWRKY transcripts may play key roles in response to the V. mali infection. Furthermore, the differentially expressed MsWRKY16, MsWRKY21 MsWRKY70, MsWRKY74 and MsWRKY85 were identified to play positive roles in resistance to canker disease. Besides, these MsWRKYs might contribute to the regulation of resistance and susceptibility of Valsa canker disease involving the SA, MeJA and ET signaling pathways. These results will provide reliable indexes of screening candidate MsWRKY transcripts and build a foundation for further investigation in M. sieversii to respond to the Valsa canker disease.
Supplementary Materials: The following are available online at https://www.mdpi.com/article/1 0.3390/f12060790/s1, Table S1: The list of the 112 MsWRKY transcripts identified in this study with properties of the predicted MsWRKY proteins.