Reference Genes for Accurate Transcript Normalization in Citrus Genotypes under Different Experimental Conditions

Real-time reverse transcription PCR (RT-qPCR) has emerged as an accurate and widely used technique for expression profiling of selected genes. However, obtaining reliable measurements depends on the selection of appropriate reference genes for gene expression normalization. The aim of this work was to assess the expression stability of 15 candidate genes to determine which set of reference genes is best suited for transcript normalization in citrus in different tissues and organs and leaves challenged with five pathogens (Alternaria alternata, Phytophthora parasitica, Xylella fastidiosa and Candidatus Liberibacter asiaticus). We tested traditional genes used for transcript normalization in citrus and orthologs of Arabidopsis thaliana genes described as superior reference genes based on transcriptome data. geNorm and NormFinder algorithms were used to find the best reference genes to normalize all samples and conditions tested. Additionally, each biotic stress was individually analyzed by geNorm. In general, FBOX (encoding a member of the F-box family) and GAPC2 (GAPDH) was the most stable candidate gene set assessed under the different conditions and subsets tested, while CYP (cyclophilin), TUB (tubulin) and CtP (cathepsin) were the least stably expressed genes found. Validation of the best suitable reference genes for normalizing the expression level of the WRKY70 transcription factor in leaves infected with Candidatus Liberibacter asiaticus showed that arbitrary use of reference genes without previous testing could lead to misinterpretation of data. Our results revealed FBOX, SAND (a SAND family protein), GAPC2 and UPL7 (ubiquitin protein ligase 7) to be superior reference genes, and we recommend their use in studies of gene expression in citrus species and relatives. This work constitutes the first systematic analysis for the selection of superior reference genes for transcript normalization in different citrus organs and under biotic stress.


