Comparative transcriptome analysis between resistant and susceptible tomato lines uncovers the response mechanism of Cf-16-mediated resistance to Cladosporium fulvum

Leaf mold disease caused by Cladosporium fulvum is a major disease in cultivated tomato plants and affects global tomato production. Some Cf genes, of which Cf-16 is an effective gene for resisting tomato leaf mold, are associated with leaf mold resistance; however, the molecular mechanism is largely unknown. we used comparative transcriptome analysis of C. fulvum-resistant (cv. Ontario7816, including the Cf-16 gene) and C. fulvum-susceptible (cv. Moneymaker) tomato lines to identify differentially expressed genes (DEGs) at 4 and 8 days postinfection with C. fulvum. We found that the number of DEGs in the Cf-16 tomatoes was significantly higher than the number of DEGs in the Moneymaker tomatoes. In addition, 1,350 DEGs were shared among Cf-16 groups at 4 and 8 dpi, suggesting the existence of a common core of DEGs in response to C. fulvum infection. Upregulated DEGs were mainly associated with defense processes and phytohormone signaling, including salicylic acid (SA) and jasmonic acid (JA), in the Cf-16 tomato. Moreover, the SA and JA contents significantly increased in the Cf-16 tomato at the early stages of C. fulvum infection. Comprehensively, more upregulated DEGs were found in the Cf-16 tomato than in the Cf-10 and Cf-12 tomatoes at the early stage of C. fulvum infection. However, the significantly enriched defense-signaling pathways involved in Cf-16 had some distinctions from those in Cf-10 and Cf-12. Our results provide new insights into the resistance response mechanism of Cf genes to C. fulvum, especially the unique characteristics of Cf-16 in response to C. fulvum infection. comparative identify differentially

2 Abstract Background Leaf mold disease caused by Cladosporium fulvum is a major disease in cultivated tomato plants and affects global tomato production. Some Cf genes, of which Cf-16 is an effective gene for resisting tomato leaf mold, are associated with leaf mold resistance; however, the molecular mechanism is largely unknown.

Results
we used comparative transcriptome analysis of C. fulvum-resistant (cv. Ontario7816, including the Cf-16 gene) and C. fulvum-susceptible (cv. Moneymaker) tomato lines to identify differentially expressed genes (DEGs) at 4 and 8 days postinfection with C.
fulvum. We found that the number of DEGs in the Cf-16 tomatoes was significantly higher than the number of DEGs in the Moneymaker tomatoes. In addition, 1,350 DEGs were shared among Cf-16 groups at 4 and 8 dpi, suggesting the existence of a common core of DEGs in response to C. fulvum infection. Upregulated DEGs were mainly associated with defense processes and phytohormone signaling, including salicylic acid (SA) and jasmonic acid (JA), in the Cf-16 tomato. Moreover, the SA and JA contents significantly increased in the Cf-16 tomato at the early stages of C. fulvum infection. Comprehensively, more upregulated DEGs were found in the Cf-16 tomato than in the Cf-10 and Cf-12 tomatoes at the early stage of C. fulvum infection. However, the significantly enriched defensesignaling pathways involved in Cf-16 had some distinctions from those in Cf-10 and Cf-12.

