Wide-ranging transcriptomic analysis of Poncirus trifoliata, Citrus sunki, Citrus sinensis and contrasting hybrids reveals HLB tolerance mechanisms

Huanglongbing (HLB), caused mainly by ‘Candidatus Liberibacter asiaticus’ (CLas), is the most devastating citrus disease because all commercial species are susceptible. HLB tolerance has been observed in Poncirus trifoliata and their hybrids. A wide-ranging transcriptomic analysis using contrasting genotypes regarding HLB severity was performed to identify the genetic mechanism associated with tolerance to HLB. The genotypes included Citrus sinensis, Citrus sunki, Poncirus trifoliata and three distinct groups of hybrids obtained from crosses between C. sunki and P. trifoliata. According to bacterial titer and symptomatology studies, the hybrids were clustered as susceptible, tolerant and resistant to HLB. In P. trifoliata and resistant hybrids, genes related to specific pathways were differentially expressed, in contrast to C. sinensis, C. sunki and susceptible hybrids, where several pathways were reprogrammed in response to CLas. Notably, a genetic tolerance mechanism was associated with the downregulation of gibberellin (GA) synthesis and the induction of cell wall strengthening. These defense mechanisms were triggered by a class of receptor-related genes and the induction of WRKY transcription factors. These results led us to build a hypothetical model to understand the genetic mechanisms involved in HLB tolerance that can be used as target guidance to develop citrus varieties or rootstocks with potential resistance to HLB.

In general, the visual symptoms were more evident in the susceptible plants, while the visual HLB symptoms were undefined in P. trifoliata and its hybrids. However, according to CLas detection, starch and callose quantification between different treatments, the hybrids were clustered into three distinct groups as follows: Susceptible Pool (S Pool), composed of three different hybrids (H109, H161, and H165) that were diagnosed as HLB-positive and presented elevated starch and callose deposition, similar to that observed for susceptible parental genotypes (Fig. 1); Tolerant Pool (T Pool), composed of three different hybrids (H113, H154, and H146) that were diagnosed as HLB-positive but did not exhibit a significant starch and callose accumulation as observed in susceptible genotypes (Fig. 1); and Resistant Pool (R Pool), composed of three different hybrids (H68, H106, and H142) that were diagnosed as HLB-negative with starch quantification similar to healthy plants (mock-inoculated plants) (Fig. 1).

Transcriptome assembly.
To elucidate the different responses to CLas infection, we studied the changes in global transcriptional level in susceptible, tolerant, and resistant genotypes infected by CLas. In this work, 36 cDNA libraries from six different genotypes of either CLas-inoculated or mock-inoculated (control) samples were evaluated. After trimming, 487 million reads were obtained, and 95% of the total was assigned (see Supple- Table 1. Detection and quantification of the bacteria by quantitative PCR (qPCR) in Citrus sinensis, C. sunki, Poncirus trifoliata and nine hybrids from an F 1 population obtained from the cross between C. sunki and P. trifoliata Raf. cv Rubidoux. Each individual is represented by five repetitions. -: Correspond to samples that were not evaluated due to the absence of leaves.  Table S1). The reads were mapped in 133,976 transcripts on the C. sinensis genome available on http:// citru s.hzau.edu.cn/. HLB-susceptible genotypes, C. sinensis and C. sunki, showed a high number of differentially expressed genes (6141 and 5624 DGEs, respectively) compared with the tolerant parental, P. trifoliata (100 DEGs) ( Table 3). A similar pattern was observed between the pool of hybrids. The S Pool showed 708 differentially expressed genes (DEGs), while the R Pool presented only 92 DGEs. The Tolerant Pool (T Pool) showed the highest number of DEGs (2027) among the hybrid pools. Most of these genes were downregulated in HLB-infected plants compared with healthy ones ( Table 3).
The principal component analysis (PCA) using the Bioconductor package (see Supplementary Fig. S1) showed the replicates of the different genotypes in general grouped according to the analyzed condition for C. sunki, C. sinensis, and the susceptible and tolerant hybrids. The resistant groups in fact presented a mixed grouping, which is not surprising if we consider that these hybrids were the ones that showed the fewer number of DEGs. The genotype grouping indicated that the global expression landscape is related more to the different genotypes and not the analyzed condition (infection by CLas). In this case, the identification of genes exclusively differentially expressed in the genotypes considered susceptible, tolerant, or resistant as well as genes that had antagonistic expression between the opposite phenotypes became important to increase our understanding of the different responses. Table 2. Candidatus Liberibacter asiaticus (CLas) quantification obtained by comparing the standard curve of the HLB primers with the standard curve of the internal control gene (GAPDH) initiators. The value quantification refers to Log 10 of the number of copies of the CLas fragment after 240 days from inoculation in each repetition per genotype included in RNAseq analysis.

