Reprogramming of a defense signaling pathway in rough lemon and sweet orange is a critical element of the early response to ‘Candidatus Liberibacter asiaticus’

Huanglongbing (HLB) in citrus infected by Candidatus Liberibacter asiaticus (CLas) has caused tremendous losses to the citrus industry. No resistant genotypes have been identified in citrus species or close relatives. Among citrus varieties, rough lemon (Citrus jambhiri) has been considered tolerant due to its ability to produce a healthy flush of new growth after infection. The difference between tolerance and susceptibility is often defined by the speed and intensity of a plant’s response to a pathogen, especially early defense responses. RNA-seq data were collected from three biological replicates of CLas- and mock-inoculated rough lemon and sweet orange at week 0 and 7 following infection. Functional analysis of the differentially expressed genes (DEGs) indicated that genes involved in the mitogen activated protein kinase (MAPK) signaling pathway were highly upregulated in rough lemon. MAPK induces the transcription of WRKY and other transcription factors which potentially turn on multiple defense-related genes. A Subnetwork Enrichment Analysis further revealed different patterns of regulation of several functional categories, suggesting DEGs with different functions were subjected to reprogramming. In general, the amplitude of the expression of defense-related genes is much greater in rough lemon than in sweet orange. A quantitative disease resistance response may contribute to the durable tolerance level to HLB observed in rough lemon.


INTRODUCTION
Huanglongbing (HLB) or citrus greening is one of the most destructive plant diseases in the world. 1 The disease is caused by the bacterium, Candidatus Liberibacter asiaticus (CLas), a phloeminhabiting, α-proteobacteria, and is transmitted by the insect vector, Asian citrus psyllid (Diaphorina citri). 2 It is generally recognized that after infection or inoculation, CLas bacteria migrate through the phloem and accumulate there, resulting in the formation of sieve plugs, which contribute to HLB symptoms. 3 Since no toxins, cell wall degrading enzymes or specialized secretion systems have been identified in the CLas genome, it is believed that the disease symptoms are the result of host metabolic imbalances brought about by nutrient depletion or interference with nutrient transport. 4 The relationship between plant hosts and their associated microbes is essentially guided by a plant innate immunity system. Plants have evolved two layers of immune systems, pathogen-associated (or microbe-associated) molecular patterns (PAMP)-triggered immunity (PTI) and effectortriggered immunity. PTI is mediated by pattern-recognition receptors (PRRs) that recognize PAMPs, whereas effectortriggered immunity is mediated by resistance (R) proteins that recognize pathogen effectors. 5 The role of these immunity systems in citrus greening disease and the transcriptomic response of citrus is poorly understood, especially in regards to the early defense response of rough lemon to HLB infection.
No source of HLB resistance has been identified in citrus species or closely-related relatives. Folimonova et al. 6 classified the response of 30 citrus genotypes to HLB from sensitive to tolerant. They classified Eureka lemon and Persian lime as tolerant, with little or no HLB symptoms and no strong correlation between bacterial titer and disease severity. Zhang et al. 7 also reported that HLB-affected scions of lemon had a higher titer of CLas, survival rate, and pathogen transmission rate than pummelo (Citrus maxima) scions. Citrus limon Burm. F., represented by 'true' lemons such as 'Eureka', contains genome contributions from three species, (C. maxima, C. reticulata, and C. medica) and is distinct from rough lemon, though related. 8 Our previous study also showed that rough lemon is tolerant of HLB. 9 Once rough lemon trees are infected and symptomatic, they can be rejuvenated by the continued growth of new shoots with few or no foliar symptoms of the disease, and they repeat this cycle for many growing seasons. In contrast, sweet orange exhibits continuous growth inhibition and eventual dieback. 9 Understanding the differences in response of the tolerant and susceptible citrus genotypes to HLB is essential for developing genetic strategies that can produce tolerant, or perhaps even resistant, varieties. Most of transcriptomic studies of citrus-CLas interactions have been conducted on different tissues (leaves, fruit and root) of sweet orange, [10][11][12][13] and roots of tangerine; 14 while only two studies have focused on the leaves of tolerant genotypes of rough lemon 9 and US-897. 10 US-897, however, is an intergeneric hybrid of Citrus spp. and Poncirus spp, and may not provide a complete understanding of HLB response in citrus varieties. 9 Using a comparative approach for understanding plant-microbe interactions can be complicated because the response and effect of candidate genes in tolerant genotypes may differ in the genetic background of susceptible trees due to epistasis. Previous studies have indicated that resistance and susceptibility can often depend on the speed and magnitude of the defense mechanism employed and the ability of the pathogen to suppress host response. 4,15,16 It has been reported that infection of sweet orange with CLas does not lead to a significant induction of defense-related genes at the early stages of the infection process. 10,11 Therefore, it is important to learn if HLB-tolerant rough lemon recognizes the presence of CLas at an early stage of infection and whether or not a defense response is elicited. Leaves of citrus trees in the field can remain asymptomatic and the titer of CLas by reverse transcriptionquantitative PCR (RT-qPCR) can remain undetectable as long as two years after the initial psyllid transmission of CLas due to the need for the presence of a high bacterial titer for detection. 17 False negatives and false positives are common due to low titers of bacteria inside the phloem, the presence of PCR amplification inhibitors in the phloem sap, and non-specific amplification. 18 CLas in greenhouse-inoculated plants can remain undetectable for 5-9 weeks 10 or 20 weeks. 9 In a time course study of CLas infection, Fan et al. 9 reported that a greater number of stress response genes were distinctively modulated in rough lemon than in sweet orange. The study, however, did not identify any genes at a significant falsediscovery rate (FDR) rate in the early stages (5 weeks) of the infection. In another study, no significant induction of defenserelated genes was observed in sweet orange inoculated to CLas over a 5-to 9-week period of time, 10 . Both of the mentioned studies utilized an Affymetrix citrus microarray for their transcriptome analysis. However, it is not known how many unique citrus genes are actually represented in the chip. Misleading interpretations of microarray results can also occur due to non-specific hybridization. 14 Next-generation RNA sequencing technology, however, can reveal rare and unknown transcripts, thus offering a more precise and accurate picture of the transcriptome. 13,19 In the present study, a whole-genome transcriptional analysis of rough lemon and sweet orange leaves was conducted at 7 weeks post-inoculation in order to identify genes that were induced as part of an early response to CLas. Leaf samples collected immediately after inoculation (week 0) were used as a baseline, which was not done in previous transcriptome studies of CLas in citrus. 4,10,11 By using week 0 as a baseline, the natural variation among biological replications could be minimized and the number of genes whose expression was significantly affected by CLas could be maximized. The response of mock-inoculated and CLas-inoculated rough lemon and sweet orange were separately analyzed and the two datasets were subsequently compared, in order to identify biological mechanisms associated with each genotype and infection versus mock-inoculation. The transcriptome analysis identified statistically significant differentially expressed genes (DEGs) between rough lemon and sweet orange in response to CLas. A distinct difference in the defense response between rough lemon and sweet orange in CLas-inoculated leaves was observed. Not only were two distinct sets of DEGs identified but greater amplitude in the profile of defense response genes in rough lemon, relative to sweet orange, was also observed. The present study provides a comprehensive overview of the early transcriptional reprogramming that occurs in rough lemon in response to CLas.