Conclusion
Our results provide new insights into the resistance response mechanism of Cf genes to C. fulvum, especially the unique characteristics of Cf-16 in response to C. fulvum infection. Background 3 Tomato is not only the second-most important vegetable crop worldwide [1,2] but also an important model plant for fleshy fruit development and plant-pathogen interactions; tomato genome sequencing was completed in 2012 [3], and breeding tomato for disease resistance is one of the current critical necessities. Leaf mold disease is considered one of the most devastating diseases in tomato (Solanum lycopersicum) and is caused by Cladosporium fulvum. C. fulvum is a nonobligate, abiotrophic pathogen. This fungus infects foliage and occasionally petioles and stems [4][5][6]. Leaf mold occurs worldwide and has long been prevalent in many countries, reducing both fruit yield and quality and causing serious economic losses, especially for covered field-grown tomatoes under hightemperature and high-humidity conditions in continuous tomato cropping fields [7]. Therefore, leaf mold is regarded as a major threat to tomato production. However, the most effective method prevent leaf mold is to breed tomato varieties containing genes resistant to leaf mold.
From the perspective of the coevolution of plants and pathogens, the plant immune system recognizes and responds to pathogens in the following stages. In the first stage of plant-pathogen interactions, PRRs (pattern recognition receptors) in the plant recognize PAMPs (pathogen associated molecular patterns) of the pathogens, thereby inducing PTI (PAMP-triggered immunity) and preventing the colonization of the pathogens. In the second stage, the pathogens successfully bypass PTI and secrete their own effectors into plant cells. At this moment, the plant cells do not have the antipathogenic proteins capable of recognizing these effectors, causing the effector-triggered susceptibility response (ETS). In the third stage, plants have gradually evolved to produce NB-LRR proteins that directly or indirectly recognize specific effectors of pathogens, inducing ETI (effector-triggered immunity). Finally, pathogens avoid plant immune responses by inhibiting or altering effectors that could be recognized by plants and generate new 4 effectors that could not be recognized by plant NB-LRR, successfully infecting plants and causing ETS. Meanwhile, the plant's immune system has evolved with changing pathogens, and new R genes could reidentify new effectors in pathogens and induce ETI again.
In the process of plant-pathogen interactions, different pathogens carry the corresponding avirulence (AVR) genes, and the AVR genes are recognized by effector proteins [8]. The proteins are secreted into the apoplastic space during infection, which induces either compatible or incompatible interactions between the fungus and infected plant [9].
Incompatible interactions (chlorosis) appear when plants resist pathogens, thus leading to the hypersensitive response (HR), whereas the compatible interaction occurs when the pathogens could grow and ramify, causing necrosis in the infected cells [10][11][12][13]. However, the tomato-C. fulvum interaction follows a typical gene-for-gene relationship, and the products of C. fulvum-resistance genes (Cf genes) in tomato could specifically recognize the products encoded by AVR genes in corresponding pathogens, thus leading to the HR [14,15]. At least 24 Cf genes have been reported since the discovery of the Cf-1 gene in the 1930s [16,17], and these genes have been introduced into cultivated tomato [18][19][20][21][22][23][24][25][26].
Recent advances in transcriptome sequencing (RNA-seq) have strongly accelerated research on host-pathogen interactions in plants, such as rice [27], maize [28], cucumber [29], watermelon [30] and strawberry [31]. Currently, Avr4/Cf-4-and Avr9/Cf-9-dependent defense gene expression has been completed by cDNA-amplified fragment length polymorphism (cDNA-AFLP) [32]; the Avr5 gene has been cloned by a combined bioinformatic and RNA-seq-based transcriptome sequencing approach [14]; Cf-19mediated, Cf-12-mediated and Cf-10-mediated resistance to C. fulvum in tomato has been characterized by cDNA-AFLP and RNA-seq, respectively [33][34][35]; however, there have been few transcriptomic studies on Cf-16-mediated resistance. In this study, comparative transcriptome analysis was used to identify differentially expressed genes (DEGs) between 5 C. fulvum-resistant and C. fulvum-susceptible tomatoes. In addition, we measured the contents of endogenous hormones of the two tomato lines in response to C. fulvum. Our study may provide a basis for cloning the Cf-16 gene, which will be useful for understanding the response mechanism of Cf resistance to C. fulvum infection in tomato and breeding tomato varieties resistant to leaf mold.

