Serial Transcriptome Analysis Reveals Genes Associated with Late Blight Resistance in Potato Cultivar Qingshu 9

: The late blight caused by Phytophthora infestans ( P. infestans ) has been a major constraint for worldwide sustainable potato production. Chinese potato cultivar Qingshu 9 has shown excellent ﬁeld performance against late blight. To understand the mechanism underlying its resistance, the transcriptomic dynamics of Qingshu 9 infected with P. infestans was systematically investigated to identify the resistance-related genes. In total, 2027 shared differentially expressed genes (DEGs) were identiﬁed when infected with P. infestans at 24, 48 and 72 h post-inoculation. Based on k-mean cluster analysis, the DEGs showed four types of expression patterns, suggesting their divergent functions in response to P. infestans infection. KEGG enrichment analysis showed that the signiﬁcantly enriched DEGs were involved in the biosynthesis of secondary metabolites, plant–pathogen interaction, and photosynthesis. Furthermore, 980 transcription factor genes belonging to 68 families were found in the DEGs, of which AP2-EREBP and MYB genes were the most abundant. Moreover, many genes related to disease resistance showed differential expression during infection. Finally, the expression of nine DEGs was validated by quantitative real-time PCR. These results provide valuable information for understanding late resistance in potato cultivar Qingshu 9. related to the resistance response of Qingshu 9 to P. infestans . These results indicate that 24 hpi is the critical period of interaction between Qingshu 9 and P. infestans . At the same time, the number of upregulated genes in Subcluster_2 increased signiﬁcantly at 48 hpi, whereas they were downregulated at 72 hpi, so it is speculated that these genes are mainly involved in the intermediate stage of disease resistance. Many genes in Subcluster_1, Subcluster_2, and Subcluster_4 were upregulated at 48 hpi, indicating that 48 hpi was also critical to the immune response of Qingshu 9. In addition, the genes in Subcluster_3 displayed downregulated expression at all four time points, indicating that they played a negative regulatory role in the interaction between Qingshu 9 and P. infestans . the transcriptome level at different times. Our results indicated that 24 h was the most critical period of interaction between Qingshu 9 and P. infestans .