Genotype
Ct value of GAPDH Ct value of HLB Quantification/Log 10 number of copies  Table S3 and S4). Among the seven genes upregulated in the R Pool, five were downregulated in C. sinensis, and one gene was downregulated in the T Pool and another one in C. sunki (see Supplementary Table S5). The study of genes with antagonistic expression between susceptible and tolerant and/or resistant genotypes may help to explain possible tolerance mechanisms as well as to identify good targets for plant resistance.

Main processes affected by CLas infection. Libraries of DEG functions assigned by Blast2GO 9 and
Gene ontology (GO) 10 analyses helped us better understand the differences in genetic responses involved in susceptibility, tolerance, or resistance (Fig. 3). Susceptible genotypes and tolerant hybrids differentially expressed many genes in comparison to resistant hybrids and P. trifoliata. These different pathways provided valuable information regarding the genetic mechanisms of CLas perception and responses activated in tolerant/resistant and susceptible hosts (Fig. 3).
Differentially expressed genes (DEGs) associated with a specific biological pathway. Signaling receptor. Plant receptors are responsible for the recognition of several external stimuli, including pathogen attack. These transmembrane proteins are directly associated with signaling pathways, which trigger a proper physiological response 11 . Several types of receptors were regulated in C. sinensis, C. sunki, the S Pool, and the T Pool, and most of them were downregulated in those genotypes (Fig. 3). In P. trifoliata and the R Pool only a few receptors were differentially expressed, and most of them were induced (Fig. 3). These receptors included G-type lectin S-receptor-like, cysteine-rich receptor kinase, and serine/threonine-protein kinase, which were upregulated in P. trifoliata, and leucine-rich repeat transmembrane kinase and leucine-rich repeat receptor-like protein kinase, which were induced in the R Pool (see Supplementary Table S6). Therefore, our results suggests that downregulation of receptors may be associated with susceptible response to CLas.
Hormones. Genes associated with auxin and ethylene pathways were barely or not affected in P. trifoliata and the R Pool, whereas many auxin and ethylene-related genes were differentially expressed in C. sunki, C. sinensis, the T Pool, and the S Pool under CLas infection. Interestingly, no important changes in the transcriptional profiles of genes related to SA and JA biosynthesis were found (Fig. 3). In addition, CLas induced key genes involved with gibberellin (GA) degradation in tolerant and resistant genotypes, while the related GA synthesis genes were downregulated. In P. trifoliata, the gibberellin-induced gene was one of the top three downregulated DEGs (log2 fold change = − 10) (see Supplementary Table S6). The opposite pattern was observed in CLas-susceptible genotypes, in which an induction of genes involved with GA synthesis and downregulation of GA degradation was observed. Thus, these findings suggests that GA plays an important role in CLas-citrus interactions, affecting plant physiology and consequently HLB symptoms.
Transcription factors. Plant responses to pathogen attack require large-scale transcriptional reprogramming. P. trifoliata showed only five transcription factor (TF)-related genes modulated by CLas infection. Only the MYB TF was downregulated. The other four TFs were upregulated, including two WRKY TFs (Fig. 3). The resistant hybrids suppressed the expression of another class of transcription factor, the SCL domain (see Supplementary  Table S6). In contrast, hundreds of TF genes showed changes at the transcription level in C. sinensis, C. sunki, the S Pool and the T Pool (Fig. 3). In this context, the large number of TFs affected in these genotypes may be directly related to the regulation of genes responsive to HLB infection. Of note, several WRKY TFs were identified in C. sinensis and C. sunki, and most of them were repressed in CLas-infected plants (Fig. 3). Therefore, these results The bar graph next to the microscopy plates show the callose quantification performed by counting fluorescent spots marked by aniline blue dye. Quantification was performed with tree replicates per genotype, inoculated plants (positive or negative HLB) and mock-inoculated plants. (c) Starch quantification. Individuals were inoculated with CLas (CLas-infected) or mock-inoculated (CLas-free) and collection was performed after 240 days, and quantification was carried by the enzymatic method. Bars represent the standard deviation between 3 biological replicates. * p_value < 0.05 (Mock-inoculated × CLas inoculated). www.nature.com/scientificreports/ indicated that the increase in transcription of WRKY TFs in P. trifoliata is associated with the genetic defense mechanism involved with HLB tolerance.
Defense-related genes. Defense-related genes are directly related to processes or production of compounds able to inhibit pathogen reproduction or to make further infection more difficult 12 . In particular, one defense-related gene, endochitinase B, was differentially expressed and highly upregulated in resistant hybrids (see Supplementary Table S6). Endochitinases have previously been reported as important bactericides, and some of them have ability to cleave peptidoglycan chains, promoting bacterial cell lysis 13 . Other defense-related genes were differentially expressed in susceptible plants by CLas. Among them, regions encoding lipid transfer, molecular factors that help the innate immune system of plants, and small lipid-transfer proteins can inhibit fungal growth and pathogenic bacteria 14 . Genes encoding these proteins were differentially expressed in C. sinensis, C. sunki, and the S Pool (see Supplementary Table S6 and Fig. S2). These results indicated the activation of defense pathways in response to CLas infection in susceptible genotypes. CDR1 also represents an important defense related gene in Poncirus and Poncirus-hybrids 15 . CDR1 showed high expression in all the Poncirus hybrids, including the S pool, but it was only induced in the R pool. Therefore, www.nature.com/scientificreports/ even though it could be associated with resistance, high CDR1 constitutive expression level seems not to be sufficient to lead to the resistance phenotype.
Secondary metabolism and cell wall composition. Secondary metabolites often play an important role in many physiological responses, such as growth, photosynthesis, reproduction, and plant defenses against pathogens 16 . The most upregulated genes in P. trifoliata included a variety of phenylpropanoids and lignin-related genes, such as caffeic acid O-methyltransferase, chalcone synthase, feruloyl ortho-hydroxylase 1, hydroxycinnamoyl transferase and laccase precursor (see Supplementary Fig. S2). In our study, the laccase precursor gene, whose protein catalyzes lignin and its derivatives 17 , was exclusive and highly induced in CLas-infected P. trifoliata (see Supplementary Table S6). Pectin hydrolysis occurs frequently in response to bacterial infection 18 . Just one pectin degradation-related gene was differentially expressed (downregulated) in P. trifoliata ( Fig. 3 and Supplementary Table S6). Many genes involved in pectin synthesis and degradation were differentially expressed in C. sinensis and C. sunki. Pectin methyltransferases are enzymes that induce pectin modification. In C. sinensis and C. sunki under stress caused by CLas infection, the pectin methyltransferase 1 gene was upregulated (see Supplementary Fig. S2).
A larger number of DEGs involved in cellulose synthesis showed mRNA levels altered in susceptible genotypes; however, P. trifoliata and the R Pool did not exhibit differentially expressed regions encoding cellulose (see Supplementary Table S6).
These results demonstrated that the cell wall is highly affected in susceptible plants even at 240 days after CLas inoculation. At the same time, genes involved in cell strengthening proved to be important in P. trifoliata.
Phloem-related genes. It is already known that callose deposition and phloem proteins (PP2) act as a physical barrier, attempting to block systemic spread of CLas; however, they also likely cause phloem disorders 19 . The current study identified DEGs coding phloem proteins that had altered expression induced by CLas in C. sinensis, C. sunki, the S Pool, and T Pool. Although P. trifoliata did not present callose-induced phloem blockage (Fig. 1), we observed modulation of PP2-B15 in response to CLas with ninefold higher expression than the control (see Supplementary Table S6). That result suggests that P. trifoliata modulates phloem genes in response to CLas without over-deposition of callose, consequently not causing important phloem function disorders. Anatomical divergences between P. trifoliata and Citrus may represent an important feature to avoid collapse of the sieve tube elements 20 . Genes are classified into nine groups (Stress response, Transporter, Carbohydrate metabolic process, Cell wall, Phenylpropanoids, Immune response, Transcription Factors Hormones and Signaling receptors) according to Blast2GO analysis and based on their expression pattern. The number of down-regulated genes in response to CLas is represented by the bars in reddish tones and upregulated in blue tones. Some bars present subdivisions and the color legend for each pathway is indicating the specific, related gene or specific pathways, which were important to illustrate the proposed tolerance mechanism to HLB. www.nature.com/scientificreports/ As shown by our phenotypic data, only susceptible plants had affected callose deposition. Different callose synthases were differentially expressed in the susceptible plants, whereas those genes were absent in P. trifoliata and the R Pool inoculated with CLas (see Supplementary Table S7).
Interestingly, genes encoding sieve element occlusion c (SEOc) and d (SEOd), which are part of a protein family that encodes specialized crystalloid phloem proteins 21 , were largely upregulated in all susceptible plants under study. Some of these genes were also upregulated in tolerant hybrids (see Supplementary Table S6 and S7).
Carbohydrate metabolism. Carbohydrate metabolism was the biological function most affected by HLB (Fig. 3). In the presence of CLas, susceptible genotypes overexpressed genes involved with starch synthesis and suppressed genes that encode enzymes for starch degradation (see Supplementary Table S8 and S9). This phenomenon was not observed for the tolerant and resistant genotypes. Several DEGs involved in the metabolism of starch were identified in C. sinensis, C. sunki, and the T Pool, especially in the former two (see Supplementary  Table S6). Genes encoding ADP-glucose pyrophosphorylase and starch branching enzyme II, which participate in the synthesis of starch and starch granules, were upregulated in C. sinensis and C. sunki (see Supplementary  Table S6 and S8). Beta and alpha-amylase, important enzymes for normal degradation of the starch in plants 22 , also had their genes expression modulated in both susceptible plants (C. sinensis and C. sunki) and the T Pool (see Supplementary Table S9). Corroborating our phenotypic data (Fig. 1), resistant and tolerant genotypes did not exhibit altered expression of the main genes involved in synthesis of starch (see Supplementary Table S6). While the R Pool had only beta-amylase-encoding gene upregulated, P. trifoliata did not have any DEGs related to synthesis and degradation of starch (see Supplementary Table S6).
Transporters. The transport of substances was also one of the main biological functions affected by CLas. The transcription levels of genes related to transporters were overwhelmingly altered by CLas infection in all genotypes and hybrids (Fig. 3). In general, susceptible plants had the greatest number of transport-related genes affected by CLas (Fig. 3). The R Pool showed few DEGs related to transport function, including ABC transporter family, phosphate transporter (PHO1-2), and amino acid transmembrane transport (Supplementary Table 2). Zinc transporter (ZIP1 and ZIP8) genes were differentially expressed in C. sinensis, C. sunki, and the T Pool (see Supplementary Table S6). Most transport family genes affected by CLas infection were involved with transport of sugars, amino acids, and ions (see Supplementary Table S6). When comparing the transporter-related DEGs in the tolerant genotypes, P. trifoliata, and the T Pool, we observed different responses among them. The T Pool exhibited 73 differentially expressed transporter-related genes. The parental P. trifoliata showed only 4 differentially expressed transporter-related genes, among which potassium transporter was exclusively differentially expressed in P. trifoliata (Fig. 3 and Supplementary Table S6).