Introduction
Real-time reverse transcription PCR (RT-qPCR) has emerged as the most widely used method to quantify changes in gene expression profiles in response to developmental transitions and environmental changes in plants. In comparison to classical methods used to measure transcript abundance, the main advantages of RT-qPCR are its higher sensitivity and specificity, even when limited amounts of RNA are available [1]. Nevertheless, to ensure reproducible and accurate quantitative expression measures, it is necessary to normalize the expression levels of target genes using suitable reference genes. An ideal reference gene should be stably expressed among samples, including those from different tissues and cell types, developmental stages, and treatment conditions [2][3][4][5]. Because there is no gene that meets all requirements for every experimental condition, a systematic validation of the stability of expression of candidate reference genes should be conducted in preliminary experiments assessing their usefulness for gene expression normalization [2,6]. Gene expression analysis in citrus in different tissues and organs and under several experimental conditions has relied on the use of traditional housekeeping genes, such as ACTIN [7][8][9][10][11][12]; EF1-a [13][14][15][16]; TUBULIN [17] GAPDH [18], and 18S rRNA [19] as reference genes, but with no previous testing of the stability of expression. It is generally assumed that housekeeping genes encoding proteins required for basal cell activities, such as central carbon metabolism, protein translation, cytoskeleton maintenance, and protein turnover, are expressed uniformly in different tissues and organs [6]. However, under many conditions, the level of transcript expressed from such genes was not stable, which may have led to the misinterpretation of results [20][21][22][23][24]. Statistical algorithms such as geNorm [25] and NormFinder [26] have been recently used to identify the best reference genes for RT-qPCR data normalization in a given set of biological samples. These algorithms have been used for assessing the expression stability of candidate reference genes across a variety of tissues and organs, developmental stages, biotic and abiotic stresses and cultivars in many plant species such as grapevine [27]; rice [28,29]; tomato [30]; soybean [31]; coffee [32]; brachiaria grass [33]; cotton [34]; eucalyptus [35]; cucumber [36] and petunia [37]. To date, only three studies relying on RT-qPCR analysis in citrus have validated candidate reference genes for transcript normalization. These studies were limited to a few test conditions such as drought [38], leaf tissues of different citrus genotypes and a few organs [39], and Phytophthora parasitica infection [40]. Citrus is one of the most important commercial and nutritional fruit crops in the world. From a scientific standpoint, citrus has proven a valuable resource for studying distinctive aspects of development and physiology such as non-climacteric fruit development, apomixis, gametophytic self-and cross-incompatibility, juvenility, deciduousness versus evergreen foliage, dormancy, seasonality, and root-shoot interaction [41]. In addition, draft genomic sequences of the sweet orange (Citrus sinensis L. Osb.) and clementine mandarin (C. clementina) are now available (http:// www.phytozome.net/clementine.php). The two reference genomes will greatly facilitate studies of functional genomics for genetic improvement in citrus and provide the opportunity to explore peculiar characteristics that cannot be easily addressed in herbaceous model plants such as Arabidopsis thaliana [42,43]. Therefore, the identification of reliable reference genes in citrus will be crucial to allow accurate measurements for gene expression analysis in functional genomics studies. In this study, we aimed to identify potential reference genes suitable for transcript normalization in different samples, tissues, and organs of citrus under different treatments and then validate them. These reference genes will enable more accurate and reliable RT-qPCR normalization for gene expression studies in citrus.

Identification of candidate citrus reference genes
In order to identify suitable citrus reference genes, 15 candidates were chosen from three sources: traditional housekeeping genes frequently used for transcript normalization in citrus; citrus homologues to superior reference genes selected from Arabidopsis transcriptome microarray data [21], and reference genes tested in Swingle citrumelo under drought stress [38]. Gene names, accession numbers, descriptions and functions according to The Arabidopsis Initiative Resource (TAIR) are listed in Table 1. To select citrus coding sequences, a BLASTN search using Arabidopsis reference genes as queries was performed in the CitEST and Harvest Citrus databases. All putative citrus homolog sequences showed very high similarities (see Table 1). Primers amplified a single PCR product as confirmed on a 2% agarose gel ( Figure S1). The stability of expression of the candidate genes was assessed by RT-qPCR in a set of 38 samples grouped into six experiments. The first experimental set was composed of different organs and flower developmental stages from healthy plants, and the remainder were composed of five biotic stresses, including some discrete infection times: two bacterial species (Xylella and Candidatus Liberibacter asiaticus); one fungus (Alternaria alternata); one oomycete (Phytophthora parasitica); and one virus (Citrus leprosis virus C) (Figure 1). In addition, five species of Citrus (C. sinensis L. Osb., C. reticulata Blanco, C. clementina, C. reshni hort. ex Tanaka and C. sunki (Hayata) hort. ex Tanaka) a related specie (Poncirus trifoliata) and a hybrid (Murcott tangor (C. sinensis6C. reticulata)), were included in the set of biological samples evaluated.

RT-qPCR analysis
RT-qPCR was optimized for each primer pair, and two or three independent biological samples under each experimental condition were evaluated in technical triplicates (see Table S1). Melting curve analysis confirmed the presence of a single PCR product from all samples with no primer-dimers ( Figure S2). Amplification efficiency was estimated using the Miner tool; the values ranged from 92 to 98%, except for GAPC2 (84.5%) and PTB1 (79.3%) (Table S1). Cycle quantification for each reaction, determined by the maximum point of the second derivative curve, was also estimated using Miner. Mean Cq values and their standard deviation are presented in Figure 1 for each transcript amplified from each biological replicate. Average Cq values ranged from 20.3 to 32.5; Ef1-a presented the highest and PTB1 transcripts the lowest expression level among all samples (Figure 1a).

Expression stability analysis
In order to find the most stably expressed genes suitable for citrus RT-qPCR normalization, we assessed the stability of expression of 15 candidate genes using the pairwise variation in expression stability implemented in geNorm v3.5 [25]. geNorm estimates two parameters to find the best-suited reference genes: the average expression stability value (M value), and the pairwise variation (V n/n+1 ). The M value is estimated by the pairwise difference between a particular reference gene and all others. At the first step, the M value for all candidate genes is calculated. At the second, the reference gene with the lowest stability of expression (highest M value) is excluded and a new M value is calculated with the remaining reference genes. Moreover, the pairwise variation (V n/n+1 ) will determine the need for inclusion of additional reference genes in the normalization factor to produce accurate and reliable normalization. Quantities (Q) of the 15 candidate reference genes calculated for each biological sample were used in geNorm to calculate M stability values. At each step, reference genes with the lowest stability of transcript accumulation (the highest M) were excluded until the two most stably expressed genes remained. Figure 2a and Table 2 display the M values of reference genes examined when all samples were considered. We found that the FBOX and SAND genes were considered the most stably expressed overall (M = 0.39), while CYP was the least (M = 1.1). In addition, all 15 genes showed acceptable expression stabilities (M#1), as observed by Hellemans and coworkers in heterogeneous samples [43]. We also calculated the optimal number of reference genes needed for a more reliable normalization in geNorm (V n/n+1 ). Taking into account the entire dataset and considering a cut-off (V n/n+1 #0.15, FBOX, SAND and GAPC2 (V 3/4 = 0.13) would be necessary for proper normalization (Figure 3a). In contrast, UPL7 was determined by NormFinder to be the most stable reference gene, whereas CYP was again ranked as the most variable. SAND, FBOX and GAPC2 were ranked in positions 6, 7 and 4, respectively, according to NormFinder (Table 3). Evaluating the six least stable reference genes in both geNorm and NormFinder, we found that TUB, ADP, UBC9, Ctp, DIM1 and CYP were ranked in the same positions. Although the results obtained by the two algorithms seem to be divergent in selecting reference genes suitable for normalizing all citrus sample sets, our results reveal that at least five more stable reference genes (FBOX, SAND, UPL7, PTB1 or GAPC2) could be selected. Besides the analysis in geNorm with the total sample sets, we divided the entire dataset into four subsets that were reanalyzed, and also analyzed by the model-based approach for estimation of expression variation proposed by NormFinder, which estimates the stability of gene expression based on the comparison between inter-and intra-group variability [26].
Analysis of the best reference genes in each experimental subset showed some differences ( Table 2). Eight different tissues or organs composed of vegetative (leaf, branch), reproductive (flower bud at two stages of development, flower and fruitlet), and meristem samples were grouped in subset 1. The DIM1/UBC21 pair was ranked as the most stable reference pair of genes by geNorm and GAPC2 as the most stable gene by NormFinder (Figure 2b and Tables 2, 3). The best combination of two genes according to NormFinder was GAPC2 and UPL7.
Despite the differences, when comparing the M value calculated for the DIM1 and UBC21 genes in geNorm, the exclusion of the GAPC2 or UPL7 genes displayed low variation (0.05 and 0.001, respectively) in the average expression stability value. Both pairs of reference genes (DIM1/UBC21 and GAPC2/UPL7) can be used to normalize the expression of target genes in different tissues or organs of citrus. Analysis of the pairwise variation revealed that the DIM1 and UBC21 genes (V 2/3 = 0.079) would be sufficient for normalizing gene expression (Figure 3b). The CtP, ACT2 and UBC9 genes were considered the most variable reference genes using both algorithms.
The second subset assessed was composed of leaves infected with A. alternata (6 and 12 h post-inoculation), and leaves collected 48 h post-inoculation with P. parasitica. The FBOX/SAND pair was selected as the least variable among all reference genes by geNorm ( Figure 2c and Table 2). In NormFinder, FBOX was the most stable, followed by GAPC2 and SAND (see Table 3). DIM1, which was selected as the most stable in subset 1, was one of the three least stably expressed in subset 2. Furthermore, when the two experimental conditions (fungus and oomycete) were analyzed by geNorm separately, and considering a cutoff of M#0.5, any reference gene except UBC9 (M = 0.55), and at least eight reference genes could be selected as good candidate reference genes for transcript normalization in citrus leaves challenged with P. parasitica or A. alternata, respectively (Table S2 and Figure S3b,e).
In bacterial stress (subset 3), the stability of expression was evaluated under two experimental conditions: in symptomatic leaves of sweet orange infected with Ca. Liberibacter asiaticus versus uninoculated controls, and in leaves 24 h and 7 days postinoculation with X. fastidiosa. FBOX and GAPC2 were calculated to be the most stable genes in geNorm and ACT2 was considered the most stable in NormFinder ( Figure 2d and Table 2). Considering only samples related to Ca. L. asiaticus treatment, the DIM1/ GAPC2 gene pair was considered the most stable in geNorm, followed by FBOX, while TUB was ranked as the worst. In general, all candidate genes except CyP and TUB presented relatively low M values (M#0.5) and could be selected as reference genes for studies of gene expression in citrus infected with Ca. L. asiaticus (Table S2 and Figure S3a). For the treatment with X. fastidiosa, FBOX and SAND were the best reference genes according to geNorm. Again, TUB showed the greatest variation among all the reference genes tested but none had a value of M greater than 0.5. In this case, all genes may therefore be candidates for normalization of gene expression levels in citrus challenged with this pathogen (Table S2 and Figure S3c).
Finally, in the fourth subset evaluated (viral stress), UBC21/ UPL7 was selected as the most stable pair by GeNorm, while DIM1/FBOX was the best combination of two genes in NormFinder. Although TUB and CyP were considered the most variable genes by both programs, none of the candidate genes differed significantly in stability of expression (Figure 2e and  Tables 2, 3). In all treatments individually evaluated, the inclusion of a third gene for more accurate normalization was not required ( Figure S4).
In summary, a comparison of geNorm and NormFinder suggested that FBOX, GAPC2, SAND and UPL7 were the most stable reference genes for all samples and subsets tested in this study. We suggest that these genes could be used as reference genes for accurate transcript normalization in citrus.

Validation of the selected reference genes
In order to validate the selected reference genes, the relative expression level of the gene encoding transcription factor WRKY70 was evaluated in plants infected with Ca. L. asiaticus. In Arabidopsis, this gene acts as an activator of salicylic aciddependent defense genes and a repressor of jasmonic acidregulated genes. We also found that citrus WRKY70 is an important gene in response to infection with Ca. L. asiaticus and americanus. According to our microarray analysis, WRKY70 was upregulated in symptomatic sweet orange plants in relation to uninoculated control plants (Mafra et al., unpublished data). Primer design, RT-qPCR and amplification efficiency calculus were performed as described above, and primer sequences are listed in Table S1. WRKY70 was normalized to the three most stable candidate reference genes (DIM1, GAPC2 and PTB11) and the two least stable (CYP and TUB) as determined by geNorm analysis. RT-qPCR analysis showed that the expression level of WRKY70 transcript significantly increased during symptoms in relation to uninoculated controls (fold change, FC = 3.19) (Figure 4). Increased expression of this transcript corroborates our microarray expression data in plants infected with Ca. L. americanus (FC = 5.13). A similar expression pattern was described by Albrecht & Bowman (2008), who found that WRKY70 transcript accumulated in sweet orange leaves infected with Ca. L. asiaticus both at 5-9 and 13-17 weeks after inoculation (FC = 1.9 and 2.3, respectively) [19]. In order to demonstrate the effect of using different normalization genes to estimate relative accumulation, we used the three most stable genes (NF 3 ) and the two least stable (NF 2 ) to calculate normalization factors. As shown in Figure 4, the use of TUB only as a reference gene or grouped with CyP to normalize the transcript level of WRKY70 led to an increase in the fold change (FC = 42 and FC = 106, respectively) compared to the values obtained with the two or three most stable reference genes. These results were expected because of the high variability of M calculated by geNorm for the Ca. L. asiaticus condition, when TUB was included among the 15 genes assessed (M = 0.876) (Table S2).

Discussion
RT-qPCR has become an important tool to understand gene expression in several biological systems. For accurate RT-qPCR measurements, endogenous reference genes are used as internal controls. An appropriate reference gene should be expressed with minimal change regardless of the experimental conditions. Because there is no reference gene that is universally stable in expression, it is necessary to identify candidate genes specifically chosen for transcript normalization for the conditions under study [2,44].
Here, we evaluated the stability of expression of eleven novel and four traditional reference genes in citrus from different tissues and under different biotic stresses. Our analysis in geNorm and NormFinder showed some differences, especially in the top ranked genes, but both programs very consistently excluded the same genes as showing unstable expression patterns. This apparent divergence probably reflects differences in the statistical algorithms. The NormFinder program employs a model-based variance estimation approach to identify genes suitable for normalization. In practice, it estimates both the intra-and inter-group variation and combines them into a stability value. This model-based approach ranks the top genes with minimal estimated inter-and intra-group variation. In contrast, the pairwise approach performed by geNorm selects two genes with the highest degree of similarity in expression profile and the lowest intra-group variation. For this reason, it is not surprising that the two algorithms differ in the ranking of the best candidate genes. This divergence in results obtained by the two methods was highlighted in the original paper describing the NormFinder strategy [26]. Discrepancies between NormFinder and geNorm were also demonstrated by other studies [4,5,32,52]. This approach could be problematic if co-regulated genes exhibit similar expression profiles and thus, might be preferentially top ranked [24]. Under viral stress, we found that the UBC21/UPL7 pair was identified as the most stable by geNorm, followed by UBC9. The UBC21 and UBC9 genes encode ubiquitin-conjugating enzymes belonging to the E2 class, whereas UPL7 encodes a ubiquitin-protein ligase grouped into the E3 class. Ubiquitin conjugation is a protein modification that occurs in a multistep reaction, sequentially involving an E1 enzyme (ubiquitinactivating enzyme), an E2 enzyme and an E3 enzyme [45]. In A. thaliana, it was estimated that there are two E1 proteins, 37 E2 proteins and more than 1,300 predicted E3 proteins [46]. Although E2 and E3 proteins participate in the same pathway, there is no evidence that UPL7 and UBC21 interact directly and may be co-regulated. Indeed, only in viral stress and when comparing different organs were these two genes ranked among the top three, while in other subsets they presented intermediate stability values. Moreover, in Arabidopsis, these genes were not top ranked by geNorm software, but occupied close positions in the ranking [21].
Our results demonstrated that FBOX, GAPC2, SAND and UPL7 were the most stably expressed reference genes in all samples and subsets studied. Nevertheless, the best combination of genes varied significantly depending on experimental condition. This observation reinforces the necessity to assay the stability of expression of candidate genes to select suitable reference genes for reliable normalization in a specific biological assay. Among the top reference genes, FBOX was identified as the most stable, followed by GAPC2 and SAND. Our results corroborate a recently published paper by Lilly et al. [47], which tested reference genes for normalization of transcripts from virus-infected A. thaliana. They found that FBOX and SAND showed the most stable transcript accumulation. Similar results in Arabidopsis were observed by Remans et al. [23], in which the same two genes, along with YLS8, were identified as the best candidates for data normalization (M,0.3) in roots and shoots in treatments with cadmium and copper. FBOX and SAND were also ranked among the top 22 most stable reference genes tested in 79 samples including different developmental stages, organs, tissues and genotypes [21]. In soybean, a possible ortholog of FBOX was identified as the most uniformly expressed gene [31]. Additionally, FBOX was considered a good reference gene for normalization of floral organs in cotton; however, when all organs were compared, this gene was ranked among the three least stable [34]. Despite slight differences found in different studies, we concluded that FBOX is a good candidate gene for normalizing a wide range of tissue and organ samples and different conditions in plants, even though the molecular function and biological process this gene is associated with remain unclear. SAND was ranked in our study as the third most stable reference gene. Similarly, SAND was revealed as one of the superior reference genes found for proper normalization in tomato development studies and a set of organs and tissues of buckwheat [30,48]. Also, SAND and RAN1 were calculated as the most stable pair when the entire dataset was evaluated in petunia, while GAPC was the most variable gene (M = 1.15). SAND was first described in the Saccharomyces cerevisiae genome. Later, with the availability of several eukaryotic genomes, homologous sequences were identified in Caenorhabditis elegans, Drosophila melanogaster and A. thaliana. In plants, only one SAND sequence was found in monocots and dicots [49]. Functional studies suggest that the SAND family proteins are involved in late steps of endocytic transport [49,50]. As suggested by Lilly et al. [47], SAND may not be a suitable reference gene in studies of gene expression in response to pathogens that could interfere with vesicle traffic, like viruses.
Genes commonly referred to as housekeeping genes, such as tubulins, actins, GAPDH, ribosomal subunits and elongation factors, have been used in several studies in citrus to normalize gene expression data. However, there is a consensus that the use of such genes arbitrarily may result in the misinterpretation of results [6,51]. In our study, GAPC2 was selected as the second most stable gene overall. In the entire dataset, GAPC2 was ranked as the third most stable (M = 0.5), and in Ca. L. asiaticus infection, GAPC2 and DIM1 were considered the best combination by geNorm (M,0.1). These results are consistent with citrus leaves subjected to drought stress, for which the EF1/ADP pair was considered by geNorm and NormFinder to be the best combination of genes, followed by  GAPC2 [38]. Similar results were observed in Brachypodium, for which NormFinder considered GAPC2 among the three most stable genes when comparing different tissues (M = 0.28), treatment with growth hormone (M = 0.18), and exposure to heat or cold and high salt or drought stress (M = 0.07) [4]. However, our results are in contrast with those of Boava et al. [35], who ranked GAPDH among the three genes least stable in all conditions assessed. In petunia, GAPDH was again considered the gene least stably expressed when assessed during leaf and flower development [37]. Unlike GAPC2, which was shown to be a good reference gene for citrus normalization in different biological contexts, TUB was ranked in the last position for different subsets analyzed in our study. Our results corroborate those obtained by Carvalho et al. [38], which considered TUB inadequate for transcript normalization in citrus under drought stress.
Actin, another reference gene frequently used in citrus, displayed an intermediate stability pattern in our analysis. Expression instability was also described by Czechowski et al. [21], who found ACT2 to be the least stably expressed gene among the 27 tested. Stability of ACT and TUB was also assessed in flax and both were considered unreliable for transcript normalization during flax development [52]. Under drought stress in citrus, ACT2 also showed unstable transcript abundance. Given these observations, we suggest that both TUB and ACT2 should be carefully evaluated before using them as reference genes for citrus gene expression studies.
Finally, to validate the suitability of the reference genes we identified in this study, we assessed the expression profile of a WRKY70 homolog in leaves of sweet orange plants infected with Ca. L. asiaticus. We demonstrated that the use of the two most variable reference genes (CYP/TUB) or TUB resulted in an increase of the relative transcript abundance of WRKY70 compared to the normalized expression data obtained using the two or three most stable ones (DIM1/GAPC2/PTB1 or DIM1/ GAPC2). These results indicate that the incorrect use of reference genes without validation may introduce bias in the analysis and lead to misinterpretation of data. Matta and collaborators [24] reported similar results in qPCR studies of Drosophila, emphasizing the need for validation of the best set of reference genes for each experimental condition tested.
In summary, we evaluated several suitable reference genes in different citrus organs and following different biotic stresses. We also identified novel reference genes that outperformed housekeeping genes commonly used in citrus and showed that some of these housekeeping genes could be inadequate for transcript normalization under particular experimental conditions. We propose FBOX, SAND, GAPC2 and UPL7 as good candidate genes to be tested as reference genes for normalization in citrus gene expression studies. In addition, we provide a list of twelve genes with the potential to be good reference genes. This work constitutes the first systematic study in citrus to identify and validate optimal reference genes for RT-qPCR normalization with consideration of different tissues, genotypes and biotic stress conditions.

Plant materials and experimental conditions
Biotic stress assays. The following citrus species and hybrids were included in the evaluation: sweet orange (C. sinensis L. Osbeck), Ponkan mandarin (C. reticulata Blanco), clementine (C. clementina hort. ex Tanaka), Sunki mandarin (C. sunki (Hayata) hort. ex Tanaka), Cleopatra mandarin (C. reshni hort. ex Tanaka), Murcot tangor (C. sinensis L. Osb.6C. reticulata Blanco), and P. trifoliata. (L.) Raf. All experiments testing a biotic stress were conducted in a greenhouse or growth chamber and are summarized in Table 4. The samples infected with systemic (X. fastidiosa, Ca. L. asiaticus, and A. alternata), or nonsystemic pathogens (CiLV-C, P. parasitica) were collected and immediately frozen in liquid nitrogen. For a detailed description of each biotic stress assay, see File S1 in supporting information.
Plant tissues, organs and developmental stages used for sampling. Three 15-year-old 'Valencia' orange (C. sinensis L. Osbeck) plants grafted onto Cleopatra mandarin (C. reticulata Blanco) were used. These trees are cultivated in an experimental field of the Centro de Citricultura Sylvio Moreira, located in Cordeirópolis, São Paulo state, Brazil. Samples of adult leaves, branches, fruitlets (8 mm length) open flowers, and flower buds (5 mm and 10 mm length) were collected during bloom. We also collected meristem samples during winter and early spring. Samples were transferred to liquid nitrogen and stored at 280uC until required.

Total RNA isolation and cDNA synthesis
About 200 mg of tissue was ground to a fine powder in liquid nitrogen using a mortar and pestle. Total RNA was extracted using an RNeasy Plant Mini Kit (Qiagen) according to the manufacturer's instructions with minor modifications. Genomic DNA contamination was removed by digestion in the RNeasy columns with recombinant DNAse I (Qiagen). Total RNA concentration and purity were determined from the ratio of absorbance readings at 260 and 280 nm using a Nanodrop ND8000 spectrophotometer (Nanodrop Technologies), and RNA integrity was tested in a denaturing agarose gel. Reverse transcription was performed with 1 mg of total RNA in a total volume of 20 mL with oligo(dT) primer using Revertaid H-Minus reverse transcriptase (Fermentas). The final cDNA products were diluted 50-fold prior to use in RT-qPCR.

Selection of potential reference genes in citrus and primer design
The 15 candidate genes evaluated in this experiment were selected from the CitEST (http://limonia.centrodecitricultura. br/blast/blast.html) and HarvEST (http://www.harvest-web. org/) citrus databases according to meeting one or more of the following criteria: (1) reference genes traditionally used in citrus for transcript normalization; (2) reference genes described in the literature for RT-qPCR normalization in Swingle citrumelo (C. paradise6P. trifoliata) under drought stress [38]; and (3) citrus homologues of reference genes tested for transcript level normalization and quantification in Arabidopsis [21]. BLASTN with a default setting was used to search for citrus coding sequences with high similarity (E-value#1e-90) to Arabidopsis genes. Primers were designed with Primer 3 (http://frodo.wi. mit.edu/primer3/) and Oligo Explorer 1.1.2 software tools (http://www.uku.fi/,kuulasma/OligoSoftware/) with the following parameters: Tm around 60uC and amplicon length of 90 to 120 bp, yielding primer sequences with a length of 19 to 23 nucleotides with an optimum at 20 nucleotides, and a GC content of 45 to 60%. Primers were also designed as much as possible to allow the amplification of transcript isoforms from all citrus genotypes. The specificity of the resulting primer pair sequences was checked against the Arabidopsis transcript database using TAIR WU-BLAST2 (www.arabidopsis.org/wublast/in-dex2.jsp). Amplicon specificity was checked by 2% (w/v) agarose gel electrophoresis and by melting-curve analysis. The sequence of the 15 amplicons was confirmed by sequencing (data not shown). PCR products were cloned into pGEM-T Easy vector and sequenced using an Applied Biosystems Model 3730 capillary DNA sequencer. RT-qPCR conditions and statistical analysis RT-qPCR was performed in a 96-well optical plate with an ABI PRISM 7500 FAST sequence detection system (Applied Biosystems). The reaction mixture contained 9 mL 2x FAST SYBR Green Master Mix reagent (Applied Biosystems), 3 mL diluted cDNA (1:50), 120 or 150 nM of each gene-specific primer pair in a final volume of 25 mL. The following standard thermal profile was used for all amplifications: 95uC for 20 sec followed by 40 cycles of 95uC for 3 sec, and 60uC for 30 sec. All assays were performed using three technical replicates and a non-template control, as well as two or three biological replicates. To analyze dissociation curve profiles, the following program was run after the 40 cycles of PCR: 95uC for 15 sec followed by a constant increase in temperature between 60 and 95uC. Primer efficiency for each experimental set was estimated using an algorithm in Real-time PCR Miner software (http://www. miner.ewindup.info/) that calculates primer efficiency and quantification cycle (Cq) values based on the kinetics of individual reactions without the need for a standard curve. Cq values, determined by the second derivative maximum for each biological sample, were converted into non-normalized relative quantities using the formula Q = E DCq , where E represents the arithmetic mean of efficiency of all samples for each gene, and DCq represents the difference between the arithmetic mean Cq value across all samples for this gene, and the Cq value of the sample in question, as recommended by Hellemans et al. [43]. These quantities were imported into geNorm v3.5 (medgen.ugent.be/,jvdesomp/geNorm/) [25] and NormFinder (www.mdl.dk/publicationsnormfinder.htm) [26] for reference gene selection. First, we performed a global analysis composed of all biological samples in geNorm. Considering the heterogeneity of treatments, we then analyzed each experimental condition individually in an attempt to identify specific reference genes according to the treatments. Finally, once NormFinder calculated both inter-and intra-group variation in the expression stability, thus identifying the best combination of reference genes, we established four subsets composed of the following treatments: tissue or organ (n = 24); viral stress (n = 18), fungal or oomycete stress (n = 48), and bacterial stress (n = 32). These subsets were then analyzed by both geNorm and NormFinder. Figure S1 RT-qPCR amplification specificity of the15 reference genes. Amplification fragments were separated by 2% agarose gel electrophoresis. UNK: unknown protein.    Figure S3 Reference genes ranked according to their expression stability as determined by geNorm for each experimental condition. A lower M value indicates more stable expression. The ranking of the reference genes is in Table  S1.  Figure S4 Pairwise variation (V) to determine the optimal number of reference genes for each experimental condition. The ranking of the reference genes is in Table S1.