RNA-Sequencing in Resistant (QL3) and Susceptible (Theis) Sorghum Cultivars Inoculated With Johnsongrass Isolates of Colletotrichum sublineola

Gene expression was analyzed at 0- and 24-h post-inoculation of two inbred sorghum cultivars known to differ in response to inoculation with Colletotrichum sublineola, the fungal pathogen that causes anthracnose. QL3 is reported to have quantitative resistance, while Theis is susceptible to most pathotypes of the pathogen; RNASeq identified over 3,000 specific genes in both cultivars as showing significant changes in expression following inoculation; in all but one gene, the changes in QL3 and Thies were in the same direction. Many other genes showed significant changes in only one of the two cultivars. Overall, more genes were downregulated than upregulated. Differences in changes in expression levels of a few genes suggested potential roles for the difference in disease response between QL3 and Theis, but did not identify known resistance genes. Gene ontology (GO) and pathway enrichment analysis identified upregulation of 23 transcription factor encoding genes as well as genes involved in the production of secondary metabolites, which are part of a typical host defense reaction.


INTRODUCTION
Sorghum [Sorghum bicolor (L.) Moench] is a crop grown for both grain and stover. Water use efficiency makes sorghum especially important in drier climates. On the other hand, the close relative of sorghum, Johnsongrass [Sorghum halepense (L.) Pers.], is considered one of the most noxious weeds and has adapted to colder climates is now found on a global scale (Fletcher et al., 2020). Estimated to have evolved separately for at least years (Paterson et al., 2020), Johnsongrass, a tetraploid that originated from a cross between the diploid species S. bicolor and S. propinquum, can cross hybridize with S. bicolor, especially if the sorghum is male sterile (Ohadi et al., 2017). The similarity of their genomes raises concern that Johnsongrass may serve as a reservoir for diseases of sorghum. This is tempered somewhat by the possibility that it may also provide sources for disease resistance.
Anthracnose, caused by Colletotrichum sublineola P. Henn, in Kabat and Bubák, is one of the most damaging diseases of sorghum in most the growing areas in the world (Prom et al., 2015). The fungus can infect the leaves, talks, and panicles, although, foliar infection is most common. Symptoms may appear as small circular to elliptical spots or elongated lesions, and as the fungus sporulates, fruiting bodies (acervuli) appear as black spots in the center of the lesions (Thakur and Mathur, 2000). Estimating grain yield losses due to foliar anthracnose infection can often be difficult, but losses as high as 50% have been reported in susceptible cultivars (Prom et al., 2015). The pathogen can be isolated from both Johnsongrass and sorghum. At least so far as sorghum is concerned, there is tremendous variability in C. sublineola, including pathogenicity. This may well be associated with breeding efforts to select resistance, which often have resulted in selecting for single genes that are soon overcome by adaptations in the pathogen. As an example, Thies, a cultivar used in this study, was released in 1974 as a source of anthracnose resistance (Broadhead et al., 1978), but in 2012 screens, it was susceptible to eight of 10 pathotypes identified by tests using 12 differentials (Prom et al., 2012). Resistance in cultivar QL3 to all isolates of the pathogen tested from the US and Puerto Rica has been found to be a quantitative trait (Vermerris et al., 2018). As would be expected, little information is available concerning the variability of pathogenic races among Johnsongrass isolates. Xavier et al. (2018) found only four of 10 isolates from S. halepense caused disease on any of 13 S. bicolor accessions, while Ahn et al. (2020) found that one of three sorghum isolates from sorghum tested via a detached leaf assay led to pathogen reproduction in all 26 Johnsongrass cultivars tested at the eight leaf stage. Ahn et al. (2020) also found that spraying Johnsongrass after the heading stage, but not at the five or eight leaf stage with that same isolate, did give disease, as demonstrated by the formation of acervuli. Personal observations of the consequences of spraying sorghum cultivars with 10 6 spores/ml of Johnsongrass isolates revealed high levels of hypersensitive response (HR), but relatively few cases of actual disease in terms of fungal reproduction were noted.
Inoculation of sorghum with C. sublineola has been shown to rapidly induce the production of pigmented phytoalexins (Lo et al., 1999), and confirmation of induction of mRNAs for phenylalanine ammonia lyase (PAL) and chalcone synthase, enzymes required for making these compounds was demonstrated by northern blots (Cui et al., 1996). Recently, the response to inoculation with C. sublineola has been measured in two large and diverse sorghum collections of accessions referred to as the minicore (Upadhyaya et al., 2009) and the Sorghum Association Panel (SAP; Casa et al., 2008). Since whole-genome sequencing data for these accessions, among others, has been made publicly available by Dr. Geoff Morris 1 (Morris et al., 2013) it has been possible to use Genome-Wide Association 1 http://www.morrislab.org/data Studies (GWAS) to tag genes by identifying Single Nucleotide Polymorphisms (SNP) alleles associated with a variety of traits. In the case of anthracnose, SNPs associated with resistance and susceptible cultivars have been identified for isolates from Puerto Rico (Cuevas et al., , 2019 and Texas Prom et al., 2019). In almost all cases, the gene nearest a high scoring SNP has been shown to have a role in host defense in other hosts/pathogen interactions, thus identifying over 30 genes that are likely to also be important in sorghum.
In this study, we have used RNASeq to investigate differences in gene expression in early defense responses of two sorghum cultivars that differ significantly in anthracnose sensitivity to try to pinpoint genes and pathways that may contribute to that difference, and to compare results with another very recent publication that did similar experiments but with emphasis on microRNAs with a different pair of anthracnose resistant and susceptible cultivars (Fu et al., 2020). Here, a mix of two C. sublineola isolates from Johnsongrass that had previously been shown lead to different visible responses on QL3 and Theis were used to examine gene expression changes in the two cultivars 24 h post-inoculation.