Discussion
The hybrids evaluated in this work and the parents, Citrus sunki and P. trifoliata, were classified as susceptible, tolerant, or resistant according to bacterial presence, callose deposition, and starch accumulation (Fig. 1). RNAseq data indicated that the genotypes responded differently under CLas infection, which was confirmed by RT-qPCR analysis. Overall, the genes showed similar patterns in the RNAseq and RT-qPCR data, but some divergent values were found, which was similar to other transcriptome studies when the results of different techniques were compared 23 .
Our findings indicated that few genes were differentially expressed according to RNAseq analysis of the tolerant and resistant plants. In contrast, RNAseq analysis of susceptible plants showed transcription modulation of many genes. Resistant and tolerant plants have a tendency to respond more rapidly and vigorously to a pathogen than susceptible plants 12 . It is possible that the resistant hybrids have an early response to CLas presence. Early molecular interactions are well-known mechanisms in plant-pathogen interactions [24][25][26] . Nevertheless, to verify that the genetic responses were due to CLas infection and to avoid false positives, the samples for transcriptomic analysis were collected 8 months after CLas infection.
P. trifoliata showed upregulation of receptor-related genes, which presented an efficient recognition of CLas and possibly an effective signaling and activation of defense response against CLas. The reprogramming of defense signaling pathways has previously been reported as a critical element of the early response to CLas in tolerant genotypes 27 , such as P. trifoliata. Previous studies have also highlighted the induction of phenylpropanoidrelated genes as a molecular mechanism of HLB tolerance 5 . Lignin-related genes and several phenylpropanoids were strongly upregulated in P. trifoliata transcriptome (Supplementary Table S6). As reorganization of plant growth and development are critical to maximize plant survival under stress 28 , cell wall reinforcement is a tolerance mechanism of P. trifoliata against CLas. When comparing P. trifoliata and resistant hybrids, we observed a distinct transcriptional response to CLas (Fig. 2). However, all replicates of the resistant hybrids did not present any detection of CLas, even after almost 1 year of the experiment (Tables 1 and 2), and probably for this reason, they exhibited few DEGs in RNAseq. Interestingly, the exclusive DEGs of the R Pool, formed by the CLas-negative hybrids, may be linked with genes and mechanisms capable of eliminating the bacteria from the plant, such as endochitinase B. Plant endochitinases cleave peptidoglycan chains, thereby promoting bacterial cell lysis 13 .
CLas infection is erratic and unpredictable, and even susceptible plants can escape from infection. Until almost 1 year, all plant replicates classified as resistant did not present CLas titer (Tables 1 and 2). Therefore, until that moment, we considered that those plants were resistant to CLas infection and that a mechanism was utilized to avoid spreading the disease.
In the transcriptome of tolerant genotypes, downregulation of GA synthesis genes and upregulation of genes involved with GA degradation were observed, and the opposite behavior was observed in the susceptible genotypes (induction of GA synthesis and repression of GA degradation). In addition, we observed upregulation of  Table S6). It is known that the GA pathway presents cross-talk with auxin and ethylene hormones, which are plant growth regulators that also have been associated with plant defense and microbial pathogenesis 29,30 . The present study showed that these regulators were strongly differentially expressed in the tolerant plants by CLas. It has been reported that auxin induces GA biosynthesis and suppresses GA degradation through modulation of several transcription factors and transporters 31,32 . In citrus-pathogen interactions, crosstalk between auxin and GA has also been reported. Inhibition of GA synthesis promotes inhibition of auxin-induced transcription, consequently reducing symptoms in the citrus-Xanthomonas citri interaction 33 . The plant tolerance mechanism is better explained by the interaction of GA and the salicylic acid (SA) hormone. The GA pathway is considered a hormone modulator of the SA signaling backbone during plant responses to pathogens [34][35][36] . In Arabidopsis thaliana, Alonso-Ramírez et al. (2009) 36 showed that GAs and the overexpression of GA-responsive genes increase not only the endogenous levels of SA but also the expression of ics1 and npr1 genes involved in SA biosynthesis and action, respectively. However, SA-related genes were almost not modulated in the present study, which might be due to the high SA level in the evaluated stage, resulting in the expression of SA synthesis-related genes no longer being necessary as shown by Oliveira et al., 2019 20 . Moreover, it is known that SA accumulation and downstream signaling events are important components of both pathogen-associated molecular pattern (PAMP)-triggered immunity (PTI) and effector-triggered immunity (ETI) 37,38 through increasing the expression of WRKY transcription factors. Many WRKY TFs were induced in the tolerant genotypes and affected in the susceptible plants (Fig. 5). WRKY TFs have been considered key regulators of plant defense against many pathogens, including CLas 27 . The function of some WRKY genes remains unexplored, but in some crop species, specific WRKYs promote tolerance or even resistance to biotic and abiotic stresses 27 . Thus, the induction of WRKY TFs may also be related to the activation of genes involved with the tolerance mechanism. For example, in P. trifoliata, the WRKY transcription factor 14-1 was induced, and its orthologue in Arabidopsis (known as WRKY22) is an essential component of MAPK-mediated plant defense responses against pathogens. MAPKs are associated with one of the earliest signaling events after plant sensing of PAMPs and pathogen effectors.
Moreover, the tolerant and susceptible genotypes had changes in the level of transcription of many callose synthases and phloem protein (PP2) genes in response to CLas infection (Supplementary Table S6 and S7). In addition, all susceptible plants showed induction of a class of genes that includes the SEOc gene (Supplementary  Table S7). This class of genes has been reported to encode P-protein subunits 21 . Overexpression of these genes could increase callose and PP2 protein synthesis in the citrus phloem sieve elements. Callose and PP2 accumulation is a crucial factor of phloem blockage in CLas-infected plants 19,39,40 . Phloem blockage causes disturbance of photoassimilate flows from source organs (leaves) to sink organs (roots), resulting in starch accumulation in the leaves as observed in this work and in previous studies 41 .
Based on the knowledge of CLas-susceptible plant interaction that culminates in HLB symptoms, a zig-zag model as illustrated previously by Jones & Dang (2006) 42 was adapted to explain such genetic molecular response to CLas (Fig. 4). During the beginning of infection, receptors from citrus plants detect the CLas PAMPs, which triggers a PTI response, resulting in the production of GA and SA as well as in the induction of several downstream genes (asymptomatic stage). In a second phase, CLas delivers effectors, such as Las5315 43 and others 44 , which interfere with PTI or enable pathogen nutrition and dispersal, resulting in effector-triggered susceptibility (ETS). In phase 3, effectors activate an ETI and an amplified version of PTI leading to induction of callose synthases and pp2 gene expression that results in callose and PP2 accumulation. Therefore, callose and PP2 accumulation and the consequent anatomical alterations of the sieve pores may lead to hypersensitive cell death (HR) of the infected plants, which spatially isolate the CLas to reduce their colonizing ability via the phloem 19,40 .
To describe the genetic mechanisms potentially involved in a susceptible, tolerant, and resistant interaction with CLas based on the data obtained in this study, we built a hypothetical model (Fig. 5). The model shows that in the susceptible plants (Fig. 5), auxin-related genes positively modulate GA synthesis, which activates response mechanisms to CLas infection, such as callose deposition, PP2 deposition, phloem dysfunction, and impaired flow transport. The impaired flow results in starch accumulation on mesophyll chloroplasts, which promotes thylakoid rupture and chlorophyll degradation, culminating in HLB typical symptoms. In the tolerant plants, including P. trifoliata (Fig. 5), the induction of signaling receptors cause a fast and efficient defense response modulated by suppression of the auxin pathway and induction of GA degradation. The suppression of these pathways prevents the events that lead to phloem dysfunction (callose deposition, starch accumulation, and transport alteration), and it activates the defense response through the synthesis of phenylpropanoids and cell wall strengthened-related genes. This transcriptional reprograming is efficient to impair the development of symptoms. In the resistant genotypes (Fig. 5), a potentially early and rapid defense may occur in response to CLas because only a few genes were differentially expressed after 240 days after inoculation. However, this response is related to induction of signaling receptors and upregulation of endochitinase B, which is associated with bacterial cell lysis.
Both hypothetical models showed that there are many pathways acting in citrus defense against CLas infection. The data acquired in this study can help to generate citrus varieties of scions or rootstocks with potential resistance to HLB based on citrus conventional breeding programs or biotechnological approaches, including the development of transgenic or cisgenic lines as well as genome editing and host-induced gene silencing.

