Morus notabilis leaves were obtained from mulberry ring rot disease-affected M. notabilis gardens in Qingshen County, Meishan City, Sichuan Province. A total of 10 diseased M. notabilis were randomly selected at the time of harvesting, and 10 mature leaves with mulberry ring rot disease spots were taken from each M. notabilis. The surface of the leaves was sterilized with 95% ethanol, then rinsed with aquae sterilization and dried with sterile filter paper. The spatial location on the leaves was divided (Figure 1) so that the area covered by all the whorls in and within the second circle from the outside to the inside was the heart zone (HT zone), which usually showed a dead yellow color. The outermost spot whorl cirri (width a) and its connected outer cirri (1.5a) were the limbic zone (LT zone), which showed both brown and green cirri. The area outside the LT zone of the lesion (3a) is the non-transmissible zone (NT), which usually shows a green color. The three zones corresponded to the pre, mid, and post periods of infection by the mulberry ring rot disease pathogen. The undeveloped, limb and core zones of each leaf were cut separately with a sterile scalpel, taking only one spot from each area. The samples were divided into three different samples according to their positions and named HT, LT, and NT. Samples were taken at the same positions of healthy leaves as control (CK). All samples were immediately frozen in liquid nitrogen and later stored in a −80°C refrigerator.

The sample was ground with liquid nitrogen, weighed 0.1 g, placed in a conical flask, precision added 75% methanol by volume 50 ml, weighed, ultrasonicated (power 250 W, frequency 40 kHz) for 20 min, cooled, then weighed, and made up for the loss of weight with 75% methanol by volume, and filtered through 0.22 μm microporous filtration membrane, and the filtrate was injected into the UHPLC-MS/MS system for analysis (Zhou, 2017). The sample was injected into the column with a 17 min linear gradient and a flow rate of 0.2 ml/min. The chromatographic gradient elution program is shown in Supplementary Table 1. The Orbitrap Q Exactive™ HF-X mass spectrometer was operated in positive and negative polarity mode, with the Spray Voltage set to 3.2 kV, the Sheath gas flow rate set to 40 arb, the Capillary Temp set to 320°C and the Auxiliary Gasflow rate was set to 10 arb.

Each sample had 1 μg of total RNA and was used as an RNA sequencing sample preparation as input material for generating sequencing libraries. A total of 12 RNA-Seq libraries were constructed with four treatments of HT, LT, NT, and CK, with three replicates for each combination. RNA integrity and total amount were assessed using the RNA Nano 6000 Assay Kit from the Bioanalyzer 2100 System (Agilent Technologies, CA, USA). Library preparations were sequenced on the Illumina Novaseq platform and 150 bp paired-end reads were generated. Raw reads in fastq format were processed through in-house perl scripts to obtain high-quality clean reads, and all subsequent analyses were based on the clean reads. The data quality of the sample sequencing outputs was assessed by calculating Q20, Q30, and GC. The clean reads after quality control were compared to the reference genome. The mulberry reference genome and gene model annotation files were downloaded directly from the genome website.¹ The index of the reference genome was built using Hisat2 (version 2.1.0),² and the locus information of the reads was obtained by using Himat2 to compare the paired-end clean reads to the reference genome (Meng, 2019). All genes were annotated for gene function based on GO (Gene Ontology) and KEGG (Kyoto Encyclopedia of Genes and Genomes) databases.

The FPKM of each gene was calculated based on the length of the gene and the reads mapped to that gene, and the relative expression of the transcripts was estimated based on the FPKM (Bray et al., 2016). Two-by-two comparisons were performed from the RNA-Seq data, i.e., HT vs. CK, LT vs. CK, and NT vs. CK. Differential expression analysis was performed for each sample based on the FPKM values of the genes in each group using the DESeq2 R package (version 1.20.0) with an adjusted p-value < 0.05 and | log2(Foldchange)| ≥ 0 as criteria for screening differential genes (DEGs), which were used for GO and KEGG pathway enrichment analysis.

Since the isolation and purification of Spondylocladium mori (Sawada) are difficult, the test organism used in the in vitro inhibition test was Corynespora cassiicola, which was isolated and identified from mulberry whorl spots, and the test strain was provided by Sichuan Academy of Agricultural Sciences. The mycelia were inoculated into potato dextrose agar (PDA) medium for 3 days with inoculation needles and left to be used. Before the experiment was carried out, the pre-test proved that the organic solvent dimethyl sulfoxide (DMSO) would not inhibit the bacteria. The bacteriostatic substance was then configured into a 3 mg/ml solution based on molecular polarity. Autoclaved 30 ml of PDA was added to the solution and poured into petri dishes (90 mm × 20 mm). The solution was cooled and inoculated with C. cassiicola grown for 72 h. Colonies were measured and photographed at 24 h intervals. Three biological replicates were performed for each treatment. In this study, 11 metabolites for in vitro antibacterial test were screened based on the enrichment of DAMs under the pathway by combining data from transcriptome and metabolome assays (Supplementary Table 2).

