Using RNAi, we previously generated CmWRKY15-1-silenced lines of ‘C029’ and CmWRKY15-1-OE lines of ‘Jinba.’ The white rust-resistant chrysanthemum cultivar ‘C029,’ the susceptible cultivar ‘Jinba,’ and various transgenic lines were used as experimental materials. These plant materials were provided by the Forestry College of Shenyang Agricultural University, Shenyang, China, and all experiments were conducted there in 2020 and 2021. Seedlings at the six- to eight-leaf stage were placed in a potting soil mixture and grown in a greenhouse under fluorescent lights for 2 weeks at 25 ± 3^°C.

Teliospores of P. horiana were collected from the abaxial sides of chrysanthemum leaves infected with white rust and placed in 1 ml sterile water. The concentration of the teliospores was adjusted to 40–60 per field of vision under a BA400 microscope (Motic, Xiamen). Teliospores were cultured according to the protocol by Takatsu (Takatsu et al., 2008). A drop of teliospore suspension was removed with a sterile straw and dropped onto a glass substrate. This was placed on a U-shaped rod in a culture dish covered with wet filter paper. The teliospores were allowed to germinate at 18–21°C for 24 h and resuspended in sterile water with 0.05% (w/v) Tween 20 (pH 4–6.5). The abaxial sides of leaves were sprayed evenly with the solution, and the plants were covered with plastic film and transferred to the dark and 90% relative humidity. After 16–24 h, infected plants were transferred to a growth chamber maintained at 17°C and 50% relative humidity. Leaves were harvested from WT, OE, and RNAi plants at 0 day (before inoculation), and at 24, 48, and 72 h after treatment with P. horiana. Leaves were wrapped in aluminum foil, immediately frozen in liquid nitrogen, and stored at –80^°C for sequencing experiments and to measure physiological indices.

H₂O₂ content in WT and CmWRKY15-1-OE lines of ‘Jinba’ leaves was measured at 0, 24, 48, and 72 h after inoculation with P. horiana. Leaves were immersed in diaminobenzidine (DAB) staining solution (1 mg/ml 3,3′-DAB, 10 mM Na₂HPO₄, and 0.05% Tween 20, pH 3.8) in the dark for 4–8 h and decolorized in 95% (w/v) ethanol. The intensity of brown coloration indicates H₂O₂ content. Quantitative H₂O₂ measurements were performed using an Amplex Red hydrogen peroxide/peroxidase assay kit (Solarbio Science and Technology, Beijing). Absorbance was measured using a Tecan infinite F200/M200 plate reader at 560 nm. H₂O₂ concentration is indicated in μmol/g fresh weight. The experiment was performed with at least three independent replicates.

Frozen leaf samples of WT and CmWRKY15-1-OE lines of ‘Jinba’ were used to measure the activity of defense enzymes, namely, PAL, POD, CAT, and SOD. SOD and POD activities were measured following the protocol by Sun et al. (2013). CAT and PAL activities were measured following the protocol by Kwon and Liu (Kwon and Nguyen, 2003; Liu et al., 2005). All experiments were performed with three biological replicates.

The Illumina HiSeq-2000 platform (Illumina Inc., San Diego, CA, United States) by Wuhan MetWare Biotechnology Co., Ltd., (¹ Wuhan, China) was used for RNA-seq. A total of 12 high-quality RNA samples [six from before inoculation (0 h) and six from after inoculation] were prepared from the leaves of WT and CmWRKY15-1-RNAi lines of ‘C029.’ First, samples were obtained from uninoculated WT and CmWRKY15-1-RNAi lines of ‘C029’ (0 h). Three samples were collected per strain for a total of six samples (labeled W0 and R0). WT and CmWRKY15-1-RNAi lines were then sampled at 24, 48, and 72 h after inoculation with P. horiana. Samples at specific time points were combined to form six inoculated samples, labeled WJ and RJ. Because most eukaryotic mRNAs have poly-A tails, mRNAs were enriched using oligo(dT) magnetic beads. Examination of RNA quality, cDNA library construction, and Illumina sequencing library construction were performed as previously described (Ana et al., 2016; Chuang et al., 2020).

