Strawberry plants (cv. Florida Beauty), which are susceptible to anthracnose fruit rot disease, were provided by Dr. Vance Whitaker at University of Florida Institute of Food and Agricultural Sciences (UF/IFAS) the Gulf Coast Research and Education Center (GCREC) in Wimauma, Florida, USA. The plants were multiplied from stolons hydroponically under greenhouse conditions at the UF/IFAS Indian River Research and Education Center (IRREC) in Fort Pierce, Florida. When the plants reached the ideal size (approximately one month after planted from stolons), strawberry leaves were freshly harvested and rinsed with sterile water and used in the detached leaf infection assays. The detached fruit infection assay was performed at GCREC with strawberry fruit at 1/3 maturity stage freshly harvested from the field.

Colletotrichum nymphaeae strain 02-179 was provided by Dr. Megan Dewdney from the UF/IFAS Citrus Research and Education Center (CREC) in Lake Alfred. C. nymphaeae strain 02-179 was collected from diseased strawberry fruit tissue in Dover, Florida, USA in 2002. Wang et al. (2019) identified strain 02-179 as C. nymphaeae within the species complex of C. acutatum (Wang et al., 2019).

The leaves were detached and placed in square bioassay dishes (Thermo Scientific, catalog No. 12-565-224). We applied the detached assay rather than the attached assay to achieve control on the environmental conditions within the laboratory, which is supported by previous studies, such as those involving ‘Camarosa’ and ‘Strawberry Festival’ cultivars, by demonstrating a correlation between detached and assays and field trial outcomes (Forcelini et al., 2017). Colletotrichum nymphaeae strain 02-179 was cultured on potato dextrose agar (PDA) plates (Oxoid, catalog No CMO139) under optimal conditions (~25°C and 24 h light). After 7-10 days, fungal spores (conidia) were harvested, the spore concentration was counted via Haemocytometer, and the final concentration adjusted to 1x10⁷ spores/mL. For the wounded leaf assay, 20 μL aliquots of the spore solution were placed as droplets on leaves previously wounded with a 10 μL sterile pipette tip. In the unwounded leaf assay, leaves were inoculated with same concentration and volume of spores without prior wounding. Both sets of leaves were placed over wet autoclaved sterile papers in the bioassay dishes maintained under optimum infection conditions (~25°C and 14 h light) for symptom development. Three to four leaves were collected 1-, 2- and 3-days post inoculation (dpi) for the wounded leaf assay, and on 1-, 2-, 3-, and 5-dpi for the unwounded leaf assay. The 5 dpi samples from the wounded leaves were excluded due to early rotting because of wounding. However, in the unwounded leaves the samples remained fresh and for that reason 5 dpi samples were included. For each leaf, samples were taken with a No. 5 cork borer (11.2 mm diameter) from areas where spore droplets were placed. All leaf samples were immediately ground to a fine powder with liquid nitrogen using a mortar and pestle and stored at -80 °C prior to RNA extraction. The entire assay, including the mycelia control, was conducted with two biological replicates. Mycelia control was harvested from strain 02-179 grown on PDA plates after 7 days and treated in the same manner as the leaf samples. For the fruit infection assay, fresh strawberry fruit at 1/3^rd maturity were collected from the fields at GCREC in Wimauma, Florida. Harvested unripe fruit were sterilized by dipping into 0.3% household bleach solution for 3 minutes and rinsing with deionized water twice. For each time point sample, 12 fruit were placed into an egg carton, and cartons were placed in clean and transparent plastic boxes with lids. Each egg carton was placed in a separate box, and 10 cartoons/boxes were used. In total 120 strawberries (n = 10 x 12 = 120) were used for 4 time points and mock control of both biological replicates. To keep high humidity conditions in the boxes, 100 mL deionized water was added to the box. C. nymphaeae strain 02-179 was grown in PDA under optimum conditions (~25°C and 24 h light), and spores were harvested after 7 days. The final concentration of the spore solution was adjusted to 1x10⁶ spores/mL. An aliquot of 20 μL of the spore concentration was placed as droplets onto the strawberry fruit placed in the egg cartons and within plastic boxes. Then the boxes were incubated at optimum infection conditions (25°C and 14 h light) for the development of disease symptoms. Fruit samples were collected 1-, 2-, 3-, and 5-days post inoculation (dpi). For each collected fruit sample, the area where the droplets were placed was taken with a scalpel. Mycelia were also used as a control for fruit infection assays. Fruit tissues and mycelia samples were immediately ground to a fine powder with liquid nitrogen using mortar and pestle. Frozen powder samples were stored at -80°C prior to RNA extraction.