Plant materials and C. fulvum inoculation
The two tomato lines, the resistant line Ontario7816 including the Cf-16 gene (kindly provided by Institute of Vegetable and Flowers, Chinese Academy of Agricultural Science) [19] and the susceptible line Moneymaker without the C. fulvum-resistance genes (kindly

Microscopic observation of C. fulvum in the Cf-16 tomato
To identify the interaction process of Cf-16-mediated HR and key time points involved in the resistance mechanism, the lactophenol trypan blue staining method was carried out according to Franco's approach [49]. The leaf samples of the resistant and susceptible lines were harvested at 0-21 dpi, immediately stained, clarified overnight in chloral hydrate solution (2.5 mg/ml) [50], and examined using an Olympus SZX10 dissecting microscope (Olympus, Japan).

Endogenous jasmonic acid (JA) and salicylic acid (SA) measurement
Leaf samples from the inoculation and control treatments from the resistant and susceptible lines were harvested at 4, 8, 12, 16 and 21 dpi, respectively. Endogenous SA and JA were extracted from leaves using the modified method of Llugany et al. (2013) [51], and their contents were measured using HPLC-MS/MS (high-performance liquid chromatography-tandem mass spectrometry), which was performed by the Public Laboratory of College of Horticulture and Landscape Architecture, Northeast Agricultural University, Harbin, China. Samples were analyzed with an AB SCIEX Triple TOF5600 + mass spectrometer [52]. Three replicates were performed.

RNA extraction, cDNA library construction and sequencing
Total RNA from each treatment was collected at 4 and 8 dpi for RNA-seq and quantitative States). Then, the library was sequenced using BGISEQ-500, and 150-bp paired-end reads were generated.

Sequencing read mapping and identification of DEGs
Raw reads in FASTQ format were generated by base calling, statistically analyzed using SOAPnuke v1.4.0 and filtered using trimmomatic v0.36 [55]. Clean reads were obtained by removing reads with adapters, reads containing more than 5% poly-N (where N represents unknown bases), and low-quality reads (whose mass value was less than 10 and the proportion of the total number of bases in the reads was greater than 20%).

Gene ontology (GO) and KEGG (Kyoto Encyclopedia of Gene and Genomes) pathway analyses
The GO and KEGG pathway enrichment analyses of DEGs were accomplished using Phyper function of R software, and the GO terms and KEGG pathways with an adjusted P-value ≤ 0.01 were regarded as significantly enriched.

Gene coexpression network analysis
Gene coexpression network analyses were completed using the WGCNA (weighted correlation network analysis) package v1.48. Gene dendrograms were constructed with colors based on the correlations between the expression of genes. These were used to build clustering trees and to divide the modules. In addition, the correlation between modules and samples was also analyzed by weighted gene coexpression network analysis.

qRT-PCR analysis
Sixteen DEGs were validated using qRT-PCR to verify the expression profiles obtained by RNA-seq. qRT-PCR was performed using AceQ ® qPCR SYBR ® Green Master Mix (Vazyme, USA) on a qTOWER 3 G Detection System (Analytik Jena, Germany). Each sample was replicated three times, and data analysis was performed using the 2 -△△CT method [60].

Microscopic analysis of C. fulvum invasion in the two tomato lines
Light microscopy was used to observe the infection process between C. fulvum and Ontario7816 or Moneymaker leaves (Fig. 1). As shown in Fig. 1a, no difference between Ontario7816 and Moneymaker was observed at 0 dpi, hyphae grew into the stomata in both Ontario7816 and Moneymaker at 4 dpi (Fig. 1e). In Ontario7816 leaves, a small number of HR areas appeared at 8 dpi (Fig. 1b), and the necrosis area gradually grew at 10 dpi (Fig. 1c) until more necrosis lesions appeared in both mesophyll cells and leaf veins between 12 and 21 dpi (Fig. 1d). Hyphae emerged through the stomata in Moneymaker at 8 dpi (Fig. 1f), and the emerged hyphae continued to increase and grow at 10 dpi, with the last few infected cells starting to undergo necrosis at 10-21 dpi (Fig. 1 g-h). Based on this observation, we collected samples from each treatment at 4 and 8 dpi for RNA-seq and qRT-PCR verification.

Analysis of SA and JA responses to C. fulvum infection
To explore the SA and JA responses to C. fulvum infection, HPLC-MS/MS was used to measure their contents. As shown in Fig. 2, compared with the SA content of CK-Cf16, the SA content of Cf16 increased rapidly to a peak at 8 dpi and was far higher than the SA content of MM between 4 and 12 dpi. Moreover, the SA content of MM was generally lower than the SA content of CK-MM. The JA content of Cf16 had the greatest value at 4 dpi, followed by a rapid reduction between 12 and 21 dpi but had a higher value than the JA content of the other samples at 4-16 dpi. These results suggested that SA and JA played important roles in regulating the plant response and enhancing plant defense in tomato plants infected with C. fulvum and had a quick response at the early stages of infection.