We assembled pathogen-sequence-free reads from the above RNA-seq data into contigs using the Trinity program. All contigs were analyzed with RNA-seq by Expectation-Maximization software to map the clean reads of each sample to the transcriptional reference sequence obtained via Trinity splicing (Michael et al., 2014). Next, we obtained the raw counts for each unigene and used the DESeq R package to identify DEGs between W0 vs. WJ, R0 vs. RJ, and WJ vs. RJ. The parameter standards included a cut-off of |log₂(fold change)| > 1 and a false discovery rate < 0.05. Finally, we analyzed unigene sequences from the NCBI non-redundant (NR) database and euKaryotic Ortholog Groups (KOG) databases and annotated DEGs between W0 vs. WJ, R0 vs. RJ, and WJ vs. RJ using the Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway databases.

To confirm the RNA-seq results, seven disease-resistance DEGs in the SA pathway were selected for verification by quantitative reverse transcription PCR (RT-qPCR). Primers used are shown in Supplementary Table 4. Leaves of the ‘C029’ WT and CmWRKY15-1-RNAi lines were harvested at 0 and 48 h after treatment with P. horiana. A PrimeScript™ RT reagent Kit with gDNA Eraser (Perfect Real Time), PrimeScript™ II 1st Strand cDNA Synthesis Kit, and SYBR^® Premix ExTaq™ II were purchased from Takara. RT-qPCR was performed using an SYBR^® Premix ExTaq™ II, following the protocols of the manufacturer. The 2^–ΔΔCT method was used to calculate the expression levels of these disease resistance DEGs (Livak and Schmittgen, 2008). All experiments were performed with three biological replicates.

One DEG, CmNPR1, was identified by sequencing. The coding sequence (CDS) of CmNPR1 was cloned from ‘C029’ by PCR with CmNPR1 forward (F) and reverse (R) primers and was cloned into pEASY^®-T1 using a Simple Cloning Kit (TransGen Biotech, Beijing). Primers are shown in Supplementary Table 4. This vector was transformed into Escherichia coli DH5α (Aidlab) and clones harboring the inserts were identified for subsequent experiments. The cDNA and deduced protein sequences from CmNPR1 were analyzed using DNAMAN software. The online tool ExPASy² was used to predict the physicochemical properties of CmNPR1, and the SOPMA³ website and Phyre² database were used to predict protein structure.⁴ Finally, MEGA 5.0 software was used for the phylogenetic analysis of CmNPR1 proteins.

To investigate the potential function of CmNPR1, the silencing vector pTRV2-CmNPR1 was constructed, and expression of the gene was silenced using the virus-induced gene silencing (VIGS) technique in ‘C029.’ Leaves of WT and TRV-NPR1 lines were sampled at 0, 24, 48, and 72 h after inoculation with P. horiana. The relative expression levels of CmPR1, CmPR2, and CmPR10 were quantified by RT-qPCR. All experiments were performed with three biological replicates. Primers are shown in Supplementary Table 4.

The primers NPR1SP1, NPR1SP2, and NPR1SP3 were synthesized by hiTAIL-PCR following the protocol of the manufacturer (Chen and Liu, 2007). Through genome walking, we obtained the 1,685-bp promoter region of CmNPR1 upstream of the ATG codon. Promoter elements of CmNPR1 were analyzed using the Plant CARE⁵ and PLACE⁶ websites. The CmNPR1 promoter was amplified by PCR with specific primers CmNPR1.1-F/R, CmNPR1.2-F/R, and CmNPR1.3-F/R (based on the predicted cis-acting elements in this promoter) using ‘C029’ DNA as a template. The 1,685-bp CmNPR1 promoter region was obtained, along with truncated sequences 850 and 450 bp long. Primers used are shown in Supplementary Table 4. These sequences were ligated into the pCAMBIA1301 vector, replacing the 35S promoter, to form the recombinant expression vectors pNPR1.1:GUS, pNPR1.2:GUS, and pNPR1.3:GUS, respectively. The vectors were introduced into Agrobacterium tumefaciens strain EHA105. Combination vectors were transformed into ‘C029’ via Agrobacterium-mediated transient transformation. Leaves of WT and transgenic plants at 0, 24, 48, and 72 h after treatment with P. horiana were subjected to β-glucuronidase (GUS) staining using a GUS staining kit (Real-Times Biotechnology Co., Ltd., Beijing, China). All experiments were performed with three biological replicates.