Plant material and experimental design
Young trees were inoculated as described by Fan et al. 9 Briefly, two-yearold seedlings of rough lemon (C. jambhiri Lush.) and 'Madam Vinous' sweet orange (C. sinensis L Osb.) were graft-inoculated with 3-4 cm bud wood from CLas-infected Carrizo citrange (C. sinensis × P. trifoliata L. Raf.) trees maintained under greenhouse conditions in order to provide an inoculum source for various HLB-related experiments; Carrizo was used as the inoculum source because it is immune to citrus tristeza virus (CTV), a common contaminant in field source trees of CLas. Control seedlings were grafted with 3-4 cm budwood from healthy Carrizo trees. All these plants were kept in a United States Department of Agriculture Animal and Plant Health Inspection Service and Center for Disease Control-approved and secured greenhouse at the University of Florida, Citrus Research and Education Center, Lake Alfred. Leaf samples were collected from the middle of the plant crown and approximately 60-80 cm from the inoculation site. Three biological replicates were produced for each species in each treatment ( Figure 1). Reverse transcription-quantitative PCR (RT-qPCR) was performed to confirm the presence of CLas in the inoculum source and in inoculated plants as previously described. 8 Four fully expanded leaves were sampled separately from CLas-inoculated plants and mock-inoculated plants (used as controls) of each species at 0 and 7 weeks after inoculation. Leaves were immediately frozen in liquid nitrogen Figure 1. Phenotype of sweet orange and rough lemon trees 18 months after CLas-inoculation.
Reprogramming of a defense signaling pathway Q Yu et al. Reprogramming of a defense signaling pathway Q Yu et al.
and stored at − 80°C until further use. Three biological replicates were analyzed for each condition ( Figure 1). In total, 12 plants and 24 leaf samples were collected (2 species × 2 treatments × 3 replications × 2 time points). Selection of week 7 as the time point representing early response to CLas was based on results of a preliminary pilot experiment conducted to determine the time point at which the maximum number of DEGs changing in expression was first initiated in response to CLas.

PCR detection of CLas
DNA from leaf midribs and petioles was extracted using the Plant DNeasyMini Kit (Qiagen, Valencia, CA, USA) according to the manufacturer's instructions. RT-qPCR assays were performed as previously described. 20 Amplifications were performed in an Agilent Mx3005P (Agilent Technology) real-time PCR system using the Brilliant III Ultra-Fast QPCR Master Mix (Agilent Technology Inc, Waldbronn, Germany). Plants were considered PCR-positive when CT (cycle threshold) values were below 30.

RNA extraction and sequencing library construction
Total RNA from the 24 samples described above was extracted using TRIzol ® Reagent following the manufacturer's protocol. RNA was further purified using the TURBO DNA-free™ kit to eliminate genomic DNA. A NanoDrop Spectrophotometer (NanoDrop Technologies, Inc., Wilmington, DE, USA) and gel electrophoresis were used to assess RNA quality and quantity. RNA quality was further assessed using an Agilent 2100 Bioanalyzer (Aglient Technologies Inc, Waldbronn, Germany). Ribosomal RNA was removed from the total RNA using a Ribo-Zero rRNA removal kit for using plant seed/plant leaf following the manufacturer's protocol. Subsequently, 50 ng of the purified RNA from each sample were used for library construction utilizing an Epicentre ScriptSeq v2 RNA-Seq library preparation kit according to the user's guide. Illumina RNA library construction was performed at the Interdisciplinary Center for Biotechnology Research (ICBR) Gene Expression Core, University of Florida (UF). Briefly, rRNA-depleted RNA was fragmented by incubation at 85°C, and then reverse transcribed using random primers containing a 5ʹ-tagging sequence. A 3ʹ-tag was added using a terminal-tagging reaction resulting in Di-tagged, single-stranded cDNA. Following purification, the di-tagged cDNA was amplified by limited-cycle PCR, in order to add the Illumina adaptor sequences. The amplified libraries were purified using Agencourt AMPure beads (Beckman Coulter, catalog # A63881). The library size and mass was assessed by analysis in the Bioanalyzer. Typically, a 200-2000 broad library peak was observed with the highest peak at 500 bp. Quantitative PCR was used to validate the library's functionality, using the KAPA library quantification kit (Kapa Biosystems, catalog number: KK4824). Finally, the libraries were pooled in equimolar concentration and sequenced on an Illumina 2 × 100 HiSeq 2000 (Illumina Inc., San Diego, CA, USA)