Materials and methods
Plant material. C. sinensis, C. sunki, P. trifoliata, and 21 hybrids obtained from a controlled cross between Citrus sunki ex Tan (female parent and susceptible to HLB) and Poncirus trifoliata Raf. cv Rubidoux (male parent and tolerant to HLB) were used in the analysis. C. sinensis was included because it is one of most important citrus scions in the world, and it can also be considered an internal control of the experiment considering that Scientific Reports | (2020) 10:20865 | https://doi.org/10.1038/s41598-020-77840-2 www.nature.com/scientificreports/ C. sinensis is characterized as a species highly susceptible to HLB 40 . The experimental design was completely randomized and consisted of five biological replicates for each inoculated genotype (CLas-infected budwoods) and mock-inoculated genotype (health budwoods). Plants were propagated using buds that were grafted onto rootstocks of Rangpur lime (C. limonia Osb.). At the end of 6 months, the plant scions were grafted using two CLas-infected buds obtained from C. sinensis (L.) Osbeck cv Pera. All plants were kept in a greenhouse at Centro de Citricultura Sylvio Moreira of the Agronomic Institute (IAC), SP with an average temperature of 25 °C for 12 months. The starch content and callose deposition were estimated only in the genotypes selected for the further analysis (C. sinensis, C. sunki, P. trifoliata, and 15 hybrids obtained from crosses between C. sunki and P. trifoliata). Leaves from inoculated and mock-inoculated plants from all evaluated genotypes were collected with three biological replicates of each genotype after 8 months of CLas infection.