Data from principal component analysis (PCA) were transformed with metaX software and then plotted with Origin 2022. Cluster heatmaps were plotted using the R package Pheatmap and metabolite data were normalized with z-score. The R package ggplot2 was used to plot bubble plots.

In this experiment, Illumina sequencing was used to transcriptionally sequence 12 M. notabilis leaf samples, and the quality control of raw reads was performed through raw data filtering, Phred score calculation (Q20 and Q30), sequencing error rate checking, and GC content distribution checking to obtain clean reads for use in subsequent analysis (Table 1). A total of 82.3G clean sequences were obtained, and the average clean sequences of CK, HT, LT, and NT were 7.05G, 6.91G, 6.69G, and 6.79G, respectively, with an error rate of 0.03%, an average of 97.65% for Q20, 93.37% for Q30, and an average of 46.19% for GC content, which indicated that the data were reliable and could be used for subsequent analysis. In addition, the relatively low Mapping rate of HT samples might be caused by too much variation of the samples relative to the reference genome at the late stage of infestation.

A total of 16,293 DEGs were detected in the six comparison groups. A total of 3,921 DEGs were detected in HT vs. CK (2,405 up and 1,516 down); a total of 4,274 DEGs were detected in LT vs. CK (2,547 up and 1,727 down); and a total of 2,626 DEGs were detected in NT vs. CK (1,619 up and 1,007 down) (Figure 2). This shows that significant changes occurred in the middle and late stages of infection compared with the pre-infection stage, which may be related to the increase of disease resistance pathways. The differential genes in HT vs. CK were slightly lower than those in LT vs. CK, which may be caused by the dieback condition in the center of the leaf spots of M. notabilis in the late stage of infection.

To further understand the function of DEGs, we performed GO functional annotation of DEGs in different partitions. We presented the 30 most significant terms in bar charts according to the three categories of biological processes (BP), cellular components (CC), and molecular functions (MF), as well as the up and down of DEGs (Supplementary Figure 1). The results showed that HT vs. CK DEGs were mainly involved in biological processes such as “response to biotic stimulus,” “defense response,” “carbohydrate derivative catabolic process,” and “photosynthesis,” etc. The cellular components are mainly enriched in processes such as “photosystem,” “photosynthetic membrane,” “thylakoid,” “thylakoid part,” etc. Mainly involved in molecular functions such as “heme binding,” “tetrapyrrole binding,” “oxidoreductase activity,” “iron ion binding,” etc. The DEGs screened by LT vs. CK were mainly involved in the following biological processes: “defense response,” “response to biotic stimulus,” “response to stress,” “carbohydrate derivative catabolic process.” It is mainly involved in “cell periphery,” “cell wall,” “external encapsulating structure,” and “tethering complex.” Mainly involved in molecular functions such as “tetrapyrrole binding,” “heme binding,” “enzyme inhibitor activity,” and “oxidoreductase activity” and so on. The DEGs screened by NT vs. CK were mainly involved in biological processes such as “multi-organism process,” “cell recognition,” “pollination,” and “pollen-pistil interaction.” It is mainly involved in “DNA packaging complex,” “protein-DNA complex,” “thylakoid,” and “thylakoid part” in cellular components. In molecular functions, it is mainly involved in “transcription regulator activity,” “DNA binding transcription factor activity,” “heme binding,” “tetrapyrrole binding,” etc. This suggests that the defense response of M. notabilis is not strong in the early stage of the disease, but in the middle and late stage of the infection, it is a direct defense through “response to biotic stimulus,” “defense response,” and “response to stress,” and at the same time, the resistance is increased through the regulation of metabolism level in the plant.