RNA-seq data analysis
Raw sequencing reads were initially cleaned as follows. Ambiguous residues were trimmed off from both sides of the sequence. Bases with a Phred quality below 20 from the 3ʹ end of the sequence were removed. Reads shorter than 40 bases or those containing 410 bases with a quality rating below 20 were also discarded. In addition, reads consisting of repetitive single bases that accounted for 460% of the length at the 5ʹ or 3ʹ end were also discarded. Approximately 50-65% of reads were retained after the initial processing, providing 45-75 million clean reads for each of the 24 samples (Table 1). The Citrus clementina genome v1.0 (JGI) genomic sequence (http://www.citrusgenomedb.org/species/clementina/genome1.0) was used as reference genome for mapping the reads. Sequencing reads from each sample were mapped independently to the reference sequences using gmap v3 (http://research-pub.gene.com/gmap/src/gmap-gsnap-2012-07-20.v3.tar.gz). 21 This step successfully mapped all of the cleaned reads to the genome; resulting in a total of~56.3% of the reads that were uniquely mapped to the genome.
Gene expression values were determined as follows. The number of mapped reads for each individual gene was counted using an in house perl script. The gene annotation file in GFF format was downloaded from Phytozome v1.0 http://www.phytozome.net (Phytozome: Cclementina 182 Reprogramming of a defense signaling pathway Q Yu et al. v1.0.gene.gff3.gz). The DESeq package in R was used for gene expression analysis. 22 DESeq uses the negative binomial distribution, with variance and mean linked by local regression, to model the null distribution of the count data. Significant up-and downregulated genes were selected using two cutoffs: an adjusted P value of 0.05 and a minimum fold-change of 2.0. Gene sets from each combined sample and treatment were annotated using Blast2GO[67] to assign Gene Ontology (GO) terms to each gene. Lists of transcripts that were differentially expressed using a false-discovery rate (FDR) o0.05 in the pairwise comparisons were used in the Fisher's Exact Test in Blast2GO in order to identify GO terms that were significantly overrepresented. The biological interpretation of the DEGs was further assessed by assigning the genes to metabolic pathways using the Kyoto Encyclopedia of Genes and Genomes (KEGG). 23 Subnetwork enrichment analyses were run on the data set using Pathway Studio 10.0. 24 Reverse transcription-quantitative PCR Nine DEGs that were identified by RNA-seq as either up-or downregulated in response of CLas were selected for validation by RT-qPCR. The RT-qPCR analyses were performed in two steps. First strand cDNA was synthesized from 0.3 μg of total RNA using an Affinityscript QPCR cDNA Synthesis Kit (Agilent Technologies), according to the manufacturer's recommendations. Primers for the twelve selected genes were designed using Primer-BLAST software at NCBI (http://www.ncbi.nlm.nih.gov/tools/primer-blast/), and are listed in Supplementary Table S1 of Supplementary File 6. RT-qPCR was performed using a Brilliant III Ultra-Fast SYBR Green QPCR Master Mix (Agilent Technologies), following the manufacturer's recommended procedures. GAPDH was used as a reference gene to normalize the expression of the other analyzed genes. PCR reactions were performed using 2 μl of cDNA, 0.25 nM of each primer, and 10 μl of 2 × SYBRGreen PCR master mix (Aglient Technologies Inc, Santa Clara, CA, USA) in a 20-μl volume. A negative control was included using water as a template for each primer pair. Specificity of the amplification reactions was checked by post-amplification dissociation curves and by sequencing the reaction products. The fluorescent intensities were used to determine relative mRNAlevels with MxPro software.

RESULTS
CLas detection and HLB symptom development Leaf samples collected at week 7 did not exhibit the characteristic blotchy mottle appearance and were RT-qPCR negative for CLas. Continued RT-qPCR analysis for CLas was conducted every two weeks after inoculation. Positive confirmation of CLas in inoculated plants was not obtained in rough lemon and sweet orange until 23 weeks post inoculation. Even after 8 months, however, rough lemon did not exhibit any signs of growth inhibition, and continued growth of new shoots with few or no symptoms was observed; however, the typical blotchy mottled appearance was commonly found on mature, older leaves ( Figure 1). In contrast, severe levels of blotchiness were observed on mature and older leaves of CLas-inoculated sweet orange and growth was significantly inhibited, with the rare production of new shoots toward the end of the experiment (Figure 1). CLas was not detected by RT-qPCR in samples from any of the mockinoculated rough lemon or sweet orange throughout the entire experiment.
Comparative RNA-seq analysis A comparison of the RNA-seq data of the tolerant rough lemon and susceptible sweet orange was conducted using samples collected in week 7. RNA-seq data obtained at week 0 were used as a baseline. independently performed using a pair-wise comparison between week 7 and the baseline week 0 for inoculated and mockinoculated rough lemon and sweet orange using DEseq for each treatment. 22 A total of 33 930 unique transcripts were identified and quantified but 1884 genes did not have any corresponding  Table S3). Among the 269 DEGs in inoculated rough lemon samples, 76 were upregulated ⩾ 3-fold (log2), abbreviated as log-fold change (LFC), while no DEGs with this magnitude of downregulation were observed ( Figure 3). CLasinoculated samples of sweet orange contained 897 DEGs, of which 119 DEGs were upregulated ⩾ 3 LFC, while 11 DEGs were found to be upregulated ⩾ 3 LFC among the 445 DEGs identified in mockinoculated samples. Genes related to the defense response, such as those in the MAPK signaling pathway, WRKY transcription factors, and other pathogenesis-related genes were identified among the highly expressed (⩾4 LFC) genes in the CLasinoculated samples of rough lemon (Table 1). Upregulated genes greatly outnumbered downregulated genes in both CLasinoculated samples of sweet orange and rough lemon ( Figure 3). The number of downregulated DEGs was slightly higher than upregulated DEGs in mock-inoculated samples of sweet orange, while the number of upregulated DEGs was much greater than the number of downregulated DEGs in CLasinoculated samples of sweet orange. A total of 29 and 25% of the DEGs in samples of CLas-inoculated rough lemon and sweet orange had a LFC ⩾ 3, respectively, whereas only 1% of the DEGs in mock-inoculated samples of rough lemon and sweet orange were observed to have an LFC ⩾ 3. These results indicate that overall gene transcription levels were higher in CLas-inoculated Expression levels of 9 defense response-related genes were analyzed by RT-qPCR in order to confirm the results obtained by RNA-seq analysis. As shown in Figure 4, the relative levels of expression of the investigated genes were consistent between the RT-qPCR and RNA-seq data, with a few exceptions but within error range, indicating that the RNA-seq data were reliable.
Functional categorization and subnetwork enrichment analysis The Fisher Extraction Test and subnetwork enrichment analysis (SNEA), as provided in Blast2GO and Pathway Studio (Elsevier/ Ariadne Genomics), were used to determine the specific Gene Ontology (GO) terms affected by CLas infection. GO terms associated with hormones, defense response, as well as stomatal and photosynthetic acclimation, were significantly overrepresented based on the Fisher's pairwise comparisons of CLasinoculated and mock-inoculated genes in both rough lemon and sweet orange samples (Supplementary File 3, Supplementary Tables S1 and S2). Only 1 GO term (FDR ⩽ 0.05), response to chitin, was identified in sweet orange, while 19 GO terms were identified in rough lemon samples ( Table 2). Three of the top GO terms in rough lemon were response to chitin, ethylene biosynthetic process, and respiratory burst involved in defense response and response to molecules of bacterial origin ( Table 2). Although a lower number of DEGs were identified in rough lemon than in sweet orange, more defense response GOs were associated with the transcriptomic response of rough lemon to CLas infection than in sweet orange, and may reflect the 'tolerance' response observed in rough lemon. The results of the Fisher's Test indicated that 'ATP binding' was the most significantly enriched GO molecular function term in rough lemon, perhaps reflecting the significant number of genes under transcriptional regulation during the early response of rough lemon to CLas infection (Supplementary File 3, Supplementary Table S1). No statistically significant GO terms in the molecular function category were found in the sweet orange RNA-seq data (Supplementary File 3,  Supplementary Table S2). SNEA uses a global expression regulatory network extracted from the entire PubMed database and full-text journals to extract regulatory networks. Using the non-parametric Mann-Whitney test, SNEA identified significant (P ⩽ 0.05), over-represented (P ⩽ 0.05) ontologies among the upand downregulated DEGS that indicate an enrichment of each sub-network. A total of 32 and 44 ontologies associated with upregulated DEGs in CLas-inoculated rough lemon and sweet orange samples were identified, while 26 and 22 ontologies associated with downregulated DEGs were found in rough lemon and sweet orange, respectively (Supplementary File 4, Supplementary Tables S1 and S3). The top 3 ontologies for upregulated genes in CLas-inoculated samples of rough lemon were associated with defense response, biotic stress, and plant defense (Table 3a). In comparison, plant immunity, jasmonate response (JA), and shade avoidance were the top ontology terms identified in CLas-inoculated samples of sweet orange (Table 3b). Among downregulated DEGs, root phototropism, seed width, and cell expansion were the top 3 ontologies identified in CLasinoculated samples of rough lemon; whereas root phototropism, response to ethylene stimulus and greening were the top ontologies identified in sweet orange (Tables 3a and b, ). In regards to mock-inoculated samples of rough lemon, photosynthesis, reductive pentose-phosphate cycle, and seed length were the top three ontologies identified for upregulated DEGs. On the other hand, plant response, nitrogen metabolism, and biotic stress were the top ontologies identified for downregulated DEGs. With respect to mock-inoculation of sweet orange, cuticle development, response to auxin stimulus, and ion homeostasis were the top ontologies identified for upregulated genes, whereas meristem development, meristem initiation, and leaf size were the top ontologies for downregulated genes (Supplementary File 4, Supplementary Tables S2 and S4).
Transcription factors A total of 68 transcription factors were identified (Table 4). Eighteen transcription factors (TFs) were identified in CLasinoculated and mock-inoculated samples of rough lemon, all of which were related to plant immunity, and 12 of which had a transcription level ⩾ 3 LFC. Only one TF ⩾ 3 LFC was observed in mock-inoculated rough lemon. ERF-1 was the only TF found in both CLas-and mock-inoculated samples of rough lemon, and was upregulated in CLas-inoculated samples and downregulated in mock-inoculated samples of rough lemon. A total of 34 TFs were identified in CLas-inoculated samples of sweet orange and 7 in   (Table 4).
Highly-expressed DEGs DEGs involved in signal transduction were overwhelmingly observed and highly expressed in CLas-inoculated leaves of both rough lemon and sweet orange ( Table 1). The AT3G47570 gene, which encodes a protein that functions in protein phosphorylation in the transmembrane receptor protein tyrosine kinase signaling pathway, was the gene most hit