Introduction
The late blight, caused by P. infestans, is regarded as the most severe disease affecting potato production worldwide. To prevent or alleviate this disease, excessive fungicides are commonly used in potato production, not only leading to increased economic input, but also posing a substantial safety threat to the environment and human health [1]. It has been reported that the economic loss caused by late blight in potato is more than 10 billion dollars annually, and P. infestans has developed tolerance to fungicide metalaxyl used in recent years [2,3]. Extensive studies have demonstrated that breeding for potato varieties with resistance to P. infestans is the safest and most effective approach to controlling this disease.
Potato resistance to P. infestans could be classified as possessing either qualitative or quantitative traits, and the former has been well characterized and widely incorporated into potato breeding [4]. For instance, a qualitative resistance trait to P. infestans has been identified in Mexican hexaploid wild species Solanum demissum (S. demissum), and numerous resistance genes (R1~R11) have been defined and even isolated [5,6]. However, those resistance genes found in S. demissum are all race-specific, and their functions in P. infestans resistance have all been lost with the emergence of new and constantly mutating races. Therefore, it is obviously inefficient to use a few single-disease-resistant R genes to improve P. infestans resistance in potato breeding [7]. Some researchers suggest that tolerant potato cultivars likely employ quantitative defense response against Pectobacterium carotovorum subsp. brasiliense [8]. Therefore, quantitative traits with broad-spectrum resistance have the potential to significantly improve the durability of potato resistance to P. infestans. Compared to qualitative resistance, quantitative resistance is much more difficult to use due to its polygenic nature of inheritance and our poor understanding of its underlying molecular and biochemical mechanisms [9]. The application of quantitative resistance genes with broad-spectrum and durable resistance to P. infestans in potato breeding is one of the effective ways to prevent disease in potato [10].
Several studies have demonstrated that many genes are involved in pathogen recognition and plant defense during the colonization process of P. infestans in potato [11][12][13][14][15][16]. Identifying candidate resistance genes induced by pathogens in plant disease resistance is the main method for studying plant disease resistance [17,18]. The sequencing of P. infestans and the potato genome provides great convenience for studying potato-pathogen interactions and candidate genes involved in potato late blight resistance [19,20]. Extensive studies have been performed to identify genes associated with potato-microbe interactions. For example, plastidic carbonic anhydrase (CA) was found to play a crucial role in regulating the compatible vs. incompatible interactions between Solanum tuberosum and P. infestans [21]. DeepSAGE transcriptome analysis on the interaction of susceptible and resistant potato genotypes with P. infestans showed that more genes associated with defense response predominate in the compatible interaction than in the incompatible interaction [22]. The SuperSAGE technique has also been employed to identify novel candidate resistance genes in potato to late blight which do not compromise late maturity [23].
With the rapid development of sequencing technologies, RNA-seq technology provides an efficient, accurate, and high-throughput method for the identification of candidate genes, as it is able to fully detect gene expression at specific times within specific tissues [24], which is also increasingly used to mine genes potentially associated with late blight resistance in potato. For example, the RNA-seq technique was employed to obtain some insight into the mechanism of potato foliage and tuber defense against P. infestans, R gene dosage and shared biochemical pathways that contribute to RB-mediated resistance in the compatible potato-P. infestans interaction [25]. The RNA-seq technique has also been employed to compare differential response in the compatible and incompatible interactions of potato with P. infestans, and 40 effector targets were identified [26]. The transcriptomic response of susceptible RG and resistant PI613981 apple fruit to wounding and inoculation with Penicillium expansum was studied using RNA-Seq, and the results showed that genes involved in the MYB-domain transcription factor family were over-represented in the resistant genotype [27]. Furthermore, 3354 differentially expressed genes (DEGs) were identified in late blight resistant potato genotype SD20 infected with P. infestans, including genes encoding transcription factors and protein kinases, and four NBS-LRR proteins [28]. By creating a new analysis method, the transcriptome analysis of three wild Solanum species compared to three potato clones identified over 400 expressed putative resistance (R-)genes [29].
Potato cultivar Qingshu 9 is highly adaptable to diverse geographic regions in China and has shown excellent field performance against late blight for many years [30]. Although it does not display full immunity to P. infestans in the field, the disease incidence is significantly delayed, largely as a result of limited lesion expansion and development. This "speed limit" disease resistance may exert less selection pressure on the P. infestans populations, potentially leading to the delayed emergence of new races and enabling a relatively persistent disease resistance. Although Qingshu 9 has been widely planted, little is known regarding the molecular mechanism underlying its late blight resistance.
In this study, we employed RNA-seq method to examine transcriptomic dynamics in Qingshu 9 in response to P. infestans infection. The DEGs were identified and were functionally enriched through Gene Ontology (GO) and Kyoto Encyclopedia of Genes and Genomes (KEGG) analysis. The expression of nine DEGs was further verified by reverse transcription-quantitative real-time polymerase chain reaction (qRT-PCR). The obtained results provide useful information and insights on the late blight resistance in Qingshu 9, an important variety highly adapted to diverse environment conditions and widely planted in China.

Plant Materials
The potato cultivar Qingshu 9 was provided by the Qinghai Academy of Agricultural and Forestry Sciences (Qinghai, China) and was grown in a greenhouse at 24 • C under a 16-h photoperiod for four-six weeks. The P. infestans strain Pk5106 was collected from the experimental field of Qinghai Academy of Agriculture and Forestry Sciences (Xining, Qinghai, China). The strain was grown in a mold incubator under the condition of 16-h photoperiod at constant 16 • C regime with a relative humidity of 95%.

Preparation of Spore Suspension and Incubation
Following collection from the field, the P. infestans strain Pk5106 was isolated and purified prior to infecting potato tuber slices 2-3 times on a Petri dish containing rye agar medium supplemented with 6% sucrose. The inoculation lasted for 8-15 days at 16 • C to induce the formation of sporangia, which were then detected by light microscopy. Precooled sterilized water was added to the Petri dishes to completely infiltrate the hyphae. The Petri dishes containing the inoculations were rinsed repeatedly with a pipette to obtain a sporangial suspension, which was kept in a freezer at 4 • C for 1-2 h to stimulate the release of zoospores. The final concentration of spores was adjusted to approximately 4 × 10 4 mL −1 .

Collection of Samples and Microscopic Examination
The spore suspension above prepared was sprayed onto the Qingshu 9 plants, which were then covered with moist transparent plastic bags for a period of 24 h to maintain high moisture conditions, as required for infection. Two leaves were collected from potato plants for the control group before spore-spraying, and two more leaves of the same size were collected at 24 h, 48 h and 72 h post-infection (hpi). Three biological replications of the samples were carried out, and immediately frozen in liquid nitrogen then stored at −80 • C.
To view infection and colonization of P. infestans in Qingshu 9, leaf samples were transferred into polypropylene tubes which were filled with diluted trypan blue solution (10 g phenol, 10 mL glycerol, 10 mL lactic acid, 10 mL water and 10 mg of trypan blue). The tubes (lid slightly unscrewed) were treated in a heated water bath and boiled for 2 min. After cooling to room temperature, tubes were boiled for another 2 min. Then, the samples were destained by replacing the staining solution with chloral hydrate solution (5 g chloral hydrate dissolving into 2 mL water) for 24 h. The samples were finally mounted in distilled water and viewed under an Olympus BX51 (Shinjuku-ku, Tokyo, Japan) microscope with differential interference contrast optics [30].