Summary of RNA-seq data
To determine the transcriptome profiles of Ontario7816 and Moneymaker following C. fulvum infection, we performed RNA-seq analysis on these two lines at 4 and 8 dpi. Three biological replicates were made at each time point for each treatment. In this study, an average of ~6.87 Gb was generated from each sample using the BGISEQ-500 platform ( Fig. S1, Table S1). The raw data were deposited in the NCBI Sequence Read Archive under the accession number GSE133678. As shown in Table S1, more than 98% of reads were ≥ 20%, and more than 91% of clean reads had a quality score of ≥ 30%. After the reads were filtered, 64.15-72.64 million clean reads were generated, and at least 93.29% of these reads were mapped to the tomato reference genome, among which more than 78.26% were aligned to unique locations. Ultimately, 18,514 novel transcripts were generated with 12,790 unknown splicing events for known genes, 2,047 novel coding transcripts without any known features, and 3,677 transcripts for long noncoding RNA.

DEGs in response to C. fulvum
A gene was defined as significantly differentially expressed when an adjusted P-value ≤ 0.001 and log 2 fold-change ≥ 2. The two standards were used to identify DEGs in the resistant and susceptible lines in response to C. fulvum at 4 and 8 dpi. All FPKM values for every gene and the fold-changes and adjusted P-values for DEGs are shown in Tables S2 and S3, respectively. As shown in Table 1, the number of DEGs was highest in CK_Cf_4dpi-vs-Cf_4dpi, indicating that Ontario7816 carries more defense-response genes than Moneymaker, and the number of DEGs at 8 dpi was markedly less than the number of DEGs at 4 dpi between Ontario7816 and Moneymaker. However, the number of upregulated genes was higher than the number of downregulated genes in the two tomato lines at 4 dpi. Overall, 8,526 and 6,938 genes were differentially expressed in the resistant and susceptible lines at 4 dpi, respectively, of which 3,382 and 2,203 genes were up-and downregulated, respectively. In addition, among three comparisons CK_Cf_4dpi-vs-Cf_4dpi, CK_Cf_8dpi-vs-Cf_8dpi and Cf_4dpi-vs-Cf_8dpi, 1,350 DEGs were shared, which further suggests that a common group of genes was activated or deactivated upon C.

GO and KEGG enrichment analyses of DEGs
To determine the functions of DEGs involved in the response to C. fulvum, we performed GO classification and KEGG functional enrichment analyses with the Phyper function of R software. GO is divided into three major functional categories: biological process, cellular component and molecular function. For the DEGs in Ontario7816, the significant GO terms were mostly enriched in biological regulation, cellular process, metabolic process and response to stimulus in the biological process category, and these terms were related to disease resistance. In the cellular component category, most of the DEGs were assigned to cell, membrane, membrane part and organelle, which were found to be specific to the resistant line. The significantly enriched terms in the molecular function category were binding, catalytic activity, transcription regulator activity and transporter activity, and among them, binding and catalytic activity terms are known to play an important role in plant hormone signal transduction (Fig. 4). Therefore, the corresponding genes of these terms might play critical roles in response to C. fulvum infection.
KEGG pathway enrichment analysis was used to investigate the biological pathways associated with DEGs. As shown in Fig. 5a, the pathways "Plant hormone signal transduction" and "Plant-pathogen interaction" were significantly enriched, and the number of genes and the rich ratios of these two pathways were significantly higher than those of other pathways. Other disease-resistance pathways were also enriched, such as "MAPK signaling pathway-plant", "Benzoxazinoid biosynthesis" and "Phosphatidylinositol signaling system". Overall, "Plant hormone signal transduction" and "Plant-pathogen interaction" may be pivotal pathways in regulating the resistance response to C. fulvum infection in tomato. In the KEGG pathway analysis based on upregulated DEGs among CK_Cf_4dpi-vs-Cf_4dpi, CK_Cf_8dpi-vs-Cf_8dpi and Cf_4dpi-vs-Cf_8dpi and shown in Fig. 5b, the pathways "Plant-pathogen interaction" and "Plant hormone signal transduction" were also significantly enriched, further confirming the existence of a common core of DEGs in response to C. fulvum infection. Meanwhile, we performed KEGG pathway analysis based on upregulated unique DEGs in CK_Cf_4dpi-vs-Cf_4dpi (Fig. 5c) and found that the "Plantpathogen interaction" pathway was also significantly enriched, and 10 DEGs in this pathway were identified (Table 4). Table 2 shows that common DEGs related to disease-resistance pathways were significantly upregulated in Ontario7816 and Moneymaker at 4 dpi with C. fulvum. Among the DEGs of the "Plant-pathogen interaction" pathway, 26 genes were significantly differentially expressed between CK_Cf_4dpi-vs-Cf_4dpi and CK_MM_4dpi-vs-MM_4dpi, while 25 DEGs in the significantly enriched pathway "Plant hormone signal transduction" were identified (Table 3). Overall, plant hormones may play a key role in the Cf-16 tomato response to C. fulvum infection.  Table 3 and Table 4 were shown at the end of this document.