Collection of Plant Materials and Inoculation
Sorghum cultivars QL3 (resistant) and Theis (susceptible to most pathotypes of C. sublineola; Prom et al., 2012) were selected for the study and grown in the greenhouse to the eight-leaf stage. The plants were then transferred to a controlled growth chamber in the Plant Pathology building for inoculation and sampling. Two aggressive isolates of C. sublineola JG16081 and JG16083 were obtained from Dr. Gary Odvody at the Texas Agrilife Research and Extension Center, Beaumont, Tx. A 1:1 mixture of spores of these two isolates collected following growth on ½ strength PDA and adjusted to 10 6 /ml of water with a drop of Tween20 as a wetting agent was used for inoculum. Spores were applied to pre-marked segments on leaves of three plants of each cultivar. Zero-time samples were collected and frozen immediately in liquid nitrogen and stored at −70°C before the RNA isolation. At 24 h, the inoculated leaf segments were collected in the same manner.

Total RNA Isolation and RNA-Seq Library Preparation
Total RNA from the leaf tissues of three biological replicates was isolated using RNeasy Plant Mini Kit (Omega Bio-tek, United States). The RNA was treated with DNAseI (Qiagen, United States) to remove co-isolated genomic DNA and subsequently purified using RNeasy MinElute Cleanup Kit (Qiagen, United States). The quality and quantity of the total RNA were analyzed using Agilent 2100 Bioanalyzer and Qubit 4 Fluorometer (Invitrogen, United States) instruments, respectively. According to the manufacturer's protocol, RNA sequencing libraries were prepared using NEBNext ® Ultra™ Frontiers in Genetics | www.frontiersin.org II RNA Library Prep Kit (NEB, United States). The mRNAs were enriched from total RNA using magnetic beads with Oligo (dT), and subsequently, they were fragmented into shorter fragments using fragmentation buffer. The first-strand cDNA was synthesized from the fragmented mRNA using random hexamer primers and converted into double-strand cDNA. The resulting double-strand cDNAs were end-repaired, and Illumina sequencing adapters were added. The adapter-ligated libraries were amplified using sequencing primers for enrichment. The library's quality and insert size were determined using a bioanalyzer (Invitrogen, United States), and the library was quantified using fluorometric analysis (Invitrogen, United States). The library was diluted to 4 nM concentration and sequenced using Illumina's NextSeq 500 platform with paired-end sequencing chemistry. The resulting image files in the bcl format were converted to FASTQ with 2 × 75 bp reads using the bcl2fastq tool (Illumina, United States).