RNA Extraction, Library Preparation, and Sequencing
Total RNAs from the P. infestans-infected Qingshu 9 leaves as well as the uninfected control leaves were extracted. Sequencing libraries were generated using the NEBNext ® Ultra™ RNA Library Prep Kit for Illumina ® (NEB, Ipswich, MA, USA) following the manufacturer's instructions. The clustering of the index-coded samples was performed on a cBot Cluster Generation System using the TruSeq PE Cluster Kit v3-cBot-HS (Illumina) following the manufacturer's instructions. The libraries were sequenced on an IlluminaHiSeq 4000 platform with PE 150, which conducted by the Novogene Bioinformatics Institute (Beijing, China).

Analysis of the Differentially Expressed Genes
Clean reads were obtained by removing reads containing adapters or poly N and reads of low quality from the raw reads using FastQC. The reference genome of the doubled  [31]. Index files of the reference genome were built using HISAT2 [32], which was used to align the clean reads against the reference genome. HTSeq v0.9.1 was used to count the read numbers mapped to each gene. The expected number of fragments per kilobase of transcript sequence per million base pairs (FPKM) was calculated for each gene [33]. The FPKM value of 0.1 or 1 was used as a threshold for determining whether a gene is expressed. Differential expression analysis was performed using the DESeq R package (1.18.0) [34], which provided statistical strategies for determining differential expression in digital gene expression data using the negative binomial distribution model. The resulting p-values were adjusted using Benjamini and Hochberg's approach for controlling the false discovery rate [35]. Genes with the fold change > 2 and adjusted p-value (padj) < 0.05 were considered as significant DEGs.