Yeast one-hybrid (Y1H) assays were conducted using the Matchmaker™Gold Yeast One-Hybrid System (Clontech, Palo Alto, CA, United States) for experimental validation. The full-length 1,685-bp CmNPR1 promoter was ligated into the pAbAi vector, generating pNPR1-AbAi. The linearized vector was digested with a single endonuclease and transfected into yeast strain Y1HGold and used as bait. The CDS of CmWRKY15-1 was cloned into pGADT7 to generate the prey vector pGAD-CmWRKY15-1; the pGADT7 vector and bait constructs were transformed into yeast cells as a negative control. Subsequently, all constructs were transformed into the Y1HGold strain. Finally, the colonies were plated on synthetic defined (SD) medium without leucine (SD/–Leu) and SD/–Leu supplemented with 800 ng/ml aureobasidin A (AbA) and allowed to grow for 3 days at 28°C. The CmNPR1 promoter was divided into two segments based on analysis of the cis-acting elements PA1 and PA2. PA1 contains the critical regulatory W-box element of the CmNPR1 promoter. Two fragments were cloned into the pAbAi vector to produce new recombinant vectors pA1-AbAi and pA2-AbAi. These two bait vectors, and the prey vector pGAD-CmWRKY15-1, were transferred into Y1HGold. The transformants were grown on medium with 800 ng/ml SD/–Leu with AbA and cultivated at 28–30°C for 3 days. All experiments were performed with three biological replicates. Primers are shown in Supplementary Table 4.

Bimolecular fluorescence complementation (BiFC) vectors were constructed for the transient transformation of Nicotiana benthamiana. The CDSs of CmNPR1 and CmWRKY15-1 were cloned into the p2YC and p2YN vectors using the BiFC-NPR1-F1/BiFC-NPR1-R1 and BiFC-WRKY15-1-F1/BiFC-WRKY15-1-R1 primers via PacI/SpeI, generating p2YC-CmNPR1 and p2YN-CmWRKY15-1, respectively. Cells of A. tumefaciens strain GV3101 carrying the BiFC constructs were infiltrated into 4- to 5-week-old N. benthamiana leaves. Yellow fluorescent protein (YFP) signals in the N. benthamiana leaf cells were observed under a laser scanning confocal microscope. All experiments were performed with three biological replicates. Primers are shown in Supplementary Table 4.

Data were analyzed using ANOVA and t-tests to determine significant differences with SPSS 24.0 software.

To determine whether ROS are involved in resistance to P. horiana, we performed DAB staining to measure the H₂O₂ content in leaves. More H₂O₂ accumulated in WT than in CmWRKY15-1-OE lines after inoculation with P. horiana (Figure 1A). We also quantified endogenous H₂O₂ content. The levels of H₂O₂ increased over time and peaked at 48 h in WT and CmWRKY15-1-OE lines after inoculation with P. horiana. Compared with the WT, CmWRKY15-1-OE lines had lower H₂O₂ levels after inoculation with P. horiana, which is consistent with the results of DAB staining (Figure 1B). This result indicated that OE CmWRKY15-1 enhanced plants’ resistance to P. horiana by inhibiting the accumulation of ROS.

We measured PAL, POD, SOD, and CAT activities in leaf samples leaves from WT and CmWRKY15-1-OE lines of ‘Jinba.’ The activities of all four enzymes showed an upward trend and then a decrease at 72 h both in WT and CmWRKY15-1-OE lines following inoculation with P. horiana. The activities of these enzymes were significantly higher in CmWRKY15-1-OE lines than in the WT. SOD and POD activities reached their highest levels at 48 h in both WT and CmWRKY15-1-OE lines. CAT and PAL activities peaked at 24 h in WT and CmWRKY15-1-OE lines (Figure 2). Our results showed that the overexpression of CmWRKY15-1 helps to activate defense systems to respond to infection with P. horiana.

Approximately 100.8 GB of clean reads were acquired from 12 samples. The average Q20 and Q30 values of the clean reads were 95.56% and 89.94%, respectively, showing that all clean reads were of high quality (Supplementary Table 1). We successfully annotated 272,669 unigenes via comparisons with databases, namely, KOG, NR, Swiss-Prot, GO, KEGG, Pfam, and TrEMBL. Results are shown in Supplementary Table 2. We performed a statistical analysis of the unigenes with homologs in the NR database to determine the species distribution of homologous sequences. The top three species with the most annotated homologous genes were sunflower (Helianthus annuus), lettuce (Lactuca sativa), and cardoon (Cynara cardunculus), accounting for 35.56, 19.09, and 16.63% of the total, respectively (Supplementary Figure 1). We also annotated the 53,463 unigenes in the KOG database. The highest number of genes (11,233) belonged to the R category: General Function Prediction Only (Supplementary Figure 2).