Kyoto Encyclopedia of Genes and Genomes pathway analysis can further understand the biological functions and interactions of related genes (Cheng et al., 2023). The results showed (Figure 3) that HT vs. CK differential genes were significantly enriched in “Galactose metabolism,” “alpha-Linolenic acid metabolism,” “Photosynthesis,” and “Glucosinolate biosynthesis” pathways. In addition, many DEGs are listed in “MAPK signaling pathway – plant,” “Tyrosine metabolism,” “Tryptophan metabolism,” “Plant hormone signal transduction,” “Cysteine and methionine metabolism,” “Flavonoid biosynthesis” and other pathways. The LT vs. CK DEGs mainly enriched in “Phenylpropanoid biosynthesis,” “Tryptophan metabolism,” “Plant hormone signal transduction,” “Cysteine and methionine,” and “Flavonoid biosynthesis” pathways. “Phenylpropanoid biosynthesis,” “Carotenoid biosynthesis,” “MAPK signaling pathway – plant,” “alpha-Linolenic acid metabolism,” “Plant hormone signal transduction,” “Glucosinolate biosynthesis” also had a large enrichment of DEGs. NT vs. CK is mainly enriched in pathways such as “Photosynthesis,” “Photosynthesis-antenna proteins,” “Carbon fixation in photosynthetic organisms,” “Plant-pathogen interaction,” “Phenylpropanoid biosynthesis,” “Flavonoid biosynthesis,” “alpha-Linolenic acid metabolism,” and “Pentose phosphate pathway.” The enrichment of the “plant-pathogen interaction” pathway indicates that plants have begun to defend against pathogenic bacteria in the early stage of infection. In response to the infection of mulberry ring rot disease, M. notabilis mainly involves phenylpropanoid biosynthesis, flavonoid biosynthesis, plant hormone signal transduction, and other pathways, which may be related to the resistance to pathogenic bacteria.

An overall clustered heat map analysis was first performed on the samples separately in the positive and negative ion modes, which revealed significant intergroup variation that could be used for the next analyses (Figure 4). PCA was then performed on each sample to determine the overall distribution trend among the different samples. As shown in Supplementary Figure 2, the contribution of the first principal component (PC1) in the positive ion mode was 59.03%, 58.05%, and 32.23% for the three regions, and the contribution of PC1 in the negative ion mode was 59.24%, 56.49%, and 39.50%, respectively; and the contribution of the second principal component (PC2) in the positive ion mode was 14.96% for the three regions, respectively, 15.86% and 29.37%, respectively, and the contribution of PC2 in the negative ion mode was 14.58%, 17.05%, and 24.87%, respectively; at the same time, all samples were significantly separated. The relationship between metabolite expression and sample species was modeled using partial least squares regression Supplementary Figure 3. The PLS-DA correlation models R2Y (cum) and Q2Y (cum) ranged from 0.99 to 1.00 and 0.59 to 0.9, respectively, among the different comparison groups under the positive and negative ion model. R2Y was larger than Q2Y, which was a better construction of the model.

Screening for differential metabolites yielded a total of 625 DAMs in positive ion mode (Figure 5A), with 493 (187 up, 306 down), 257 (84 up, 173 down), and 87 (25 up, 62 down) obtained from the CK vs. HT, CK vs. LT, and CK vs. NT comparator groups, respectively. The negative ion mode screened 331 DAMs, and the three comparison groups of CK vs. HT, CK vs. LT, and CK vs. NT screened 234 (66 up, 168 down), 123 (37 up, 86 down), and 71 (13 up, 58 down) species, respectively (Figure 5B). As can be seen from the Wayne diagram, there were 17 differential metabolites common to the three regions in the positive ion mode, 303 species accumulated only in the HT region, 77 species accumulated only in the LT region, and 50 species accumulated only in the NT region. There were four differential metabolites common to the three regions in the negative ion mode, 149 species accumulated only in the HT region, 43 species accumulated only in the LT region, and 46 species accumulated only in the NT region. In both positive and negative ion modes, DAMs were most abundant in CK vs. HT, followed by CK vs. LT, which was consistent with the transcriptome results. It suggests that M. notabilis has the strongest defense response against pathogens in the middle and late stages of infection.

Kyoto Encyclopedia of Genes and Genomes enrichment analysis of DAMs (Figure 6) showed that HT vs. CK pos enriched a total of 43 metabolite pathways, including “Metabolic pathways,” “Biosynthesis of secondary metabolites,” “Flavonoid biosynthesis,” “Phenylalanine metabolism,” “Flavone and flavonol biosynthesis,” “Phenylpropanoid biosynthesis” and other pathways were more abundant. LT vs. CK pos enriched a total of 33 metabolite pathways, “Biosynthesis of secondary metabolites,” “Biosynthesis of amino acids,” “Flavonoid biosynthesis,” “Tryptophan metabolism,” “ABC transporters,” “Tropane, piperidine, and pyridine alkaloid biosynthesis” were enriched in large quantities. NT vs. CK pos had the least enrichment pathway, and a total of 12 metabolite pathways were enriched, including “Tryptophan metabolism,” “Biosynthesis of amino acids,” “Phenylpropanoid biosynthesis,” “Flavone and flavonol biosynthesis,” “Purine metabolism,” and other pathways. In negative ion mode HT vs. CK neg, LT vs. CK neg, and NT vs. CK neg were enriched in 35, 24, and 17 pathways, respectively. The flavonoid biosynthesis pathway was significantly enriched in both the mid- and late stages of infection, which is consistent with the results of transcriptomics analysis.