CLas quantification.
CLas presence and HLB symptoms were evaluated according to previously described methodology 40 . Briefly, 30, 90, 180, 240, and 360 days after inoculation, to confirm HLB infection, leaves above the inoculation point were collected and tested by qPCR using 16S ribosomal DNA primer sets and FAM/Iowa Black FQ label probe (IDT Inc., Coralville, IA) probes as described by Li et al. (2006) 45

Phenotypic analysis. Starch and callose quantification of CLas-inoculated and mock-inoculated plants
was performed after 240 days of infection. Callose quantification was performed following the methodology reported previously 40 . Leaf petioles were fixed in FAA solution (50 mL of formaldehyde, 50 mL of glacial acetic acid, and 900 mL of 70% ethanol) for 72 h and then kept in 70% ethanol. Transversal sections of 10 μm were generated using an automatic slide microtome (Leica SM2010R). The sections were stained with blue aniline, and the stained samples were examined on an Olympus BX61 fluorescence microscope using 355-375 nm excitation filter, 400-nm dichromatic mirror, and 435-490 nm emission filter. Callose quantification was performed by counting fluorescent spots in the total phloem area in 10 fields of view for each sample. The starch measurement was performed using leaves dried in an oven at 60 °C for 48 h and ground. Starch content was estimated by enzymatic analysis using 10 mg of dried leaves according to 46 . Absorbance was measured in 96-well microtiter