DISCUSSION
The current study examined transcriptomic changes associated with HLB during the early stage of CLas infection in HLB-tolerant rough lemon and HLB-susceptible sweet orange. Detection of CLas in citrus is delayed after initial inoculation; trees remain visibly asymptomatic until the bacterial titer reaches levels sufficient for reliable detection by RT-qPCR. 17 Infection of CLas-inoculated plants of rough lemon and sweet orange could not be confirmed by RT-qPCR in our study until 23 weeks after the plants had been inoculated. A previous study reported that HLB could be detected in sweet orange at 5-9 weeks after inoculation under greenhouse conditions. 10 The discrepancy in the time needed for confirmation in the two studies is likely associated with the use of different inoculum sources with different titers of bacteria, and possibly  differences in plant age and growing conditions. Stover et al. 26 found that inoculation with Kuharske citrange results in lower CLas titers than other inoculum sources such as lemon and sweet oranges. Since both Kuharske and Carrizo are citranges, the inoculum used in our study could take a longer time to be detected by qPCR than 'Lisbon' lemon in the earlier study. Identifying the proper time point to collect samples to characterize an early response to infection can also be problematic. In our study, bulk samples were collected at 0, 5, 7, 9, 12 and 17 weeks in order to determine the appropriate time point to represent an early response. Preliminary RNA-seq data indicated that the number of DEGs changing in expression was maximized at week 7; as a result, this time point was selected for detailed sequencing and analysis (data not shown). An appropriate statistical test is needed to determine whether or not observed differences in the number of obtained sequences (read counts) of a specific gene is significant, that is, whether the difference is greater than what would be expected due to natural random variation. 22 In addition, it is also critical to design the experiment so that natural variation (not due to a treatment effect) in gene expression is minimized and the number of identified DEGs is maximized. Previous studies directly comparing CLas-and mock-inoculated plants, or tolerant and susceptible genotypes, may have not sufficiently reduced biological variation after inoculation at week 0. As a result, a low number of earlystage DEGs were identified. 4,10,11,27 In the present study, gene expression in week 0 was used as a baseline in order to reduce the level of random variation between individual plants. To identify the DEGs induced by CLas, only DEGs that were uniquely expressed in CLas-inoculated and not in mock-inoculated samples at week 7 were selected for comparative analysis. Using this approach, we identified a greater number of DEGs in response to CLas inoculation than by direct pairwise comparisons of CLasinoculated rough lemon and sweet orange at week 7. This approach has not been used in previous transcriptome studies of CLas infected citrus. However, we would not identify as DEGs those genes which had high absolute levels of expression in both week 0 and week 7 after inoculation. We speculate that some of these genes might play important roles against CLas, although these genes were not necessarily induced by CLas.
Gene expression associated with tolerance and susceptibility to CLas Plant disease resistance is generally divided into two categories, qualitative resistance, determined by major R-genes, and quantitative resistance, determined by multiple genes with minor effects. 28 The first tier of defense, PTI, is triggered by the perception of PAMP/DAMPs by membrane-anchored PRRs, which is then followed by the activation of a MAPK cascade and downstream transcription factors, leading to immune responses. The second tier of defense is elicited by pathogen effectors via an interaction with an R protein (effector-triggered immunity), where the interaction between R protein and a pathogen effector oscillates between compatible and incompatible reactions over time. Plant pathogens are broadly divided into biotrophs and necrotrophs. Plant defense against biotrophic pathogens is largely due to major gene resistance. 29 CLas is an obligate biotrophic pathogen and currently no resistant citrus varieties has been identified. R-gene-mediated resistance usually induces a hypersensitive response (HR), which is thought to combat biotrophic pathogens by restricting their access to water and nutrients. 29 R gene-mediated resistance also activates SA-dependent signaling, leading to an activation of a string of presumed defense effector genes. The activation of SA signaling occurs throughout the plant to establish systemic acquired resistance (SAR) against subsequent pathogen infections. 29 During SAR, deposition of callose and lignin occurs in plant cell walls, and plants acquire the ability to mount a rapid HR. Fan et al. 9 reported that callose-plugged phloem sieve elements and inhibition of phloem loading were observed in CLasinoculated leaves of both rough lemon and sweet orange, but that phloem transport was less affected in rough lemon than in sweet orange. In our study, phloem protein 2 (PP2-15) was identified as one of the most highly upregulated genes in CLas-infected sweet orange, but not in CLas-infected rough lemon. Several other studies also reported that phloem protein 2 (PP2) genes are upregulated in HLB-infected leaves. 4,[9][10][11]14 In addition to being involved in the differentiation of vascular tissue, a lectin-like protein (PP2) is associated with the plugging of sieve plates in response to wounding and as a defense against pathogens and insects. 30 Accumulation of PP2 in sieve plates, in conjunction with phloem necrosis and blockage of the translocation stream, appears to be a major factor in disease symptom development of citrus greening. 11 PP2-B15 was the most highly upregulated gene in CLas-inoculated roots of 'Sanhu' red tangerine at 50 days post inoculation. Zhong et al. 14 speculated that PP2-like genes haves an active role in defense against the invading bacteria after a plant has been infected with CLas. In our study, SNEA specifically identified 'defense response' (P = 1.77E − 06) and 'plant immunity' (P = 0.0003) as overrepresented ontologies among the upregulated genes (Table 3a and 3b, ) in CLas-inoculated leaves of rough lemon and sweet orange, respectively. Ontologies for 'jasmonate response' and 'stomatal movement' were also significantly over-represented in CLas-inoculated rough lemon. The accumulation of callose, which is synthesized between the cell wall and the plasma membrane, as well as stomatal closure, are classic markers of PTI. 31 In addition, SA, JA, and ethylene hormones are induced during PTI. DEGs related to R gene-mediated resistance, which activates an SAdependent signaling pathway, were not found in our study. EDS1and PAD4, which play important roles in SA signaling, and NPR1, a master regulator of SA, 32 were also not significantly induced in any of the samples analyzed in the current study. Wang 33,41 indicated that CLas contains CLIBASIA_00255, which encodes a salicylate hydroxylase that can convert salicylic acid (SA) into catechol, a metabolite that does not induce a resistance response. 'AT1G73805', Systemic Acquired Resistance Deficient 1 (SARD1,), however, was highly expressed in CLas-inoculated samples of rough lemon (LFC = 4.35) and sweet orange (LFC = 9.00). SARD1, a plant-specific DNA-binding protein, is a key positive regulator of SA synthesis, and was induced after exposure to SA and JA. 34 Knocking out SARD1 compromises both basal resistance and SAR. 31 Further study is needed to understand the mechanism of SARD1 in citrus defense system. The SNEA analysis of the DEGs identified in the present study, specifically identified 'hypersensitive response' and 'cell death' as over-represented ontologies (P ⩽ 0.05) among upregulated DEGs in CLas-inoculated leaves of both rough lemon and sweet orange, but these ontologies were not identified in mock-inoculated samples of rough lemon and sweet orange (Supplementary File 4). A total of 4 and 7 DEGs associated with hypersensitive response and cell death were identified in CLas-inoculated leaves of rough lemon, respectively, whereas 8 and 16 DEGs associated with hypersensitive response and cell death were identified in CLas-inoculated leaves of sweet orange, respectively. Among the identified DEGs, BAP2, GLIP1 and ACL5, associated with inhibition of cell death, were found in CLasinoculated samples of sweet orange, whereas only BAP2 was found in CLas-inoculated samples of rough lemon. Arabidopsis C2 domain proteins, BAP1, and its homologue BAP2, negatively regulate biotic and abiotic cell death. 35 ACL5 controls Arabidopsis xylem morphogenesis through the prevention of premature cell death. 36 Treatment of Arabidopsis plants with GLIP1 protein systemically inhibited cell death in distant leaves inoculated with A. brassicicola, where cell death would be otherwise strongly induced. 37 MPK3, PLA2A, RBOHD, RPM1, SAUL1, SOBIR1, VTC2, and WRKY22, all of which are associated with hypersensitive response and cell death, were upregulated in CLas-inoculated leaves of sweet orange, whereas SOBIR1, MKK9, RBOHD, and MPK3, also associated with hypersensitive response and cell death, were upregulated in CLas-inoculated leaves of rough lemon.
Starch, which accumulates extensively in photosynthetic cells, as well as in phloem elements and vascular parenchyma cells of leaf blades and petioles, has also been observed to accumulate in the xylem parenchyma and phelloderm of HLB-affected 'Valencia' orange trees but not in control (non-infected) samples. 38 Albrecht and Bowman 4 found that glucose-6-phosphate/phosphate transporter (GPT2), which mediates the import of glucose-6-phosphate, an essential substrate for starch biosynthesis, was more highly upregulated in infected 'Cleopatra' mandarin than in sweet orange. In the present study, GPT2 expression was significantly induced to a high degree in infected rough lemon but not in sweet orange. A total of 19 enzymes belonging to the starch and sugar metabolism pathway were identified. Fourteen were identified in CLas-inoculated sweet orange and 5 in mockinoculated sweet orange, whereas 6 sugar and starch metabolism enzymes were identified in mock-inoculated rough lemon and 2 in CLas-inoculated rough lemon.  Tables S1 and S2). Although no symptoms of HLB were observed in inoculated leaves of rough lemon and sweet orange, a greater number of DEGs related to sugar and starch metabolism were found in inoculated leaves of sweet orange than in inoculated rough lemon leaves (Supplementary File 5, Supplementary Tables  S1 and S2). These data indicate that transcriptomic reprogramming in carbohydrate metabolism and defense is occurring in response to CLas even before symptoms are visibly observed or positive RT-qPCR results for HLB were found. Previous studies have observed that many genes involved in photosynthesis are repressed in response to HLB, most likely due to increased accumulation of sucrose/glucose levels in leaves. 10,11,39 In our study, SNEA analysis indicated that photosynthesis (P = 0.001), as an ontology term, was the most significant in DEGs obtained from mock-inoculated rough lemon, whereas defense response (P = 1.77 E − 06), as an ontology term, was the most significant term among the DEGs obtained for CLas-inoculated rough lemon (Supplementary File 4, Supplementary Table S2). These data suggest that a dramatic shift from photosynthesis-to defenserelated gene expression occurs in rough lemon as an early response to CLas. In previous proteome studies, it was found that accumulation of starch synthase increased, but the production of photosynthesis-related proteins decreased in infected rough lemon, 40 and most significantly, upregulated proteins of infected sweet orange were involved in stress/defense response. 39