GO and KEGG Enrichment Analyses of DEGs
The DEGs were sequence-aligned with the GO and KEGG databases. GO enrichment analysis of DEGs was implemented by the GOseq R package, in which gene length bias was corrected [36]. GO terms with corrected p-values less than 0.05 were considered to be significantly enriched for DEGs. KEGG pathway enrichment (http://www.genome. jp/kegg/, accessed on 15 October 2018) was used to identify significantly enriched signal transduction pathways or metabolic pathways among the DEGs [37]. KOBAS (2.0) was used for the calculation and analysis of KEGG pathway enrichment factor (KS) [38].

Transcription Factor (TF) Analysis of DEGs
Plant TFs were predicted with iTAK software, the TFs were identified and classified according to the method as described by Perez-Rodriguez et al. [39]. The hmmscan program was used to identify TFs among the DEGs [40]. The potential role of the TFs in the interactions between Qingshu 9 and P. infestans was further analyzed based on the relevant biological processes in which those differentially expressed TFs might be involved in.

Comparative Analysis of Plant R Gene Database (PRGDB)
The DEGs were compared with the plant R gene database (http://prgdb.crg.eu/wiki/ Main_Page, accessed on 21 November 2018). On the basis of the comparison results, we selected the annotated genes related to potato disease resistance.

Validation of the Expression of DEGs by qRT-PCR
Nine DEGs closely related to plant disease resistance were selected to validate the expression patterns using qRT-PCR. The qRT-PCR primers as provided were synthesized by Shanghai Shenggong Bioengineering Co., Ltd. (Table 1). The potato Actin gene (GenBank accession no. X83206) was used as the reference gene. PrimeScript™ RT reagent Kit with gDNA Eraser (Perfect Real Time) was used for the first-strand reverse transcription.

Leaves at 0, 24, 48, 72 hpi Were Chosen to Prepare in Planta Materials
Inoculated potato leaf materials were continuously harvested at different time points. To choose the right time of infection, the infection process of Qingshu 9 sprayed with P. infestans was monitored by microscopic observation. At 0 hpi, the potato leaves were not infected by P. infestans, and normal stomas could be seen clearly ( Figure 1A,E). At 24 hpi, sporadic necrotic spots could be observed in the potato leaves ( Figure 1B,F). With increasing infection time, necrosis spots gradually increased and the area of necrosis gradually increased at 48 hpi and 72 hpi ( Figure 1C,D,G,H). On the basis of these results, the time points of 24 hpi, 48 hpi and 72 hpi were chosen to prepare infestation material to perform deep sequencing, using 0 hpi as a control.

A Large Amount of Credible Data Were Obtained through Transcriptome Sequencing
To investigate the transcriptional response of Qingshu 9 to P. infestans, the inoculated leaf tissues were sampled at 0 hpi, 24 hpi, 48 hpi and 72 hpi, and total RNA was isolated for RNA sequencing, which were performed with three biological replicates. The sample symbols were Pi0h_1, Pi0h_2, Pi0h_3, Pi24h_1, Pi24h_2, Pi24h_3, Pi48h_1, Pi48h_2, Pi48h_3, and Pi72h_1, Pi72h_2, Pi72h_3. In total, 588 million raw paired-end reads were generated, and an average of 47,753,443, 50,630,987, 49,302,732 and 48,372,645 raw reads from Pi0h, Pi24h, Pi48h and Pi72h were produced ( Table 2). The correlation coefficient (R 2 ) between different biological replicates was very close to 1 (Table 3), indicating the high reproducibility between replications and the reliability of sequencing data. Clean reads were obtained after quality control and used for subsequent mapping, assembly and gene expression calculation. When the clean reads were mapped onto the potato DM reference genome by HISAT 2.0, approximately 71.02-84.9% reads were mapped (Table 2), which is higher than the values of 78.7% [25] and 57.29-72.13% previously reported [28]. Finally, 24,251, 23,579, 23,412 and 23,736 genes were detected on average in the sampling at the four time points, including a large number of plant R genes.

A Large Amount of Credible Data Were Obtained through Transcriptome Sequencing
To investigate the transcriptional response of Qingshu 9 to P. infestans, the inoculated leaf tissues were sampled at 0 hpi, 24 hpi, 48 hpi and 72 hpi, and total RNA was isolated for RNA sequencing, which were performed with three biological replicates. The sample symbols were Pi0h_1, Pi0h_2, Pi0h_3, Pi24h_1, Pi24h_2, Pi24h_3, Pi48h_1, Pi48h_2, Pi48h_3, and Pi72h_1, Pi72h_2, Pi72h_3. In total, 588 million raw paired-end reads were generated, and an average of 47,753,443, 50,630,987, 49,302,732 and 48,372,645 raw reads from Pi0h, Pi24h, Pi48h and Pi72h were produced ( Table 2). The correlation coefficient (R 2 ) between different biological replicates was very close to 1 (Table 3), indicating the high reproducibility between replications and the reliability of sequencing data. Clean reads were obtained after quality control and used for subsequent mapping, assembly and gene expression calculation. When the clean reads were mapped onto the potato DM reference genome by HISAT 2.0, approximately 71.02-84.9% reads were mapped (Table 2), which is higher than the values of 78.7% [25] and 57.29-72.13% previously reported [28]. Finally, 24,251, 23,579, 23,412 and 23,736 genes were detected on average in the sampling at the four time points, including a large number of plant R genes.    Because the disease incidence in Qingshu 9 is significantly delayed, we proposed exploring the resistance molecular mechanism of Qingshu 9 by screening the differentially expressed genes during the infection process. DEGs were identified using DESeq software with the corrected Padj < 0.05 and fold change > 2 as a threshold for significant differential expression. We identified 5035, 6093, 9419 and 2027 shared DEGs in the three groups Pi24h vs. Pi0h, Pi48h vs. Pi0h and Pi72h vs. Pi0h (Figure 2A), of which 6289 were upregulated and 7023 were downregulated ( Figure 2B). Obviously, the number of DEGs gradually increased with the extension of infection time. Although different infection time points contained many specific DEGs, there were still 371 and 1459 DEGs that were commonly up-and downregulated among the groups ( Figure 2B). These DEGs can be induced by infection with P. infestans, and also play an important role in the resistance of Qingshu 9 to late blight.

24 hpi and 48 hpi Are Important in the Immune Response of Qingshu 9
Because the disease incidence in Qingshu 9 is significantly delayed, we proposed exploring the resistance molecular mechanism of Qingshu 9 by screening the differentially expressed genes during the infection process. DEGs were identified using DESeq software with the corrected Padj < 0.05 and fold change > 2 as a threshold for significant differential expression. We identified 5035, 6093, 9419 and 2027 shared DEGs in the three groups Pi24h vs. Pi0h, Pi48h vs. Pi0h and Pi72h vs. Pi0h (Figure 2A), of which 6289 were upregulated and 7023 were downregulated ( Figure 2B). Obviously, the number of DEGs gradually increased with the extension of infection time. Although different infection time points contained many specific DEGs, there were still 371 and 1459 DEGs that were commonly upand downregulated among the groups ( Figure 2B). These DEGs can be induced by infection with P. infestans, and also play an important role in the resistance of Qingshu 9 to late blight. The Venn diagram on the left (A) represents total differential genes, and the one on the right (B) represents differential genes that are upregulated and downregulated.
To further analyze the expression patterns of DEGs, k-mean cluster analysis was performed. There were four different expression patterns obtained, with 5366, 1581, 3265 and 434 genes, respectively ( Figure 3). The genes in Subcluster_1 and Subcluster_4 were up- To further analyze the expression patterns of DEGs, k-mean cluster analysis was performed. There were four different expression patterns obtained, with 5366, 1581, 3265 and 434 genes, respectively ( Figure 3). The genes in Subcluster_1 and Subcluster_4 were upregulated on the whole, and the number of upregulated genes at 24 hpi was greater, especially in Subcluster_4. Among them, the expression of genes in Subcluster_1 was slightly upregulated at the four time points, and Subcluster_4 showed the highest expression level at 24 hpi. Due to the high expression of the gene Osmotin (PGSC0003DMG400003044) and the significant upregulation of many genes, including Auxin repressed/dormancyassociated gene (PGSC0003DMG400029732) and PR1 (PGSC0003DMG400005113), at 24 hpi, we believe that these genes are induced by early plant immunity and are closely related to the resistance response of Qingshu 9 to P. infestans. These results indicate that 24 hpi is the critical period of interaction between Qingshu 9 and P. infestans. At the same time, the number of upregulated genes in Subcluster_2 increased significantly at 48 hpi, whereas they were downregulated at 72 hpi, so it is speculated that these genes are mainly involved in the intermediate stage of disease resistance. Many genes in Subcluster_1, Subcluster_2, and Subcluster_4 were upregulated at 48 hpi, indicating that 48 hpi was also critical to the immune response of Qingshu 9. In addition, the genes in Subcluster_3 displayed downregulated expression at all four time points, indicating that they played a negative regulatory role in the interaction between Qingshu 9 and P. infestans. is speculated that these genes are mainly involved in the intermediate stage of disease resistance. Many genes in Subcluster_1, Subcluster_2, and Subcluster_4 were upregulated at 48 hpi, indicating that 48 hpi was also critical to the immune response of Qingshu 9. In addition, the genes in Subcluster_3 displayed downregulated expression at all four time points, indicating that they played a negative regulatory role in the interaction between Qingshu 9 and P. infestans.

DEGs among the Groups Are Enriched in Different GO Terms and KEGG Pathways
To better understand the functions of DEGs, GO enrichment analysis was performed. In addition to the phosphorus metabolic process and the phosphate-containing compound metabolic process, which had the largest numbers of enriched genes, kinase activity and protein kinase activity also enrich many DEGs in the group Pi24h vs. Pi0h ( Figure  4A). The DEGs in group Pi48h vs. Pi0h were mainly enriched with respect to biological

DEGs among the Groups Are Enriched in Different GO Terms and KEGG Pathways
To better understand the functions of DEGs, GO enrichment analysis was performed. In addition to the phosphorus metabolic process and the phosphate-containing compound metabolic process, which had the largest numbers of enriched genes, kinase activity and protein kinase activity also enrich many DEGs in the group Pi24h vs. Pi0h ( Figure 4A). The DEGs in group Pi48h vs. Pi0h were mainly enriched with respect to biological processes and metabolic processes in the biological processes category ( Figure S1A). Additionally, cytoplasm and protein kinase activity were also significantly enriched. Similarly, DEGs in group Pi72h vs. Pi0h were mainly enriched in BP category and some DEGs were also enriched with respect to kinase activity and protein kinase activity ( Figure S2A). The results indicated that many DEGs could be annotated to kinase activity (GO:0016301; 700DEGs) and protein kinase activity (GO:0004672; 616DEGs). As intermediate molecules in signal transduction pathways, protein kinases (PKs) played an important role in the interaction between Qingshu 9 and P. infestans. Combined with the analysis of expression patterns, it was found that except for group Pi24h vs. Pi0h, the genes enriched in the other two groups were mainly of the subcluster_1 type ( Figures 4A, S1A and S2A). GO enrichment results showed that the co-upregulated DEGs were mainly enriched in the biological process category and the cellular component category ( Figure 5A). The co-downregulated differential genes were mainly enriched in the molecular function category, which also includes kinase activity and protein kinase activity ( Figure 5B).   The KEGG pathway enrichment analysis was further used to identify the metabolic pathways that were involved in Qingshu 9 response to P. infestans. It was found that DEGs of Pi24h vs. Pi0h were enriched in numerous metabolic pathways, including the biosynthesis of secondary metabolites (ID: sot01110, 274 DEGs), plant-pathogen interaction (ID: sot04626, 69 DEGs), photosynthesis, and zeatin biosynthesis ( Figure 4B). The enriched pathways of DEGs of Pi48h vs. Pi0h included biosynthesis of secondary metabolites (ID: sot01110, 270 DEGs), plant-pathogen interaction (ID: sot04626, 71 DEGs), ribosome and so on ( Figure S1B), while that of Pi72h vs. Pi0h included plant-pathogen interaction (ID: sot04626, 116 DEGs), ribosome, ribosome biogenesis in eukaryotesis ( Figure S2B). These results suggested that biosynthesis of secondary metabolites, plant-pathogen interaction pathways were mainly involved in the resistance response of Qingshu 9 to P. infestans. Therefore, the DEGs involved in these pathways were considered to participate in the Qingshu 9-P. infestans interaction.
By comparing the groups of Pi24h vs. Pi0h, Pi48h vs. Pi0h and Pi72h vs. Pi0h, we found that the abundance of genes related to disease resistance, such as PR1, NBS-LRR, WRKY transcription factor and cell wall protein, was gradually increasing. Additionally, the gene cytochrome P450 and ethylene response factors associated with hypersensitive response and many lignin synthesis-related genes were significantly induced by P. infestans infection. These results suggest that Qingshu 9 activates a series of regulated and functional genes to respond to P. infestans infection. The resistance of Qingshu 9 is controlled by multiple genes.

Transcription Factors Play Crucial Roles in Qingshu 9's Resistance Response
Previous studies have shown that TFs, as key components of the plant defense system, activate the expression of stress-related genes and participate in biological processes [42][43][44][45]. In this study, a total of 980 differentially expressed TFs belonging to 68 TF families

Transcription Factors Play Crucial Roles in Qingshu 9's Resistance Response
Previous studies have shown that TFs, as key components of the plant defense system, activate the expression of stress-related genes and participate in biological processes [42][43][44][45]. In this study, a total of 980 differentially expressed TFs belonging to 68 TF families were identified. Among them, the AP2-EREBP and MYB families were the most abundant group, with 86 and 74 members, respectively, accounting for 8.78% and 7.55% of the total TFs. Furthermore, more than 30 members of bHLH (55)

Comparison of Results between DEGs and Plant R Gene Database
To explore whether the selected DEGs contain plant R genes, we scanned the DEGs into the plant R gene database. The results showed that 2532 DEGs were plant R genes, of which 656 genes were further annotated in detail. These annotated genes can be divided into 12 R gene types, namely CN type (18), CNL type (169), Mlo-like type (21), N type (76), NL type (164), RLK type (29), RLK-GNK2 type (11), RLP type (117), RPW8-NL type (1), T type (11), TNL type (34), and other types (5). Further screening results showed that 311 genes were annotated as potato resistance genes, thus providing candidate genes for the selection and use of potato disease resistance genes.

The Results of qRT-PCR and RNA-Seq Are Consistent
To validate the results of the RNA-Seq, qRT-PCR was performed on nine randomly selected DEGs. The proteins encoded by these genes were usually involved in plant disease resistance, such as zinc finger protein 593, PR1 protein, heat shock cognate protein 80, WRKY transcription factor-30, and salicylic acid-inducing protein 19. As illustrated in Figure 7, the expression patterns of these nine genes were consistent with the results of the RNA-Seq, which validated the reliability of RNA-Seq data.

Comparison of Results between DEGs and Plant R Gene Database
To explore whether the selected DEGs contain plant R genes, we scanned the DEGs into the plant R gene database. The results showed that 2532 DEGs were plant R genes, of which 656 genes were further annotated in detail. These annotated genes can be divided into 12 R gene types, namely CN type (18), CNL type (169), Mlo-like type (21), N type (76), NL type (164), RLK type (29), RLK-GNK2 type (11), RLP type (117), RPW8-NL type (1), T type (11), TNL type (34), and other types (5). Further screening results showed that 311 genes were annotated as potato resistance genes, thus providing candidate genes for the selection and use of potato disease resistance genes.

The Results of qRT-PCR and RNA-Seq Are Consistent
To validate the results of the RNA-Seq, qRT-PCR was performed on nine randomly selected DEGs. The proteins encoded by these genes were usually involved in plant disease resistance, such as zinc finger protein 593, PR1 protein, heat shock cognate protein 80, WRKY transcription factor-30, and salicylic acid-inducing protein 19. As illustrated in Figure 7, the expression patterns of these nine genes were consistent with the results of the RNA-Seq, which validated the reliability of RNA-Seq data.

Potato Cultivar Qingshu 9 Contains Special Resistance to P. infestans
Qingshu 9 is an elite cultivar with high adaptation to diverse geographic regions in China. It has shown excellent field performance against late blight in different regions for many years. Its "speed limit" disease resistance model may alleviate the selection pressure on the P. infestans population, which might delay the emergence of new races with virulence. This special resistance is urgently required for the selection of disease-resistant potato varieties. However, the underlying molecular and biochemical mechanisms of Qingshu 9's resistance remain elusive.

Transcription Analysis of Qingshu 9 Infected by P. infestans
A comprehensive understanding of the physiological and biochemical processes of potato interaction with P. infestans at the molecular level would help in the further development of control strategies, leading to the effective control of late blight in potato. Hostpathogen interaction studies have been carried out in plants. Many studies have applied RNA-Seq technology to identify genes involved in the response to potato late blight at the genome level [25,28,29,46,47]. However, few studies identifying special resistance genes to potato late blight have been reported. In this study, a comprehensive transcriptome analysis of the interaction process between Qingshu 9 and P. infestans was performed, which could contribute to better understanding the resistance response mechanism in Qingshu 9. In the resistance reaction of Qingshu 9 to P. infestans, 0 h was used as the control group, and the treatment groups were 24 h, 48 h and 72 h after inoculation. This experiment was able to distinguish the treatment groups and the control group at the transcriptome level at different times. Our results indicated that 24 h was the most critical period of interaction between Qingshu 9 and P. infestans.

Potato Cultivar Qingshu 9 Contains Special Resistance to P. infestans
Qingshu 9 is an elite cultivar with high adaptation to diverse geographic regions in China. It has shown excellent field performance against late blight in different regions for many years. Its "speed limit" disease resistance model may alleviate the selection pressure on the P. infestans population, which might delay the emergence of new races with virulence. This special resistance is urgently required for the selection of disease-resistant potato varieties. However, the underlying molecular and biochemical mechanisms of Qingshu 9's resistance remain elusive.

Transcription Analysis of Qingshu 9 Infected by P. infestans
A comprehensive understanding of the physiological and biochemical processes of potato interaction with P. infestans at the molecular level would help in the further development of control strategies, leading to the effective control of late blight in potato. Host-pathogen interaction studies have been carried out in plants. Many studies have applied RNA-Seq technology to identify genes involved in the response to potato late blight at the genome level [25,28,29,46,47]. However, few studies identifying special resistance genes to potato late blight have been reported. In this study, a comprehensive transcriptome analysis of the interaction process between Qingshu 9 and P. infestans was performed, which could contribute to better understanding the resistance response mechanism in Qingshu 9.
In the resistance reaction of Qingshu 9 to P. infestans, 0 h was used as the control group, and the treatment groups were 24 h, 48 h and 72 h after inoculation. This experiment was able to distinguish the treatment groups and the control group at the transcriptome level at different times. Our results indicated that 24 h was the most critical period of interaction between Qingshu 9 and P. infestans.

Transcription Factors May Play an Important Role in the Resistance Response to P. infestans in Qingshu 9
Transcription factors (TFs) play important roles in the interaction between plants and pathogens. The AP2/ERF family, WRKY family, MYB family and bZIP family are involved in the immune response to external stimuli and biotic stress in plants [48][49][50]. The AP2/ERF is the largest family of TFs in higher plants, and many of its members have been reported to act as activators or suppressors in response to biotic stress in plants [51]. In the interaction between plants and pathogens, ERF1 and ERF2 transcriptionally activated the expression of related disease resistance genes, while ERF3 and ERF4 inhibited gene expression in plants during the disease resistance process [52]. It has been reported that the ERF transcription factor StERF3 negatively regulates resistance to P. infestans and salt tolerance in potato [53]. In addition, the MYB, WRKY, bZIP and NAC TF families are well known for their roles in plant defense responses [54]. Genes of transcriptional regulation play a regulatory role in multiple pathways of potato against P. infestans. For example, StWRKY1 enhances resistance to P. infestans by regulating the expression of genes in the phenylpropionic acid biosynthetic pathway [55]. Similarly, StWRKY8 is involved in the regulation of the biosynthesis and regulation of the benzyl isoquinoline alkaloid, an antibacterial agent precursor, which could inhibit the early infection of P. infestans [56]. Moreover, potato NAC43 and MYB8 are able to regulate the secondary cell wall biosynthesis to prevent P. infestans infection [57]. In this study, there were 980 differentially expressed TFs, belonging to 68 TF families, mainly including the AP2-EREBP, MYB, bHLH, WRKY, HB, C2H2, NAC, GRAS, bZIP, and Orphans families. Among them, the AP2-EREBP family, with 86 members were the most numerous. Different TF genes were found to display differential expression during P. infestans infection in Qingshu 9, and these genes may positively or negatively regulate the response to P. infestans and co-regulate the resistance of Qingshu 9. These results provide a resource for future research into resistance genes, and the cloning of disease-resistant TF genes is of great significance for genetic engineering to improve disease resistance in potato.

Protein Kinases May Participate in the Interaction between Qingshu 9 and P. infestans
Previous studies have shown that plant receptor protein kinases recognize the pathogenassociated molecular patterns of pathogens as well as activated and transmitted immune signals [58]. In this study, we found that the molecular functional ontology of group Pi24h vs. Pi0h and group Pi72h vs. Pi0h both exhibited enriched kinase activity, phosphotransferase activity, alcohol group as acceptor, and protein kinase activity. Many DEGs could be annotated to kinase activity and protein kinase activity. Calcium-dependent protein kinases (CDPKs or CPKs) have been found to be involved in growth and developmental processes as well as in defense against diverse environmental stresses [59]. Reactive oxygen species (ROS) have been implicated in plant innate immunity and calcium-dependent protein kinases with regulation of the production of ROS by potato NADPH oxidase [60]. The CDPK family members are particularly involved in plant-pathogen interaction, and StCDPK7 has emerged as a relevant player that could be manipulated to deploy disease resistance in potato [61]. StCDPK5 has been shown to phosphorylate the N-terminal region of plasma membrane RBOH proteins, and participate in StRBOHB-mediated ROS burst [62]. Our results suggest that at the 24 hpi and 72 hpi of infection, protein kinase may play an important role in signal recognition and signal transduction, contributing to Qingshu 9's resistance.

Disease Resistance-Related Genes May Enhance the Disease Resistance of Qingshu 9
In this study, DEGs were significantly enriched in the categories of biosynthesis of secondary metabolites, plant-pathogen interaction pathways, and transcriptional regulation, etc., which are considered to play important roles in Qingshu 9's resistance to late blight and is similar to the results found for plant-pathogen interaction in previous reports [63][64][65][66][67]. In addition, a large number of disease resistance-related genes, such as PR1, NBS-LRR, cell wall protein, and many lignin synthesis-related genes were significantly induced by P. infestans infection. Previous studies have demonstrated that in the defense response process of potato against P. infestans, secondary metabolites such as phenylpropanoid and cell wall synthesis-related biosynthesis participated in the plant's disease resistance defense by synthesizing lignin and strengthening the cell wall [57]. The lignin content of the resistant potato variety Orion increased significantly compared to the susceptible variety Majestic after infection [68]. The increase in the enzyme activity of phenylalanine ammonia lyase, and lignin content were able to enhance the resistance of potato tubers to P. infestans [69]. After inoculation with P. infestans, the Phytophthora-resistant variety Kennebec showed increased lignin content and early expression of defense-related genes PAL-1, HMG-2, PR-1 and PR-5 compare to the susceptible variety Russet Burbank [70]. Similarly, the lignin content in tomato roots after 16 days of Verticillium dahliae infection was 43% higher in resistant varieties than the susceptible varieties [71]. In this study, many lignin synthesis-related genes were induced, we speculated that the high expression of lignin biosynthesis-related genes could strengthen the cell wall of Qingshu 9 and then enhance its resistance against cellular invasion by P. infestans. Additionally, the gene cytochrome P450 and ethylene response factors associated with hypersensitive response were also induced by P. infestans infection. These results suggest that the infection of P. infestans activates a series of resistance-related genes, and that the resistance of Qingshu 9 is controlled by multiple genes and may owed the quantity resistance. In addition, the 311 potato R genes were obtained by comparison with the plant R gene database. Among them are known genes, such as R1 and R3a, but most of them are unnamed new genes, which will provide a large number of candidate genes for researching the disease resistance of Qingshu 9.

Conclusions
In summary, through transcriptomic analyses of Qingshu 9 infected by P. infestans at four time points, we systematically identified genes associated with host-pathogen interaction. GO and KEGG analysis of the DEGs provided some insights into molecular mechanisms underlying late blight disease resistance of Qingshu 9. This study not only provides important candidate genes for further functional study, but also has potential to contribute to genetic improvement of the tolerance to late blight in potato and beyond.