To achieve a comprehensive insight into the functions and mechanisms of DEGs in ‘C029’ after inoculation with P. horiana, we analyzed the DEGs among three comparison pairs: W0 vs. WJ, R0 vs. RJ, and WJ vs. RJ. The expression levels of 5,449 and 6,128 genes significantly differed in W0 vs. WJ and in R0 vs. RJ, respectively. In addition, we also found 1,684 DEGs in WJ vs. RJ. To explore the functions of these DEGs, they were annotated using GOseq. The annotations were grouped into three categories: cellular components, molecular functions, and biological processes. In the cellular components category, the “cell” and “cell part” subclasses were the most enriched in W0 vs. WJ, R0 vs. RJ, and WJ vs. RJ. In the molecular functions category, the most highly enriched GO terms were “catalytic activity” and “binding” in W0 vs. WJ, R0 vs. RJ, and WJ vs. RJ (Supplementary Figure 3). Analysis of the biological process category showed that the “cellular process” and “metabolic process” subclasses were the most abundant in W0 vs. WJ, R0 vs. RJ, and WJ vs. RJ.

We used KEGG enrichment pathway analysis to explore the molecular mechanisms of these DEGs. Most of the DEGs in W0 vs. WJ, R0 vs. RJ, and WJ vs. RJ belonged to the KEGG pathways “metabolic pathway” and “biosynthesis of secondary metabolites” (Supplementary Figure 5). A Venn diagram created to visualize the number of DEGs between W0 vs. WJ, R0 vs. RJ, and WJ vs. RJ revealed an overlap of 65 DEGs (Supplementary Figure 4).

To validate the reliability of DEGs identified by RNA-seq analysis, we analyzed the expression levels of seven candidate genes in the SA signaling pathway from the 65 DEGs between W0 vs. WJ, R0 vs. RJ, and WJ vs. RJ via RT-qPCR (Supplementary Table 3). These genes included two SA synthetic signaling pathway genes (CmPAL1 and CmPAL2), three SA-response-related genes (CmPR1, CmPR5, and CmNPR1), and two TGA transcription factor genes (CmTGA1 and CmTGA3). Our results showed that the RT-qPCR data were closely related to the RNA-seq data. In the silenced lines, expression of DEGs was inhibited, as their expression levels were always lower in these lines than in the WT (Figure 3). One DEG, CmNPR1, is an important regulator and a core point in the downstream signal transduction network of plant disease resistance. Expression of this gene significantly increased in both the WT and CmWRKY15-1-RNAi lines after treatment with P. horiana.

We cloned the CDS of CmNPR1 from ‘C029,’ showing that it is 1,670 bp long and encodes a 558-amino-acid protein (Figure 4A). Analysis of the CmNPR1 protein sequence revealed that it has 97, 77, and 78% sequence similarity to the NPR protein sequences from sweet wormwood (Artemisia annua), H. annuus, and C. cardunculus, respectively (Figure 4B), confirming that this gene is CmNPR1. According to the predicted results, the CmNPR1 protein N-terminal region includes an evident BTB/POZ domain between amino acids 137 and 552; this domain is a potential target of the Cullin3-E3 ligase ubiquitin degradation pathway. An ANK domain was present in the middle region of this protein at amino acids 802–1,059. This domain is responsible for mediating the interaction between NPR1 and WRKY or TGA. The C-terminal NPR1/NIM1 region (between amino acids 1,102 and 1,632) and the N-terminal BTB/POZ domain participate in the binding of NPR1 to components of the SA signaling pathway.

Analysis using ExPASy revealed that NPR1 protein primarily comprises leucine, serine, valine, glutamic acid, aspartic acid, and lysine residues, which account for 12.2, 9.1, 8.2, 7.5, 7.2, and 7.0% of the amino acid residues, respectively. The molecular formula of the protein is C₂₇₇₄H₄₅₀₁N₇₇₃O₈₄₇S₂₉. Based on an analysis of the predicted protein structure, we submitted the CmNPR1 protein sequences to SOPMA to predict its secondary structure. Of the 310 alpha-helices present in CmNPR1, 176 were part of a random coil, 52 were part of an extension strand, and 21 were part of a beta turn, accounting for 55.46, 31.48, 9.30, and 3.76% of the protein, respectively (Figure 4C). We predicted its tertiary structure using the Phyre² database (Figure 4D). In addition, a comparison of CmNPR1 with other proteins from various species indicated that it is most similar to AaNPR from A. annua and RcNPR1 from the castor bean (Ricinus communis) (Figure 4E).