Transcriptome KEGG enrichment revealed that 36 DEGs were injected into the flavonoid biosynthesis pathway (Figure 7), including 4 PAL genes, 3 4CL genes, 1 C4H gene, 16 chalcone synthase (CHS) genes, 1 CHI gene, 1 flavanone-3-hydroxylase (F3H) gene, 4 DFR genes, 2 flavonol synthase (FLS) genes, 2 HCT genes, and 2 F3′H genes. The expression of most genes was increased compared to CK, especially CHS, F3H, and FLS, and they increased more in HT and LT. HCT and C4H genes, however, showed down. Among the DAMs, flavonoids included quercetin, kaempferol, eriodictyol, naringenin, hesperetin, apigenin-8-C-glucoside, myricetin, and their corresponding glycosides.

The DEGs and DAMs associated with the flavonoid biosynthesis pathway were mapped to the KEGG pathway (Figure 8). Based on the expression of metabolites it can be seen that most of the metabolites showed upregulation compared to CK. L-phenylalanine, p-coumaric acid, naringenin, hesperetin, dihydrokaempferol, kaempferol-3-O-β-glucopyranosyl-7-O-α-rhamnopyranoside, apigenin, vitexin, pelargonin chloride, quercetin-O-glucoside, quercetin 3-alpha-L-arabinofuranoside (Avicularin), C-hexosyl-luteolin O-feruloylhexoside, C-pentosyl-luteolin-C-hexoside, luteolin, (+)-catechin, (−)-catechin gallate, gallocatechin gallate, procyanidin A3, procyanidin B1, procyanidin B2, dihydromyricetin, delphinidin 3-sophoroside-5-rhamnoside, (−)-epigallocatechin, myricetin, myricetin 3-O-galactoside, and epicatechin, metabolites that showed up in LT or HT. Shikimic acid, neohesperidin, kaempferol, apigeninidin 5-glucoside, apigenin C-pentoside, eriodictyol, cyanidin-3-O-glucoside, cyanidin-3-O-glucoside chloride, cyanidin O-rutinoside, and quercetin showed up in NT compared to CK. Only naringin dihydrochalcone, 1-caffeoylquinic acid, prunin, apigenin C-glucoside, hesperidin, naringin, di-C, C-pentosyl-apigenin, eriodictyol O-malonylhexoside, and eriodictyol-7-O-glucoside, all of these metabolites showed down compared to CK. Most of the metabolites were significantly upregulated in the middle and late stages of infection, and some of them were also significantly upregulated in the early stages of infection, suggesting that flavonoid metabolites play a role in the defense against mulberry ring rot disease.

Metabolomics analysis revealed significant accumulation of M. notabilis flavonoid compounds after infection by the mulberry ring rot disease pathogen, therefore we selected seven metabolites related to flavonoid biosynthesis pathway, L-phenylalanine, proanthocyanidins, (−)-epicatechin (EC), (+)-catechin, quercetin, naringenin, and kaempferol for in vitro antimicrobial test in C. cassiicola. In addition literature suggests that phenylpropanoid biosynthesis and plant hormone signal transduction also play an important role in plant disease resistance (Bao et al., 2020), also in the present study the enrichment of both pathways was found, therefore L-lysine and caffeic acid associated with phenylpropanoid biosynthesis, and salicylic acid and methyl jasmonate associated with plant hormone signal transduction were selected to be used for the in vitro antibacterial test along with the phenylpropanoid biosynthesis. After 48 h of inoculation, all nine metabolites showed strong inhibitory effects on C. cassiicola, except for two amino acids, L-phenylalanine, and L-lysine, which had low antimicrobial properties (Figure 9). Significance analysis of colony diameter at 24, 48, and 72 h showed (Figure 10) that kaempferol and (+)-catechin had the best inhibitory effect at 24 h of inoculation, and proanthocyanidins had the worst inhibitory effect. L-lysine and caffeic acid showed strong inhibitory effect at 24 h, and weaker inhibitory effect at invasive naringenin and salicylic acid showed the opposite inhibitory effect, which was significantly greater than most of the other substances at 48 and 72 h, and the inhibitory effect also showed a trend of gradual enhancement; methyl jasmonate showed strong inhibitory effect throughout the whole process; moreover, the inhibitory effect of 24 and 48 h (−)-epicatechin (EC) was the best. The (−)-epicatechin (EC) and (+)-catechin inhibited C. cassiicola better than proanthocyanidins at 24 and 48 h, and proanthocyanidins produced higher inhibition than (−)-epicatechin (EC) and (+)-catechin at 72 h.