RNA-Seq Data Analysis
The sequencing adapters and the low-quality reads (Phred score QV < 30) were filtered using the read trimming tool Trimmomatic v.0.39 (Bolger and Giorgi, 2014). The quality-filtered reads were mapped to the S. bicolor reference genome build "assembly Sorghum_bicolor_NCBIv3" (Cooper et al., 2019) using STAR RNA-Seq aligner v.2.7.9 (Dobin et al., 2013) to generate alignments. The read count table was created using BAM alignment, general feature format (GFF) of genome annotation, and HTSeq v.0.13.5 (Anders et al., 2015). The differentially expressed genes (DEGs) among different experimental pair-wise combinations were identified using the DESeq2 R package (Love et al., 2014). The DEGs were filtered based on the minimum log 2 Fold Change of 1 and false discovery rate (FDR) of 0.05. Gene Annotation and GO Enrichment analysis were performed using Omicsbox. 2 Pathway mapping was performed using the KEGG database, MapMan v.3.6.0 (Thimm et al., 2004;Schwacke et al., 2019) and Pathview web tool (Luo et al., 2017). In these analyses, q values were used to minimize the number of potential false positives associated with the use of p values (Storey and Tibshirani, 2003). Plant disease resistance genes were identified using the plant disease resistance genes database (Osuna-Cruz et al., 2018).

Physiological Responses Observed
In order to determine susceptibility to the mixture of two isolates, greenhouse inoculation was conducted in 2017 and 2019 with 19 sorghum lines at eight leaf stage as described in Prom et al. (2012). Among the tested lines, QL3 was evaluated to be one of the most resistant lines, while Theis was scored to be the most susceptible line. Even though, we describe sorghum lines as susceptible based on necrosis on leaves, acervuli formation was not seen during greenhouse inoculations.
At 7 days post-inoculation (dpi), small necrotic lesions were shown on a few leaves in QL3, while Theis showed severe necrotic lesions. By 28 days, the differences were dramatic, as shown in Supplementary Figure 1.

RNA Sequencing and Identification of DEGs
RNA sequencing was used to identify the transcriptome response to the anthracnose infection in resistant (QL3) and susceptible (Theis) cultivars of Sorghum. The transcriptome analysis was performed with the samples at 0 and 24 h after infection (hpi). A total of 321, 992, 578 reads were obtained, which included 281,724,736 quality-filtered reads. The filtered reads were aligned to the Sorghum reference genome with an average mapping percentage of 88.4%. A summary of the RNA-Sequencing and genome mapping is given in Table 1.
The number of genes showing DEGs between the 0 and 24-h samples is shown in Table 2, where the cutoff used was a Log 2 fold change > 1. More DEGs for both the Resistant and Susceptible cultivars had decreased rather than increased expression. Comparison of the DEGs for the two hosts for samples taken 24 hpi showed that 3,160 were shared, while approximately 2,000 differed ( Figure 1A). As shown in the volcano plot' in Figure 1B, while most of the expression differences were close to the cutoff, a significant number had more than 5-fold changes, and even some of the more minor changes were highly significant, indicating consistency among replicates. Because exposure to pathogens in sorghum is known to involve the induction of defense-response genes such as those encoding enzymes, including those needed to make PAL-pathway derived phytoalexins and chitinase that can degrade fungal cell walls, emphasis was placed on identifying DEGs that are expressed to a higher degree in QL3 than Theis. A list of the top 25 of the highest scoring genes with increased expression 24 hpi in QL3 and the relative increase, if any in Theis, is shown in Table 3. A complete list of significantly up or downregulated DEGs in QL3 (resistant) and Theis are available in Supplementary Tables 1 and 2. DEGs common to both the hosts is available in Supplementary Table 3.