To verify the function of CmNPR1, we obtained CmNPR1-silenced lines (TRV-NPR1) of ‘C029’ generated by VIGS. We then performed RT-qPCR to identify any changes in the expression of these downstream PR genes associated with plant defense in WT and TRV-NPR1 lines. Three genes, CmPR1, CmPR2, and CmPR10, showed lower expression in the silenced plants than in the WT controls at 0, 24, 48, and 72 h. Transcript levels of all three genes initially increased, peaked at 48 h, and then decreased in both WT and TRV-NPR1 lines (Figure 5). The results showed that CmNPR1 can respond to infection with P. horiana and induce the expression of PR genes.

We cloned the 1,685-bp promoter of CmNPR1 by hiTAIL-PCR (Figure 6A). We analyzed cis-acting elements in the CmNPR1 promoter using the Plant CARE and PLACE databases. Numerous motifs were discovered in the CmNPR1 promoter, including core cis-acting elements: 13 transcription start site TATA-boxes; 19 CAAT-boxes; one MYB element (GA-responsive); one SARE sequence; seven as-1 (SA-responsive) elements; one TGACG motif [methyljasmonate (MeJA)-responsive]; one AT-rich sequence (stress-responsive); three GT1 motifs (light responsive); and five W-boxes (Figure 6B). The W-box is a binding site for WRKY transcription factors, suggesting that the promoter of CmNPR1 may be regulated by the transcription factor CmWRKY15-1. Based on analysis of cis-acting elements, using ‘C029’ DNA as a template, we constructed three expression vectors with the GUS gene driven by three promoter sequences and obtained by PCR (in place of the Cauliflower mosaic virus 35S promoter): the full-length CmNPR1 promoter sequence and two truncated sequences (Figure 6C). To explore the mechanism of CmNPR1 promoter fragments on P. horiana response, we transiently transformed three promoter fragments into ‘C029.’ We sampled leaves at different time points after inoculation with P. horiana and evaluated GUS activity. Little GUS staining was observed in any of the transgenic leaves at 0 h. Compared with the mock control, the three transgenic leaves (harboring pNPR1.1:GUS, pNPR1.2:GUS, and pNPR1.3:GUS) gradually turned blue with increasing time after inoculation, with the deepest staining observed at 48 h, after which staining decreased again. GUS activity driven by the CmNPR1 promoter was significantly upregulated after treatment with P. horiana, and GUS activity driven by the CmNPR1.1 promoter (pNPR1.1: GUS) was significantly higher than that of pNPR1.2:GUS and pNPR1.3:GUS (Figure 6D). We propose that CmNPR1 has a pathogen-inducible promoter and can be induced by P. horiana. In addition, the change in CmNPR1 expression levels in ‘C029’ leaves after treatment with P. horiana is closely related to its change in promoter activity.

The interaction between the promoter of full-length CmNPR1 and the CmWRKY15-1 protein was verified with a Y1H assay. We transferred the bait (pNPR1-AbAi.) and prey (pGAD-CmWRKY15-1) vectors into Y1HGold yeast cells; pGADT7 and p53 were used as a positive control. The strain was able to grow on plates of SD/-Leu medium supplemented with 800 ng/ml AbA when pNPR1-AbAi was coexpressed with pGAD-CmWRKY15-1 in yeast (Figure 7B). The CmNPR1 promoter was divided into two fragments, PA1 and PA2, and the promoter region PA1 contains W-box elements (Figure 7A). We transformed Y1HGold cells with the pPA1-AbAi or pPA2-AbAi bait vector and the pGAD-CmWRKY15-1 prey vector and plated the cells on SD/-Leu medium containing 800 ng/ml AbA. Y1H assays showed that CmWRKY15-1 interacted with fragment PA1 but not PA2 (Figure 7A). Therefore, CmWRKY15-1 can bind to the W-box elements in the promoter of its target gene CmNPR1. These results show that CmWRKY15-1 has DNA binding activity and specifically binds to the CmNPR1 promoter.

Finally, to examine the interaction between CmWRKY15-1 and CmNPR1, we performed a BiFC analysis. We cloned the CDS of CmWRKY15-1 and CmNPR1 into vectors containing the N-terminal of p2YN and C-terminal regions of p2YC and then cotransfected N. benthamiana leaves with these constructs (p2YC-CmNPR1 and p2YN-CmWRKY15-1). The p2YC vector coexpressed with p2YN-CmWRKY15-1 vectors served as the negative controls. YFP signal was observed only in N. benthamiana leaves coexpressing p2YC-CmNPR1 and p2YN-CmWRKY15-1, but not in cells coexpressing the negative control constructs (Figure 7C). This indicates that CmWRKY15-1 and CmNPR1 interact in vivo.