RNA extraction and sequencing (RNA-seq).
Leaves from three biological replicates of the three genotypes (C. sinensis, C. sunki, and P. trifoliata) and the three hybrid pools (S Pool: H109, H161, and H165; T Pool: H113, H154, and H146; and R Pool: H68, H106, and H142), either CLas-infected inoculated or mock-inoculated plants, were collected for transcriptomic analysis after 240 days of infection. It is difficult to establish the ideal time for studying the first responses and stages of infection because it is difficult to confirm that the plant tissue is colonized by bacteria. Thus, to verify that the genetic responses were due to CLas infection, we performed RNAseq analysis at 8 months. Total RNA was isolated with the MasterPure Plant RNA Purification Kit (EPICENTRE Biotechnologies, Madison, WI, USA) according to the manufacturer's instructions. A total of 10 µg of RNA from each sample was sent for sequencing at the Centro de Genômica Funcional in Centro de Biotecnologia Agricola in ESALQ/USP (http://www.esalq .usp.br/genom icafu ncion al/). RNA-seq was performed using the Illumina HiSeq 2500 platform. All procedures were performed according to Illumina's protocols. RNA-seq was performed in triplicate with a total of 36 samples. Interference on substance transport along with callose deposition causes phloem dysfunction resulting in starch accumulation on photosynthetic tissues. Starch accumulation promotes thylakoid rupture and chlorophyll degradation culminating in HLB classical symptoms. Tolerant plants, the induction of signaling receptors causes a fast and efficient defense response modulated by suppression of auxin pathway and induction of GA degradation. The suppression of these pathways prevents the events that lead to the phloem dysfunction (callose deposition, starch accumulation and transport alteration) and activates defense response through the synthesis of phenylpropanoids and cell wall-strengthened related genes. This transcriptional reprograming is efficient to impair the development of symptoms. Resistant genotypes, a possibly early and fast defense may occur in response to CLas, since low numbers of the genes are modulated after 240 days post inoculation. Nonetheless, this response is related to induction of signaling receptors and upregulation of Endochitinase B, which might be associated with bacterial cell lysis (Created with BioRender.com).
Real time PCR (RT-qPCR) validation. To ensure reproducibility of the biological phenomenon observed by transcriptomic analysis, we performed a second experiment with other plants following the same design used for RNA-seq. We sampled one hybrid of each pool to represent the susceptible, tolerant, and resistant pools. We used only one hybrid from each pool because it represents the hybrids that comprise each pool regarding CLas infection behavior. Total RNA was extracted using the protocol described by Chang et al. (1993) 50 . Traces of genomic DNA were eliminated using the DNase RNase-Free Ket (QIAGEN, Valencia, CA, USA) according to the manufacturer's instructions. cDNAs were synthesized from 1.0 μg of total RNA using Superscript III (200 U/μL) (INVITROGEN) with an oligo (dT) primer (dT12-18, INVITROGEN) according to the manufacturer's instructions.
Ten genes that showed the opposite expression profile between the genotypes with different responses were selected, including chalcone synthase, lipid transfer, cytochrome P450, gibberellin-regulated 9, sieve element occlusion c, cinnamoyl-reductase, pectin methylesterase 1, starch branching enzyme II, PRR response regulator, and choline transporter like-protein 2 (see Supplementary Table S10). Primers were designed using Primer3Plus 51 , and the Primer-BLAST tool 52 was used to check the specificity of the primers. Two endogenous genes, GAPDH and FBOX, were used for normalization of the data. Relative gene expression was calculated with the 2 −ΔΔCt method 53 . www.nature.com/scientificreports/ Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creat iveco mmons .org/licen ses/by/4.0/.