Characteristics of DEGs
Scanning the RNASeq defined DEGs that are above the statistical cutoff, including those that increased or decreased in both or only in one cultivar or the other, turned up several interesting observations. With a single exception, in all 3,160 cases where the same gene was identified as being differentially expressed in both cultivars, both were either upregulated or both showed downregulation; the exception was a chaperone protein dnaJ C76, chloroplastic isoform X2, which was upregulated 1.03-fold in QL3 and downregulated 1.19-fold in Theis. As in this case, the vast majority of relative differences for the cultivars, up or down, were within 2-fold. The most extreme cases include a lichenase-2 precursor, which was up 10.8-fold in QL3 and only 3.9-fold in Theis. Lichenases are involved in the hydrolysis of D glucans in plant cell walls during growth, and a lichenase has been shown to inhibit the growth of Fusarium in vitro (Mirakhorli et al., 2019). A lichenase has also been shown to be differentially upregulated in potatoes in response to Alternaria alternata (Brouwer et al., 2020). An opposite effect was seen for a peroxidase mRNA that was 9.1-fold increase in Theis but only 3.8-fold increase in QL3, implying a greater level of the oxidative burst in the susceptible cultivar, which suggests a greater degree of pathogen ingress may occur. Searching for genes that based on previous results, were predicted to show induced levels of expression within 24 hpi gave negative or inconsistent results. For example, expression of a chitinase-2 gene was up 2.6-fold in QL3 and 1.9 in Theis but mRNA levels for chitinase-6 were down 1.8 and 2.0-fold, respectively in the two cultivars. QL3 also had a chitinase-1 that was up 2.7 and a chitinase-11 that was down-1.1, neither of which were among the DEGs for Thies. While this could indicate activation of a defense-activated chitinase concomitant with inactivation of a different chitinase, there is not sufficient information on the specific genes to assign those roles. For phenylalanine ammonia lyase, the first step required to make the sorghum phytoalexin apegeninidin, mRNA levels actually decreased 2.8 and 3.5-fold in QL3 and Theis, respectively. Flavone synthase, required to make luteolinidin, another antifungal that originates from a shared precursor (Poloni and Schirawski, 2014) did not appear, but the expression of a downstream gene that encodes isoflavone-3-hydroxylase was downregulated by factors of 3.3 (QL3) and 6.0 (Theis). No evidence for altered expression of pathogenesis response gene PR10, which is often found in sorghum soon after inoculation, was detected. However, PR-4 was up three in both as was PR3, while three PR1 genes upregulated between 2 and 3-fold in Theis were not among the DEGs seen in QL-3. Of the genes that showed dramatic decreases in specific gene expression in one or both cultivars, the only genes recognized as having a potential role in host defense were the previously listed flavonol synthase/flavanone 3-hydroxylase and a probable cinnamyl alcohol dehydrogenase 5 that functions in the PAL pathway. Its expression was down 9.3-fold in QLE and 5.2-fold in Theis.
Many of the 3,160 DEGs showing significant changes between the zero time and 24-hpi samples for both cultivars have potential roles in host defense, where even small changes in expression can have enormous downstream consequences. This includes genes with LRR domains that are characteristic of cloned resistance genes, genes that function in signal transduction pathways, and transcription factors, including zinc-finger proteins or WRKY motifs, among others ( Table 4). As a consequence, we elected to take advantage of various gene function classification systems and pathway enrichment analysis to seek further understanding.

Pathway Enrichment Analysis
After identifying DEGs via sequence alignment to the annotated sorghum genome, pathway enrichment analysis was used to compare the relative numbers of genes from multiple biological pathways among the genes identified as having differential expression in QL3 compared to the Theis (Figure 2). MapMan analysis revealed several genes for biotic stress differentially expressed in the two hosts 24 hpi, including genes related to cell wall metabolism, redox function, and hormone responses (Figure 3). Mapping the DEGs to the pathways of plant-pathogen interactions reveled important genes differentially expressed in QL3 vs. Theis in 24 hpi (Figure 4).

GO Enrichment Analysis of DEGs
Gene Ontology (GO) enrichment analysis helps to understand DEGs regulation under the three important categories: biological processes, molecular functions, and cellular components. GO enrichment analysis of upregulated DEGs from resistant genotype (QL3) revealed the functional categories activated with the Johnsongrass anthracnose infection (Figure 6). The significantly enriched biological processes upregulated in QL3 upon infection included several stress responses processes and plant development processes. This includes response to stress (103 DEGs), response to abiotic stimulus (54 DEGs), cellular response to stimulus (58 DEGs), signal transduction (58 DEGs), cell communications (64 DEGs), anatomical structural development (53 DEGs), and post-embryonic development (31 DEGs).