Regulation of WRKY expression by MAPKs
The peptide flg22, a bacterial flagella protein recognized by most plants, activates a MAPK cascade which then leads the activation of WRKY-type transcription factors, key regulators of plant defense. 41 Flagella in CLas, however, have not been observed in any of the numerous electron micrographs of these bacteria infecting plants and psyllids. 2 The reduced genome of CLas and their transmission by psyllids may allow it to avoid PTI. CLas, however, still possesses 57 genes coding for products that function in cell envelope biogenesis and the outer membrane, including lipopolysaccharides, and most flagellar genes, which might function as PAMPs. 33 Elongation factor Tu (EF-Tu) is one of the most abundant bacterial proteins and is recognized as a PAMP by Arabidopsis. 42 The plant PRR for EF-Tu is the LRR-RLK EF-Tu receptor (EFR), which belongs to the same subfamily (LRRXII) as FLS2. 42 In this study, upregulation (LFC = 2.84) of EFR was only found in inoculated sweet orange. No expression of FLS2 was observed in CLas-inoculated leaves of either rough lemon or sweet orange. During PTI, activation of the MAPK cascade leads to the activation of WRKY-type transcription factors and other key regulators of plant immunity. 43 It is thought that the MAPK cascade regulates plant immunity through the activation of defense-related genes via direct phosphorylation of downstream transcription factors, such as WRKYs and ERFs. 44,45 MKK9 is upstream of MPK3, and MAPKKK19 is upstream of MKK9 ( Figure 5). In our study rough lemon exhibited a stronger MAPK response than sweet orange (Table 1). After the perception of flg22, MPK6 activates ethylene biosynthesis through the phosphorylation of ACS6. 46 ACS6 was highly upregulated in CLasinoculated rough lemon but not sweet orange. We also observed the strong upregulation of ERF1 (LFC = 3.18) and ERF9 (LFC = 3.31) in CLas-inoculated leaves of rough lemon (Tables 3a and b, Figure 5). In response to pathogen attack, ET and JA cooperate through transcriptional induction of ET response factor 1 (ERF1). 47 Activation of the MAPK cascade induces members of the WRKY family of transcription factors and defense-related genes 48 in tobacco. 49 For example, WRKY22 is activated via a MAPKcascade induced by flg22. 50 In our study, WRKY22 expression was only found in sweet orange, however, the upregulation was only moderate (LFC = 1.96). These data were much lower than the upregulation of group 1 WRKY (WRKY33, WRKY70 and WRKY40) TFs in rough lemon. Emerging evidence has indicated that group I WRKY transcription factors, which contain a conserved motif in the N-terminal region, are also activated by MAPK-dependent phosphorylation, underlining their importance in plant immunity. 51 Our study indicates that WRKY transcription factor genes may play an important role in the tolerance to HLB exhibited by rough lemon. A total of eight WRKYs were identified in our RNA-seq analysis; with five in CLas-inoculated rough lemon and three in sweet orange. The magnitude of upregulation of all the WRKY genes in CLas-inoculated rough lemon was highly significant (LFC ⩾ 3), while only WRKY41 approached an increase of LFC = 3.00 in CLas-inoculated sweet orange. Higashi et al. 52 indicated that WRKY41 may be a key regulator in the crosstalk between the salicylic acid and jasmonic acid pathways. 52 In our study, WRKY41 was upregulated in inoculated samples of both rough lemon and sweet orange. In response to mock-inoculation, however, only rough lemon exhibited downregulation (LFC = − 1.62) of WRKY42 and no significant difference was observed in sweet orange. The magnitude of the changes in the expression level of all of the WRKYs was higher in CLas-inoculated rough lemon than in CLasinoculated sweet orange (Tables 3a and b). In Arabidopsis, WRKY46 was specifically induced by salicylic acid (SA) and the biotrophic pathogen Pseudomonas syringae. 30,53 The WRKY33 transcription factor in Arabidopsis is required for resistance to necrotrophic fungal pathogens. 30 Previous studies have presented evidence indicating that WRKY70 regulates plant disease resistance, 53 and identified WRKY70 as a node of convergence for integrating salicylic acid (SA)-and jasmonic acid (JA)-mediated signaling events during plant response to bacterial pathogens. 30 Over-expression of WRKY70 in Arabidopsis resulted in the suppression of several JA responses, including expression of a subset of JA-responsive and Alternaria brassicicola responsive genes. SA-induced expression of WRKY46, WRKY70, and WRKY53 is mainly dependent on NPR1. 30 In our study, expression of NPR1 did not respond significantly to either mock-or CLas-inoculation, however, WRKY46 and WRKY70 were strongly upregulated in CLasinoculated rough lemon. Martinelli et al. 54 reported that several WRKY members were more highly upregulated in fully ripe fruits than in young leaves of Clas infected sweet orange. 54 ABA-induced activation of MAPKs in rough lemon MAPK cascades ( Figure 5) have also been shown to be involved in ABA signaling. 55 In the present study, MAK 3, MAPKKK19 and MKK9 were all significantly upregulated in CLas-inoculated rough lemon (Supplementary Files, Supplementary Figures S1 and S2), and may play a role in the ABA signaling pathway as suggested in a previous Arabidopsis study. Blast2GO also indicated strong over representation of an ABA signaling GO term in CLas-inoculated rough lemon (Supplementary File 6). MPK3 is activated by both H 2 O 2 and ABA in Arabidopsis seedlings, and overexpression of MPK3 increased ABA sensitivity to ABA-induced post-germination arrest of growth. 56 MPK9, which is preferentially expressed in guard cells, is also activated by ABA and has been shown to mediate ABA signaling in guard cells. 57 MPK3/MPK6 negatively regulates ABA signaling in plants. 58 The WRKY superfamily of TFs is the major regulator of plant defense and SA-mediated signaling, but also participates in ABA-mediated signaling. 59 AtWRKY40 has been reported to be a negative regulator of ABA signaling during seed germination and also interacts with AtWRKY18 and AtWRKY60 to inhibit the expression of stress-responsive genes. 53 In our study, WRKY40 was highly upregulated (LFC = 4.9) in CLas-inoculated rough lemon. In addition, it was the most highly expressed among all of the upregulated WRKYs in CLas-inoculated rough lemon, but was not expressed in inoculated sweet orange. An and Mou 32 stated that although ABA is known to play a crucial role in adaptation to abiotic stress, its role in biotic stress responses is less understood. In general, however, ABA is considered to be a negative regulator of disease resistance. 55 Rapid production of reactive oxygen species (ROS) is one of the early detectable events following ABA perception. 60 Plant recognition of PAMPs also induces the rapid and transient production of ROS in an oxidative burst. 55 Interestingly, components of an ABA-activated MAPK cascade are also activated by ROS, suggesting that ABA and ROS may converge at the MAPK level in regulating stomatal closure. 55 Previous studies of the hypersensitive response (HR) in pepper revealed that the activation of HR-induced cell death and defense responses is linked to signal transduction pathways that are coordinated by the activity of defense signaling molecules, such as ROS and plant defense hormones. 61 In our current study, PUB22 (respiration burst in defense) was highly upregulated in both CLas-inoculated rough lemon (LFC = 5.08) and sweet orange (LFC = 5.42). RBOHD (respiratory burst oxidase-D) was also upregulated in CLas-inoculated rough lemon (LFC = 3.15) and sweet orange (LFC = 2.88). ROS mediate ABA signaling in guard cells. 62 Although most of the links between ABA and MAPKs are poorly understood, it is evident that these pathways are part of the complex cellular signaling network that exists in plants to integrate various environmental cues, such as pathogen challenges, nutrient status or developmental programs. 55 AtRbohD and AtRbohF were initially described in Arabidopsis as key components of plant defense. RBOH-dependent ROS generation is associated with pathogen recognition during the oxidative burst linked to the perception of MAMPs/PAMPs and during the HR, which is coupled to the recognition of pathogen avirulence factors through resistance proteins. 37 Quantitative defense signaling In a review, Eulgem 16 stated that R-mediated resistance, basal defense, and SAR are related defense systems that share similar regulatory components and are effective against biotrophic pathogens. Differences between transcript profiles associated with R-mediated resistance and basal defense are quantitative rather than qualitative. 16 Our results indicate that defenseresponse-related DEGs may be mediated by both basal defense and R-mediated resistance mechanisms. Distinct sets of DEGs were induced in CLas-inoculated rough lemon and sweet orange. We observed that the overall global mRNA expression Reprogramming of a defense signaling pathway Q Yu et al.
(upregulated and downregulated) profiles in CLas-inoculated and mock-inoculated rough lemon and sweet orange were similar (Supplementary File 2, Supplementary Table S1). However, the amplitude of expression of the defense response in susceptible sweet orange was lower than in the tolerant rough lemon ( Figure 3, Table 5). Resistance is regulated by changes in the levels of response of a variety of defense mechanisms and by quantitative enhancements that make the defense response more effective. 15,16,63 The signaling mechanism that controls the activation of defense mechanisms consists of a highly interconnected network. Differences in the level of resistance often lie in the speed and intensity of specific defense responses. Resistant plants often respond more rapidly and vigorously to pathogens than susceptible plants. Hence, it is important to understand how plants sense the presence of a pathogen and initiate a subsequent response. Although the expression profile of up and downregulated genes in our study was similar in tolerant rough lemon and susceptible sweet orange, the proportion of downregulated DEGs with an LFC ⩽ − 3 was 0.0% and 5.8% in rough lemon and sweet orange, respectively, while the proportion of upregulated DEGs with an LFC ⩾ 3 was 29.0% and 25.1% in CLas-inoculated rough lemon and sweet orange, respectively. In mock-inoculated rough lemon and sweet orange, however, the proportion of DEGs that were down regulated was 3.8% and 1.2% in rough lemon and sweet orange, while the proportion of upregulated genes was 1.0% and 1.1%, respectively. RT-qPCR demonstrated that the transcription amplitude of 9 selected genes was higher in CLasinoculated than in mock-inoculated treatments. In general, these data were consistent with the expression levels obtained from the RNA-seq data. The LFC values of T29M8.8, AGP16, PUB22, CBF, AT1G32928, JAZ1, MPK3 and WRKY70 were also higher in CLasinoculated rough lemon than in the CLas-inoculated sweet orange. It should be noted, however, that there was a slight discrepancy between the expression levels of some genes in the RNA-seq data versus the RT-qPCR data, including DIC2, EXL2, and HSP90.1. Our results demonstrated that a higher number of DEGs with an LFC ⩾ 3 were present in CLas-inoculated rough lemon sampled at week 7 post-inoculation. Activation of resistance and defense responses is stronger and faster in CLas-inoculated rough lemon than in CLas-inoculated sweet orange. Mock-inoculated plants of both rough lemon and sweet orange exhibited fewer WRKYs and defense-related genes than CLas-inoculated plants.
These findings suggest that sweet orange exhibits an inability to suppress effectively the pathogen, resulting in the compatibility of the interaction between the host and pathogen. Maleck et al. 64 found that defined sets of up-regulated genes are more highly expressed in an accelerated manner during incompatible interactions. Tao et al. 15 proposed a mechanism that is commonly used in both R-gene and basal defense responses, where signal input to gene expression output occurs in a quantitatively determined manner. This mechanism is likely to involve regulators commonly used by both defense systems, such as WRKYs and defense response genes (Tables 4 and 5).

CONCLUSIONS
The current study presents a comparative transcriptional analysis of the response of rough lemon and sweet orange to mockinoculation and CLas-inoculation at 7 weeks post-inoculation in order to identify the early response to CLas infection in tolerant (rough lemon) and susceptible (sweet orange) trees. Rough lemon and sweet orange may recognize CLas via a pathogen recognition receptor, and activate plant immune systems even before positive RT-qPCR results can be obtained for CLas. Our collective results indicate that quantitative disease resistance may contribute to the durable tolerance to HLB exhibited by rough lemon, characterized by distinct transcriptional regulation of genes in various functional categories. Signaling pathways used by different defense systems appear to converge and target overlapping gene sets. The defense response DEGs mediated by basal defense and R-mediated resistance mechanisms were accelerated and amplified by CLas inoculation. Some transcription factors have been proven to play a pivotal role in disease resistance. Regulatory circuits linking signaling and gene regulation are emerging which suggest that a complex interplay between transcriptional activators and repressors fine-tunes the expression of the defense-related transcriptome. The current study provides new insights into the complex network of transcriptional regulation that occurs during the early response of rough lemon and sweet orange to CLas. Reprogramming of a defense signaling pathway Q Yu et al.