Gene coexpression network analysis
Weighted gene coexpression network analysis is a common algorithm for constructing gene coexpression networks [36]. A total of 13 different modules were obtained using a gene dendrogram colored according to the correlations between gene expression levels ( Fig. 6a). Among them, genes in MEred and MEgreenyellow were highly expressed in Cf_4dpi, and genes in MEpurple had a relatively high expression in Cf_4dpi and MM_4dpi (Fig. 6b). We performed KEGG analysis for the three modules. For the MEred module, pathways related to "Plant-pathogen interaction", "Oxidative phosphorylation" and "Phenylalanine, tyrosine and tryptophan biosynthesis" were enriched, whereas for MEgreenyellow, pathways related to "Pentose phosphate pathway", "Flavonoid biosynthesis", "Phenylpropanoid biosynthesis" and "Plant hormone signal transduction" were enriched (Fig. S2).

Validation of RNA-seq data by qRT-PCR
To verify the RNA-seq data, 16 DEGs were chosen for qRT-PCR using three biological replicates. These 16 genes were selected from significantly enriched KEGG pathways (such as "Plant hormone signal transduction", "Plant-pathogen interaction" and "Metabolic pathways"). The expression data of qRT-PCR were consistent with the RNA-seq results, indicating a similar trend between the transcriptome analysis and qRT-PCR data (Fig. 7).