Total RNA extraction of strawberry infected tissues was performed with RNeasy^® Plant Mini Kit for each time point of the wounded leaf assay and the fruit assay according to the manufacturer’s protocol (Qiagen, Catalog No. 74903). To increase the RNA concentration, during the QIAshredder step, 4-6 tubes of the lysate samples were eluted into the same lilac collection tube. For the unwounded leaf assay, Spectrum™ Plant Total RNA Kit was applied for each time point according to the manufacturer’s protocol (Sigma-Aldrich, Catalog No. STRN50). RNA samples were digested with DNase enzyme following the manufacturer’s protocol to remove traces of DNA (Invitrogen TURBO-DNA-free kit, Catalog No. AM1907). We performed quality control on the total RNA samples by quantifying and measuring RNA purity (absorbance 260/280 ratios) with Nanodrop Lite and Qubit Fluorometer (v2.0). A total of 2 μg of total RNA for each plant tissue (wounded, unwounded and fruit infected samples) and time points were sent to Novogene Company (Sacramento, California) for cDNA library preparation (250-300 bp insert cDNA library size) and RNA-sequencing. Novogene used Illumina NextSeq500 platform for RNA-sequencing and generated paired-end (PE) read data of 150 bp in length. For each sample, we obtained around 6 Gb of read data. The raw read sequence data was deposited in SRA under accession PRJNA1074080.