DISCUSSION
As indicated previously, past work has shown that initial responses to exposure to Colletotrichum induce defense response reactions in both resistant and susceptible hosts. Examples include visual production of PAL-derived phytoalexins (Lo et al., 1999) and increased levels of mRNAs for phytoalexin pathway enzymes such as PAL and chalcone synthase (Cui et al., 1996). In addition, cytological studies showed that spore penetration occurs within 24 h in both the resistant and susceptible hosts, but that spread from the initially penetrated cells was not seen in the resistant host (Wharton and Julian, 1996). Thus, it has been assumed that resistant cultivars must either trigger an earlier or higher level of activation of genes critical for the defense response or that resistance involves genes other than those for which similar levels of induction occur in both the resistant and susceptible hosts. Consequently, it was a surprise to see that several of the genes such as PAL that were expected to be induced in both cultivars actually decreased in expression. A potential explanation besides the fact that different cultivars with different forms of resistance were used is that the earlier studies were based on northern blots or qPCR primers that were not specific to gene family members, some of which may respond to different signals. Since fold level changes do not reflect actual mRNA levels, only changes from the control, plus, and minus DEG values may not reflect actual levels of expression. There may also be cases when short sequences from conserved regions of gene family members get assigned to different genes in mapping to the genome. Even so, RNAseq data provide an excellent opportunity to explore the potential molecular basis for resistance further. For example, as shown in Table 4, 18 of 23 genes identified as transcription factors were more highly expressed in QL3 24 hpi than in Thies. Also, mRNAs for a DNA cytosine methylase, an enzyme that leads to cytosine methylation in the di-or trinucleotide sequences CG or CXG, were highly enriched only in QL3. The epigenetic effect of 5 C methylation in promoters is associated with turning off transcription and is also a downstream consequence of PTGS, a mechanism of combating RNA viruses in plants (Deleris et al., 2016). Several of the significantly upregulated genes detected are protein kinases that typically activate proteins that function in signal transduction pathways that function between recognizing pathogen elicitors (Antolín-Llovera et al., 2012) and are reflected in Figure 3 as map kinases. Basic blue proteins in plants are copper-containing plantacyanins that are regulated via conserved microRNA408 (miR408), which is present in sorghum (Zhang et al., 2011). Although, miR408 was not identified among the upregulated miRNAs detected by RNAseq in anthracnose-resistant sorghum cultivar SC383 following inoculation with C. sublineolum (Fu et al., 2020), a miR408 induced basic blue protein in wheat has been demonstrated to have a role in resistance to stripe rust (Puccinia striiformis f. sp. tritici; Feng et al., 2013). The ent-cassadiene C2-hydroxylases function in pathways for converting geranylgeranyl-PP into a variety of compounds, including terpene-based phytoalexins.
Frontiers in Genetics | www.frontiersin.org In rice, the expression of these enzymes is regulated by a microRNA encoded by DCL1; inactivation of DCL1 altered the production of terpene-based phytoalexins in rice, thus increasing sensitivity to rice blast caused by Magnaporthe oryzae and another fungal pathogen Fusarium fujikuroi (Salvador-Guirao et al., 2019). The presence of these genes in the upregulated transcripts implies that these genes may also play a role in multigenic pathogen resistance in QL3. Table 3) may also reveal new avenues important to active host defense. As an example, lichenase 2 precursor was at the top of the list of upregulated genes in both hosts. While lichenase [an endo-(1−>4)-beta-glucanase] has not been reported as a defense-related gene in sorghum, it has also been found to be highly induced in banana 27 hpi with Fusarium oxysporum f. sp. cubense in a race-specific manner . Chitinase 2 was upregulated in both hosts, but chitinase six expressions were reduced by nearly the same fold level, implicating chitinase 2 as more likely to be active host defense. Furthermore, there are a number of LRR and F box proteins, some of which are upregulated and others downregulated that may be important: leucine-rich repeat (LRR) regions of a protein allow interactions between proteins, including those of host and pathogen origin. Most known resistance genes in plants have an LRR domain (Gururani et al., 2012). F-box proteins are best known for functioning in growth regulation and development guided by hormone signals, but at least one is known to be involved in cell death triggered by pathogen recognition as part of the (HR; van den Burg et al., 2008). We found 3,150 DEGs downregulated in susceptible cultivar (Theis), and these genes are important to study the functions affected by Johnsongrass anthracnose infection in Theis. Complete list of up/downregulated DEGs from Theis available in Supplementary Table 2. Interestingly, early nodulin-93 was the most downregulated gene in the susceptible cultivar upon Johnsongrass anthracnose infection. This gene was downregulated more than 10-fold in Theis. This gene is essential for maintaining cell structure and nitrogen-use efficiency in plants (Okubara et al., 2000;Bi et al., 2009). The downregulation of this gene might have affected the growth and health of susceptible cultivar.