Discussion
In this study, we identified the interaction process in the Cf-16 tomato by microscopic observation. An HR was observed in the Cf-16 tomato at 8 dpi. By contrast, hyphae emerged through the stomata and continued to increase and grow in Moneymaker at 8 dpi (Fig. 1). The systemic defense response of the resistance genes was activated at the early stage of C. fulvum infection. This result is consistent with the results of previous studies of other Cf genes [33][34][35]. We demonstrated that the Cf-16 tomato triggered a resistance response to the leaf pathogen C. fulvum by comparative transcriptome analysis between resistant and susceptible tomato lines using RNA-seq, and the reliability of RNA-seq data 13 was verified by qRT-PCR analysis. In summary, we identified many significant DEGs between the resistant and susceptible lines in response to C. fulvum infection, and the number of DEGs in the Cf-16 tomato was significantly higher than the number of DEGs in the Moneymaker tomato. The significant GO terms were mostly enriched in biological regulation, cellular process, metabolic process and response to stimulus in the biological process category, and these terms were related to disease resistance. Meanwhile, KEGG enrichment analysis indicated that most DEGs were classified into "Plant hormone signal transduction" and "Plant-pathogen interaction". In future work, we will be committed to a comprehensive comparative analysis of Cf-19, Cf-12, Cf-10 and Cf-16 to further explore the mechanism of Cf gene-mediated resistance response to C. fulvum infection.
Two main mechanisms of the plant innate immune system include PTI and ETI, which identify pathogens and form a defense response. First, the PRRs recognize a series of PAMPs and activate PTI in the early stages. At present, many PAMPs have been found in plant pathogens, and the corresponding PRR proteins have also been found in plants.
Among them, the most comprehensive PRR explored is FLS2 in Arabidopsis. In our study, the pattern recognition protein CERK1 (chitin elicitor receptor kinase 1, BGI_novel_G000519 and BGI_novel_G000515) was significantly upregulated in the Cf-16 tomato at 4 dpi. This result was in agreement with the research of Xue et al. (2017) in Cf- 12. In the future, we could conduct a more in-depth study on whether a higher expression level of CERK1 is involved in the activation of chitin signaling and affects the tomato-C. fulvum interaction process.
After recognition of infection, the Cf-16 tomato quickly established complex signal defense pathways. CNCGs (plant cyclic nucleotide-gated ion channels) are crucial in plant hypersensitivity, and these channels can regulate Ca 2+ influx. Ca 2+ can activate CDPKs (calcium-dependent protein kinases), which play important roles in plant responses to 14 both abiotic stress and pathogens [37,38]. In our study, CDPKs (100316879, 101055527 and 101255379) were expressed at higher levels during the early stage of infection in the Cf-16 tomato. This is consistent with previous studies suggesting that these genes play crucial roles in the response of the Cf-16 tomato to C. fulvum infection [39]. On the other hand, Ca 2+ combined with a CML (calcium binding protein) produces NO, which further promotes plant hypersensitivity responses or autoimmune reactions [40]. Interestingly, the results of our study showed that eleven CML genes were significantly expressed at 4 dpi in the Cf-16 tomato compared with the Moneymaker tomato (Fig. 8a). Particularly, 543942 and 101245711 were upregulated approximately eight-fold in the Cf-16 tomato compared with the Moneymaker tomato. This is consistent with previous studies suggesting that these genes are very important in the response of the Cf-16 tomato to C. fulvum infection. Ca 2+ is involved in signal transduction, especially in the early response of C. fulvum infection. At the same time, the signal triggers a series of protein kinases that rely on endocytosis [41,42]. FLS2 recognizes flg22 and subsequently activates downstream signaling pathways involving WRKY transcription factors to promote defense responses against bacterial and fungal pathogens and nematodes [43,44]. In our study, a total of 12 WRKY genes were specifically upregulated at 4 dpi, as shown by the hierarchical clustering of DEGs in both tomato lines (Fig. 8b) insects. The first NB-LRR protein that was found to interact directly with Avr was the rice CC-NB-LRR protein Pi-to, whose LRR domain was able to directly recognize the effector Avrpita of Magnaporthe oryzae and induce ETI [45]. However, it has been demonstrated that RPM1 is an NBS-LRR protein from Arabidopsis thaliana that confers resistance to Pseudomonas syringae. RPM1 is also related to the onset of the HR [46,12]. Our results showed that RPM1 (100736444, 101246473, 101246761 and 109120689) was significantly upregulated in the Cf-16 tomato at 4 dpi, and the expression levels were higher in Cf-16 than in Moneymaker. These genes may be very important in the response of the Cf-16 tomato to C. fulvum infection. Meanwhile, they may be candidate genes of Cf-16. With deepening research on the plant immune system, it has been suggested that an increasing number of disease-resistance mechanisms will be discovered, which will help us to more deeply understand the recognition mechanism of plant pathogen effectors.

Conclusions
In summary, microscopic, hormone content and RNA-seq analyses were performed to explore the interactions between the Cf-16 tomato and C. fulvum. The upregulation of CDPKs, CMLs, LRRs and other elicitor-response proteins indicated that PTI and ETI immunity systems may exist in the Cf-16 tomato to combat C. fulvum. When the mycelia of C. fulvum grew into the interspace of the stomata and mesophyll cells, the effector proteins secreted by C. fulvum were rapidly recognized by the Cf-16 tomato. This triggered downstream defense-signaling transduction associated with Ca 2+ channels and other pathways, including pathways involving SA and JA. Meanwhile, specific defense-related transcription factors, WRKYs, were activated as an early plant response, which actively regulated downstream resistance pathways. Finally, the HR was induced, and the C.

Consent for publication
Not applicable.

Availability of data and materials
All the data pertaining to the present study have been included in the tables and figures of the manuscript, and the authors are pleased to share all the data upon reasonable request.

Competing interests
The authors declare that they have no competing interests.

Funding
The publication costs of the article was funded by the China Agriculture Research System

Supplementary Files
This is a list of supplementary files associated with the primary manuscript. Click to download. Table S4.xls   Table S1.xls