Colletotrichum nymphaeae SA-01 was used in this study as the reference species for secretome prediction and effector annotation. The C. nymphaeae SA-01 proteome sequences, and functional annotations were downloaded from the Joint Genome Institute (JGI) (https://genome.jgi.doe.gov/portal/Colny1/) and the Ensembl Fungi Repository (https://fungi.ensembl.org/Colletotrichum_nymphaeae_sa_01_gca_001563115/). Secreted proteins from C. nymphaeae SA-01 were identified by filtering the proteome for characteristics such as the presence of signal peptides, absence of transmembrane domains, and absence of Glycosylphosphatidylinositol (GPI) membrane anchor domains. We used SignalP v2.0 and v5.0 software to predict proteins carrying an N-terminal signal peptide leader sequence (http://cbs.dtu.dk/services/SignalP/) (Almagro Armenteros et al., 2019; Nielsen et al., 1997). We filtered out proteins with transmembrane domains using TMHMM v2.0 (https://services.healthtech.dtu.dk/service.php?TMHMM-2.0) and proteins with predicted GPI anchor domains using PredGPI (http://gpcr.biocomp.unibo.it/predgpi/pred.htm) and GPI-Som (http://gpi.unibe.ch/) (Pierleoni et al., 2008). We carry out the prediction of effector types (apoplastic or cytoplasmic) with EffectorP v2.0 (http://effectorp.csiro.au/) (Sperschneider et al., 2018). In addition, we run predictions for the subcellular localization of the secreted proteins with LOCALIZER (https://localizer.csiro.au/) which is based on the transit peptides (associated with chloroplast and mitochondria) and predictions of nuclear localization signals (NLS) (Sperschneider et al., 2017). CAZymes were predicted using CAZYme database on dbCAN2 server (Zhang et al., 2018a).

Raw reads were processed with FASTX-Toolkit for quality check (Hannon, 2010). Colletotrichum nymphaeae SA-01 coding gene sequences (CDS) fasta files containing 14,404 sequences were used as the reference database for alignments. First, an index for C. nymphaeae SA-01 fasta CDS file was created with bowtie2-build function from Bowtie2 v2.2.6 software (Langmead and Salzberg, 2012). RNAseq reads from the three conditions (wounded and unwounded leaves, and fruit samples) were aligned to C. nymphaeae SA-01 using TopHat2 v2.1.1 software (Kim et al., 2013). The number and percentage of aligned reads for each condition are reported in Table 1 . The output of the alignments in BAM (Binary Alignment Map) format was sorted by the read name with samtools-sort function (with –n parameter) from Samtools package (https://sourceforge.net/projects/samtools/files/) for downstream expression quantitation analysis. To count reads in features, we used htseq-count function (stranded parameter=no) from HTSeq v2.2.0 software (Putri et al., 2022). The raw counts from HTSeq were processed with DESeq2 package installed in RStudio to identify differentially expressed genes (DEG) (Love et al., 2014) by performing normalization using the median of gene count ratios relative to geometric mean per gene. Genes with Log2FC fold change (LFC) values greater than or equal to 1 and p-adjusted values less than 0.05 were considered significantly up-regulated. Genes with (LFC) values less than or equal to -1 and p-adjusted values less than 0.05 were considered significantly down-regulated ( Supplementary Table S1 ). We searched for up-regulated genes among the differentially expressed gene list with known relevant functions associated with pathogenicity ( Table 2 ). The expression of some selected genes based on the significant upregulation was further validated with reverse-transcription quantitative polymerase chain reaction RT-qPCR analysis (Supplementary Figure S3). The visualization of significant up-regulated and down-regulated genes was done with Volcano Plots using ggplot2 package in R ( Figure 1 ) (Wilkinson, 2011). Variations between the time points and the control were checked with Principal Component Analysis (PCA) Plots using ggplot2 package in R ( Figure 2 ). Comparisons of upregulated gene list among infected tissues were performed with E-venn (http://www.ehbio.com/test/venn/#/) and OrthoVenn2 tools (https://orthovenn2.bioinfotoolkits.net/home; Figure 3 ). Heatmap analyses were performed by publicly available Morpheus software (https://software.broadinstitute.org/morpheus) with a hierarchical tree analysis using one minus Pearson correlation method ( Supplementary Figure S1 ).

The cDNA synthesis was performed using the SuperScript™ IV kit (Invitrogen, Catalog No.18091050). Total RNA samples were treated with TURBO DNA-free kit prior to cDNA synthesis (Invitrogen, Catalog No. AM1907). We first prepared a PCR tube per sample using: 1 µl 50 µM Oligo d(T), 1 µl 10 mM dNTP mix and 11 µl of template total RNA (250 ng). cDNA synthesis was carried out in a Thermal Cycler (Bio-Rad T100), the mixture was heated at 65°C for 5 minutes, and after that incubated on ice for at least 1 minute. In a second PCR tube, the reverse transcriptase reaction mix was prepared with 1 µl 100 nM DTT, 1 µl Ribonuclease inhibitor, 1 µl SuperScript™ IV Reverse Transcriptase (200 U/µl), and 4µl 5x SSIV Buffer. Then, the first and second tubes were combined, and the reaction mixture was incubated at 50°C for 10 min followed by inactivating the reaction by incubating at 80°C for 10 min. A quality check was performed on the single stranded cDNA samples. Quantification was performed with Nanodrop Lite and the absorbance A260/A280 ratios were measured with Qubit Fluorometer v2.0. RT-qPCR was performed on a StepOnePlus™ Real-time PCR system (Thermo Scientific, Catalog No. 4376600). For C. nymphaeae detection in this study, we designed specific primers for a total of eight C. nymphaeae candidate target genes: CNYM01_08295 (KXH55186) Tannase, CNYM01_13686 (KXH64601) Tannase, CNYM01_02731 (KXH43592) Multicopper oxidase, CNYM01_07521 (KXH38331) Glycosyl hydrolase family 7, CNYM01_06500 (KXH30597) Endochitinase, CNYM01_11346 (KXH38857) Cytochrome P450, CNYM01_03788 (KXH42765) Pectate Lyase, and CNYM01_07125 (KXH45695) NUDIX domain-containing protein. For C. nymphaeae housekeeping reference gene we designed specific primers for the Actin gene CNYM01_04906 (KXH42418) ( Supplementary Table S2 ). The expression levels of the target genes of the fungal pathogen C. nymphaeae in the infection process were normalized to the constitutively expressed Actin reference gene CNYM01_04906 KXH42418) and were calibrated compared with the levels recorded in the control mycelia samples (set as 1) using the 2^–ΔΔCt method (Livak and Schmittgen, 2001) (Supplementary Figure S3). The quantitative qPCR assay was carried out in a final volume of 10 μL, consisting of 1 μL of cDNA template (from a 50 ng/μL stock to a final amount of 5 ng per 10 ul of SYBR green reaction), 0.5 μL of primer forward, and 0.5 μL of primer reverse (from 10 uM stock to a final concentration of 500 nM for each primer), 5 μL 2× Fast SYBR, and 3 μL of nuclease-free water. Three technical replicates were run for each sample.

The gene ontology GO terms from C. nymphaeae strain SA-01 was downloaded from the Joint Genome Institute (JGI) (https://genome.jgi.doe.gov/portal/Colny1/). For the GO enrichment analysis, we used as input the list of up-regulated genes during infection and encoding for secreted proteins and load it into the SEA custom option from AgriGO Gene Ontology Analysis Toolkit http://systemsbiology.cpolar.cn/agriGOv2/c_SEA.php. We ran the GO enrichment analysis for three gene lists associated with the strawberry infected tissues ( Supplementary Table S3 ). The GO enrichment statistical method used was Fisher test with a significance level threshold of 0.05.

Protein sequences of a set of 14 annotated tannase family members were retrieved from C. nymphaeae SA-01 proteome and combined in a multi FASTA format file. Additionally, we obtained protein sequences of two tannase genes, sscle_08g067140 and sscle_10g076570, from Sclerotinia sclerotiorum reported to be up-regulated during infection in different host species (L. angustifolius and B. napus respectively) from JGI https://mycocosm.jgi.doe.gov/Sclsc1/ (Allan et al., 2019). These two sequences were used as references then combined with the 14 tannase genes of C. nymphaeae SA01. The combined sequences underwent alignment using MUSCLE alignment function in MacVector 17.0.10 (Olson, 1994). Subsequently, a phylogenetic tree was constructed using the neighbor-joining method in MacVector. For the prediction of active sites in tannase enzymes, the I-TASSER server (https://zhanggroup.org/I-TASSER/) was employed, taking individual protein sequences as input (Zheng et al., 2021) ( Supplementary Figure S2 ).

Transcriptomic analysis was performed from the RNAs extracted from the wounded and unwounded leaves, and fruit tissues of susceptible strawberry cultivar (Fragaria x ananassa cv. ‘Florida Beauty’) infected with Colletotrichum nymphaeae strain 02-179. RNAs were extracted from mycelia tissues grown in vitro conditions for the pathogen control and from the infected host tissues at 1, 2, 3, and 5 dpi. In wounded leaf samples, the percentage of paired reads mapped to C. nymphaeae SA-01 ranged between 1.35% to 1.65% on average. Although the alignment rate between time points only slightly changed, there was a significant difference between the biological replicates ( Table 1 ). In the unwounded leaf samples, the alignment rate changed from 0.11% in the 1 dpi samples to 0.22% in the 3 dpi samples ( Table 1 ). In the fruit samples, the paired read percentage dramatically increased from 0.02% in the 1 dpi samples to 35.35% in the 5 dpi samples. Moreover, unlike the wounded or unwounded leaf samples, a stable increase between each time point from the first time point (1dpi) to the last time point (5 dpi) was observed ( Table 1 ). For mycelia controls, the alignment rate ranged between 44.49% to 46.12%. The paired reads were kept and applied for downstream analysis.

Out of 14,404 proteins in the C. nymphaeae proteome, we identified 2,034 (14.12%) with secretion signals predicted with SignalPv2.0 and SignalPv5.0. After removing the proteins with predicted transmembrane domains (TM) and glycosylphosphatidylinositol (GPI) anchors, we identified a total of 1,517 (10.5%) secreted proteins ( Supplementary Table S1 ). It is known that hemibiotrophic pathogens often possess a significant portion of secreted proteins in their proteomes, typically exceeding 10%. Colletotrichum spp. proteomes have been reported to exhibit a range of 13-16% secreted proteins (Crouch et al., 2014), which aligns with our findings. During the infection, we identified 268, 534 and 555 upregulated genes encoding for secreted proteins in the wounded leaves, unwounded leaves and fruits, respectively ( Supplementary Tables S1 , S4 ). Only 184 upregulated genes encoding for secreted proteins were commonly upregulated in wounded leaves, unwounded leaves and fruit samples ( Supplementary Tables S1 , S4 ). Out of these 1,517 secreted proteins, we found a total of 515 predicted as effectors with EffectorPv2.0. Out of the 515, 62 secreted effectors were commonly upregulated in wounded leaves, unwounded leaves and fruit samples ( Supplementary Tables S1 , S4 ).

For the wounded leaf samples, downregulated and upregulated gene numbers and percentage (in order) were 315 (2.19%) and 347 genes (2.41%) at 1 dpi, 407 (2.83%) and 561 genes (3.89%) at 2 dpi, and 584 (4.05%) and 716 genes (4.97%) at 3 dpi ( Figure 1 ; Supplementary Table S4 ). In the unwounded leaf samples, downregulated and upregulated gene numbers and percentage (in order) were 481 (3.34%) and 796 genes (5.53%) at 1 dpi, 415 (2.88%) and 774 genes (5.37%) at 2 dpi, 804 (5.58%) and 1,338 genes (9.29%) at 3 dpi, and 680 (4.72%) and 1,298 (9.01%) genes at 5 dpi ( Figure 1 ; Supplementary Table S4 ). For the fruit samples, downregulated and upregulated gene numbers and percentage (in order) were 444 (3.08%) and 292 genes (2.03%) at 1 dpi, 350 (2.43%) and 482 genes (3.35%) at 2 dpi, 1,943 (13.49%) and 1,465 genes (10.17%) at 3 dpi, and 2,431 (16.88%) and 2,687 genes (18.65%) at 5 dpi ( Figure 1 ; Supplementary Table S4 ). The percentage of the upregulated and downregulated genes in wounded leaves gradually increased from 1 dpi to 3 dpi. However, in unwounded leaves, the percentage of both upregulated and downregulated genes did not show a steady increase. In the fruit samples, on the other hand, while the percentage of upregulated genes increased from 1 dpi to 5 dpi, the percentage of downregulated genes began to increase continuously only after 2 dpi.

After identifying upregulated genes from the samples at each time point, volcano plots were applied to visualize the data with the adjusted p-value and the Log2 Fold Change (LFC) value. Adjusted p-value were converted into –log₁₀ value and plotted on the y-axis. On the x-axis, LFC values from the expression result are plotted and the threshold value was chosen as (-1) for downregulated genes or (+1) for upregulated genes. The greater expression values are represented closer to the left or right ends; the downregulated genes are shown in green color on the left side of the origin while the upregulated genes are shown in red color on the right side.

To check variation between the time points and the control, and assess quality, Principal Component Analysis (PCA) was applied using the read counts provided by HTSeq. In this analysis, closer distance between the points demonstrates more similar gene expression profiles of the samples (Sun et al., 2019). PCA results showed an accumulation of the points belonging to the samples from each time point of wounded leaves, unwounded leaves and fruits, away from the points belonging to the control group ( Figure 2 ). This separation indicated the big differences between the control and sample groups in each assay. However, when it came to the variation within the time points in each assay, a different profile was observed. In the wounded leaves, despite their proximity, the same time points of the replicates were not necessarily grouped together which showed a variation between the biological replicates, which might be resulted from the wounding process. In the unwounded leaf and fruit samples, lower variation was observed which indicated the reproducibility of the replicates.

Considering all time points and all samples, we identified 4,877 upregulated DEGs, among which 514 (10.5%) core upregulated DEGs were common to fruits, wounded and unwounded leaves ( Figure 3A ). 68.01% of upregulated DEGs (3,317 genes) were specific to the infection of one plant organ with 1,449 upregulated genes during leaf infection only (wounded and unwounded), and 1,868 upregulated genes during fruit infection only, highlighting the prevalence of host colonization programs dedicated to each host organ. However, when the 5 dpi results (which were missing for the wounded leaf) were excluded from the diagram, significant changes were seen. When 5dpi time point was excluded, the number of upregulated DEGs belonging all remaining time points was 3,730 and the number of upregulated DEGs only in the wounded leaves, unwounded leaves, and fruits were 1,036 (27.77%), 2,124 (56.94%), and 1941 (52.04%), respectively. Therefore, compared to unwounded leaves and fruits, wounded leaves had lower number of upregulated genes ( Figure 3B ).

In the upregulated gene profile of all the samples, upregulation was observed for the genes that have previously been shown to be involved in pathogenicity in other Colletotrichum species ( Table 2 ; Supplementary Figure S1 ). These genes included CAZymes involving enzymes with Auxiliary Activities (AAs) and different CWDEs such as pectinases, esterases and glycosyl (glucoside) hydrolases; SM related enzymes including cytochrome P450s, transporters and oxidoreductases; effectors associated with the switch to the necrotrophy such as NUDIX or NEP-1; and LysM proteins.

We identified 784 CAZyme genes in the C. nymphaeae genome through CAZyme prediction. The numbers of the upregulated CAZyme genes were as follows: 40 (5.10%), 65 (8.29%), and 83 (10.59%) on the 1, 2, and 3 dpi of the wounded leaves; 70 (8.93%), 102 (13.01%), 240 (30.61%), and 226 (28.83%) on the 1,2,3, and 5 dpi of the unwounded leaves: and 23 (2.93%), 41 (5.23%), 134 (17.09%) and 289 (36.86%) on 1, 2, 3, and 5 dpi of the fruits, respectively ( Supplementary Table S1 ). The increase in the percentage of upregulated CAZymes from the beginning until the end of the infection in all three samples might be explained by their association with necrotrophy.

Multicopper oxidases are a subfamily of AA enzymes that are more recently integrated into the CAZyme database (Sützl et al., 2018). Two multicopper oxidase genes, KXH30257 and KXH43592, were highly upregulated in all the samples ( Supplementary Table S1 ). KXH43592 exhibited a significant upregulation at the later stages for the wounded and unwounded leaves; respectively, by 7.4-fold at 2 dpi and 8.4-fold at 3 dpi, by 8.5-fold at 3 dpi and 8.3-fold at 5 dpi ( Table 2 ). However, in the fruits it was upregulated at all time points; by 8.6-fold at 1 dpi, 8.7-fold at 2 dpi, 11.5-fold at 3 dpi and 10.3-fold at 5 dpi ( Table 2 ; Supplementary Figure S1 ).

Another important group of CAZymes, CWDEs, were remarkably represented in the DEGs of all the samples. Carboxylesterase gene KXH42185 and glycosyl hydrolase family 7 gene KXH38331 were highly upregulated at the later stages in all samples. For KXH38331, LFC = 2.30 at 1dpi, LFC = 3.48 at 2 dpi, and LFC = 4.01 at 3 dpi in the wounded leaves, LFC = 4.41 at 2 dpi, LFC = 5.17 at 3 dpi and LFC = 5.89 at 5 dpi in the unwounded leaves, and LFC = 4.96 at 1 dpi, LFC = 6.20 at 2 dpi, LFC = 6.55 at 3 dpi and LFC = 3.99 at 5 dpi in the fruits ( Table 2 ; Supplementary Figure S1 ). On the other hand, despite the similarities for the stage being expressed, leaf samples and fruits diverged with the different types of pectinases; pectin esterase gene KXH39301 was upregulated in the wounded (LFC = 10.8 at 2 dpi and LFC = 11.0 at 3 dpi) and unwounded leaves (LFC = 11.7 on 3 dpi and LFC = 12.0 at 5 dpi), meanwhile, in the fruits pectin lyase gene KXH42765 was upregulated on the late stage (LFC = 4.99 at 5 dpi; Table 2 and Supplementary Figure S1 ).

Among the up-regulated genes, one particular tannase gene, KXH55186, had very distinct and overlapping expression between all the samples and different time points: LFC = 8.8 at 1 dpi, LFC = 8.6 at 2 dpi and LFC = 7.1 at 3 dpi in the wounded leaves; LFC = 10.6 at 3 dpi and LFC = 9.1 at 5 dpi in the unwounded leaves; and LFC = 13.2 at 1 dpi, LFC = 12.8 at 2 dpi, LFC = 8.7 at 3 dpi and LFC = 5.0 at 5 dpi in the fruits ( Table 2 ; Supplementary Figure S1 ). Alongside KXH55186, another tannase gene, KXH64601, was found to be significantly differentially expressed and had a different expression pattern compared to the other 12 tannase genes in the same tannase family ( Figure 4 ).

In differential expression analysis, upregulation of several secondary metabolite SM related genes, such as Cytochrome P450s, oxidoreductases and transporters, during the infection were found up-regulated at the later stages.

From the large family of Cytochrome P450 genes with 143 members, one of the most significantly upregulated ones in all the samples was KXH38857. Its expression was increased by 6.40-fold at 2dpi, 5.49-fold at 3 dpi in the wounded leaves, by 4.59-fold at 1 dpi, 5.74-fold at 2 dpi, 5.21-fold at 3 dpi and 6.98-fold at 5 dpi in the unwounded leaves, and by 4.40-fold at 3 dpi and 3.24-fold at 5 dpi in the fruits ( Table 2 ; Supplementary Figure S1 ).

While ABC transporters show a significant upregulation only in the unwounded leaves ( Supplementary Table S1 ), major facilitator superfamily (MFS) transporters were widely upregulated. MFS transporter gene KXH44404 was considerably upregulated: LFC= 6.9 at 3 dpi in the wounded leaves, LFC = 7.4 at 3 dpi and LFC = 9.8 in the unwounded leaves, and LFC = 7.2 at 3 dpi and LFC = 11.1 at 5 dpi in the fruits ( Table 2 ). Among the upregulated SM related genes that overlapped between all the tissues were oxidoreductase genes. NADPH2:quinone reductase gene KXH57053 was one of the oxidoreductase genes that was immensely expressed in all the samples; LFC = 11.0 at 2 dpi and LFC = 11.1 at 3 dpi in the wounded leaves; LFC = 10.8 at 1 dpi, LFC = 9.9 at 3 dpi and LFC = 11.7 at 5 dpi in the unwounded leaves, and LFC = 15.3 at 3 dpi and LFC = 12.5 at 5 dpi, in the fruits ( Table 2 ).

LysM proteins in DEG showed big differences at different stages of infection and diverged between the samples. Intracellular hyphae protein 1, encoded by KXH41657, showed very high upregulation in all the samples: LFC = 8.1 at 1 dpi, LFC = 12.5 at 2 dpi and LFC = 11.9 at 3 dpi in the wounded leaves; LFC = 10.1 at 1 dpi, LFC = 12.9 at 3 dpi and LFC = 12.1 at 5 dpi in the unwounded leaves; and LFC = 15.31 at 2 dpi, LFC = 17.6 at 3 dpi and LFC = 14.3 at 5 dpi, in the fruits ( Table 2 ). Meanwhile, chitinase genes diverged between the samples; endochitinase gene KXH30597 was not differentially expressed in the fruits while chitinase-1 gene KXH64700 was only expressed in the fruit samples at all time points. LFC values of KXH30597 were 5.2 and 3.1 at 1 dpi and 2 dpi, respectively in the wounded leaves, and 5.8 and 3.7 at 1 dpi and 2 dpi, respectively, in the unwounded leaves. LFC values of KXH64700 were 9.1 at 1 dpi, 8.1 at 2 dpi, 5.5 at 3 dpi and 3.5 at 5 dpi in the fruits ( Table 2 ).

Necrotrophic effectors in DEG like necrosis-inducing proteins and NUDIX domain-containing proteins were generally seen in the later stages of the infection. The most significantly expressed genes of these groups were, respectively, KXH60658 and KXH45695. In the wounded leaves, KXH60658 had increased expression at 2 dpi (LFC = 9.6 and at 3 dpi (LFC = 11.0); in the unwounded leaves the upregulation started at 5 dpi (LFC = 11.8), and in the fruits, it started at 3 dpi (LFC = 14.0) and continued at 5 dpi (LFC = 11.5; Table 2 ). KXH45695 had an increased expression at 3 dpi (LFC = 4.4) in the wounded leaves, and at 3 dpi (LFC = 7.0) and 5 dpi (LFC = 6.2) in the unwounded leaves. On the other hand, in the fruits general expressions including all time points were seen; LFC = 7.1 at 1 dpi, LFC = 6.5 at 2 dpi, LFC = 8.5 at 3 dpi and LFC = 4.4 at 5 dpi ( Table 2 ).

I-TASSER generated structure predictions for the protein sequences of two differentially expressed tannase genes (KXH64601 and KXH55186). These predictions were utilized to analyze the active site and ligand binding site residues of the enzymes. Serine and Histidine residues of the catalytic triad, Ser-Asp-His, associated with tannase enzymes identified in the ligand binding site of KXH64601 ( Figure 5A ) and KXH55186 ( Figure 5B ). In the case of the active site belonging to KXH64601, all three residues of the catalytic triad were predicted within the active site, but no prediction was obtained for the active site of KXH55186 ( Supplementary Figure S2 ).

Multiple sequence alignment results showed the presence of the conserved pentapeptide motif “GXSXG”, where the serine residue of the catalytic triad is located, in all tannase sequences ( Figure 6 ). However, the constructed phylogenetic tree, based on this alignment, formed a clade exclusively including tannase enzymes that were differentially expressed during infections. Moreover, the sequences in this clade had a different pattern in their pentapeptide motifs, having either “GCSDG” or “GCSEG” phenotypes while most of the remaining sequences had “GCSTG” phenotype.

Validation of RNAseq gene expression results via RT-qPCR has been widely used in various transcriptomic analyses to confirm technical reproducibility. RNAseq analysis has a superior sensitivity than previous transcriptomic analysis techniques and does not necessarily apply probes which can cause bias (Wang et al., 2009). Nonetheless, RT-qPCR analysis can still be useful to validate RNAseq results as they can perfectly correlate. In our results, we were able to confirm the expression of two tannase genes (KXH64601 and KXH55186) ( Figure 7 ), and other candidate virulence genes encoding for CAZYmes including multicopper oxidase (KXH43592), glycosyl hydrolase family 7 (KXH38331), endochitinase (KXH30597) and pectate lyase (KXH42765), an SM-related enzyme cytochrome P450 (KXH38857), in wounded and unwounded leaves, and fruits ( Supplementary Figure S3 ).

Gene ontology GO enrichment analysis was performed to determine the functional roles of in planta expressed genes in three groups: biological processes, molecular function, and cellular component ( Figure 8 ). Total numbers of enriched GO terms in C. nymphaeae during infection in the wounded leaves, unwounded leaves and fruit were, respectively, 33, 62 and 61 ( Supplementary Table S3 ). The most enriched GO term was in the molecular function domain in the wounded leaves, unwounded leaves and fruits. However, while the hydrolase activity was the most enriched in the wounded leaves, it was the second most enriched GO term in the unwounded leaves and fruits. The most enriched GO term was the catalytic activity for both unwounded leaves and fruits.