Enrichment Analysis
Results gained from the various bioinformatics resources based on the proteins encoded by the enriched mRNAs are also informative. The databases used classify functions in cellular processes (GO), pathway components (KEGG and BioCys), and evolutionary relationships (Panther). As seen in Figure 2, the largest category of upregulated genes in QL3 involves secondary metabolites, which would include the flavonoid and PAL pathway genes for making sorghums' antimicrobial compounds and coumarins that are listed. A surprise in the list includes camalexin, a defense compound made by Arabidopsis thaliana, but it turns out that the pathway components are analogous to the enzymes sorghum uses to make a similar defense compound, dhurrin (Klein et al., 2013). It is also of interest that the highest relative expression in QL3 (Supplementary Table 1) was an enzyme involved in benzoate biosynthesis. While benzoate is widely used as a preservative, some sorghum fungal have a Cyt P450 enzyme that can detoxify them (Korošec et al., 2014). However, sorghum also makes cinnamic acid derivatives that can inhibit the detoxification (Nguyen et al., 2015).

Notable DEGs Based on the MapMan Biotic Stress Pathway
A sorghum equivalent of an Arabidopsis gamma-Tubulin ring complex (gamma-TuRC) gene was found to be 1.2-fold upregulated in resistant cultivar QL3 in response to the infection. This protein is important for the plant cytoskeleton, including microtubule organization and elongation (Zeng et al., 2009). Plant cytoskeleton provides a major defense against invading fungal pathogens and is the major regulator and target for biotic interaction in plants (Takemoto and Hardham, 2004). Additionally, another protein, MT plus-end-tracking protein (EB1), is expressed at 2.3-fold higher 24 hpi than immediately after inoculation in QL3. This protein is also important for the microtubule bundling and organization in plants, which is essential for the plants biotic stress response and plant growth (Molines et al., 2018). The results infer that microtubule organization played an important role in QL3 in biotic stress response. Calcium cation channel (DMI1/Pollux Castor) expressed 1.3-fold only in QL3 with anthracnose infection. DMI1 has been shown to interact with Ca 2+ −permeable cyclic nucleotidegated (CNGC) channels at the nuclear envelopes. DM1proteins function as ligand-gated K + channels that can modulate the membrane potential across the nuclear envelope and work together with CNGC channels to modulate nuclear Ca 2+ release. CNGCs play a crucial role in pathways related to cellular ion homeostasis, development, and as a "guard" in defense against biotic stress (Jha et al., 2016). Sorghum cDNAs that mapped to Arabidopsis CNGCs were expressed about 3.8-fold higher in QL3 in response to infection. CNGCs play an important role in the innate immune response in plants (Ma et al., 2009). Three CNGCs (CNGC1, CNGC 4, and CNGC 9) were found to be upregulated (3.8, 2.3, and 1.3-fold, respectively) in QL3 in response to infection when compared to Theis (Supplementary Table 1).
RPM1 is a plant immune receptor that explicitly recognizes pathogen-released effectors to activate effector-triggered immunity (ETI) in Arabidopsis thaliana (Grant et al., 1995). Sorghum genes equivalent to Arabidopsis RPM1 were more highly expressed in resistant cultivar QL3 than the susceptible Theis (Figure 3). Three different sorghum homologs of Arabidopsis RPM1 were observed to be uniquely upregulated in the QL3 file (1.6, 1.6, and 1.3-fold, respectively), and two others for different Arabidopsis RPM genes were up 1.7 and 1.0. Expression of an RPM homolog found in Theis decreased 6.1-fold. Overall, these observations clarify the value of grouping genes by GO and functional categories based on RNASeq data to reveal genes that may play a role in hast diseases that would be difficult to detect on a one gene at a time basis.

CONCLUSION
We have analyzed the transcriptome changes associated with the inoculation of C. sublineola in two inbred sorghum cultivars (QL3 and Theis) known to differ in their anthracnose disease development. Differential gene expression analysis identified over 3,000 specific genes in both cultivars as showing significant changes in expression following inoculation. Gene ontology and pathway enrichment analysis identified upregulation of 23 transcription factor encoding genes as well as genes involved in the production of secondary metabolites, which are part of a typical host defense reaction.

DATA AVAILABILITY STATEMENT
The raw RNA sequencing reads generated from the current study have been deposited at NCBI under the BioProject accession number PRJNA736032.

AUTHOR CONTRIBUTIONS
CM and UR conceived the study and designed the research experiments. PN, EA, UR, RP, CM, and LP conducted the experiments, analyzed the data, and wrote the manuscript. All authors contributed to the article and approved the submitted version.