key: cord-0911981-whm2ce4u
authors: Li, Qingdi Quentin; Skinner, Jeff; Bennett, John E
title: Evaluation of reference genes for real-time quantitative PCR studies in Candida glabrata following azole treatment
date: 2012-06-29
journal: BMC Mol Biol
DOI: 10.1186/1471-2199-13-22
sha: 69d2e69e949dab8909be8c7b764a0d5469546b31
doc_id: 911981
cord_uid: whm2ce4u

BACKGROUND: The selection of stable and suitable reference genes for real-time quantitative PCR (RT-qPCR) is a crucial prerequisite for reliable gene expression analysis under different experimental conditions. The present study aimed to identify reference genes as internal controls for gene expression studies by RT-qPCR in azole-stimulated Candida glabrata. RESULTS: The expression stability of 16 reference genes under fluconazole stress was evaluated using fold change and standard deviation computations with the hkgFinder tool. Our data revealed that the mRNA expression levels of three ribosomal RNAs (RDN5.8, RDN18, and RDN25) remained stable in response to fluconazole, while PGK1, UBC7, and UBC13 mRNAs showed only approximately 2.9-, 3.0-, and 2.5-fold induction by azole, respectively. By contrast, mRNA levels of the other 10 reference genes (ACT1, EF1α, GAPDH, PPIA, RPL2A, RPL10, RPL13A, SDHA, TUB1, and UBC4) were dramatically increased in C. glabrata following antifungal treatment, exhibiting changes ranging from 4.5- to 32.7-fold. We also assessed the expression stability of these reference genes using the 2(-ΔΔCT) method and three other software packages. The stability rankings of the reference genes by geNorm and the 2(-ΔΔCT) method were identical to those by hkgFinder, whereas the stability rankings by BestKeeper and NormFinder were notably different. We then validated the suitability of six candidate reference genes (ACT1, PGK1, RDN5.8, RDN18, UBC7, and UBC13) as internal controls for ten target genes in this system using the comparative C(T) method. Our validation experiments passed for all six reference genes analyzed except RDN18, where the amplification efficiency of RDN18 was different from that of the ten target genes. Finally, we demonstrated that the relative quantification of target gene expression varied according to the endogenous control used, highlighting the importance of the choice of internal controls in such experiments. CONCLUSIONS: We recommend the use of RDN5.8, UBC13, and PGK1 alone or the combination of RDN5.8 plus UBC13 or PGK1 as reference genes for RT-qPCR analysis of gene expression in C. glabrata following azole treatment. In contrast, we show that ACT1 and other commonly used reference genes (GAPDH, PPIA, RPL13A, TUB1, etc.) were not validated as good internal controls in the current model.

The investigation of gene expression has become increasingly prevalent in numerous animal, human, microorganism, and plant studies [1] [2] [3] [4] [5] . The quantitation of gene expression requires sensitive, precise, and reproducible measurements for specific mRNA sequences. Generally, gene expression levels can be determined by a variety of techniques, including Northern blotting, RNase protection assay, semi-quantitative reverse-transcription PCR, and real-time quantitative PCR (RT-qPCR) [4] . RT-qPCR has gained favor as it is a highly sensitive, accurate, and fast technique that offers highthroughput and the ability to detect low-abundance mRNAs [6] and quantify mRNA copy number [7] . Thus, RT-qPCR has been used for countless different applications [1] [2] [3] [4] [5] .

One of the main uses of RT-qPCR, when coupled with reverse transcription, is to measure gene expression at the mRNA level in various biological samples. However, there is substantial technical variability associated with RT-qPCR, arising from inherent differences in samples, sample collection, RNA degradation and extraction efficiency, quantity and quality of input RNA, reverse transcription and PCR efficiency, and pipetting accuracy or error. Researchers have employed a number of strategies to normalize their data, including normalization to (i) genomic DNA, (ii) total RNA, (iii) an external standard, and (iv) a reference gene. The most common practice is to normalize to an internal control gene termed a reference gene. A reference gene is subject to the same errors in cDNA preparation as the gene of interest, making it an excellent normalizing control. However, selection of an inappropriate reference gene can add large unpredictable error to the analysis and result in incorrect estimates [8] . The ideal reference gene should have a stable RNA transcription level under different experimental conditions and be sufficiently abundant across different tissues and cell types. However, it has become apparent that such an ideal reference gene has not yet been identified [9] . The most commonly used reference genes, including β-actin, cyclophilin, GAPDH, tubulin, and 18S and 28S ribosomal RNAs, have shown variable expression levels in different cells and tissues under different conditions, and therefore they are unsuitable for normalization purposes owing to large measurement error [6, . Hence, it is no longer acceptable to arbitrarily select any reference gene for normalization; it must be demonstrated that the reference gene of choice is suitable for the experiment in question.

In recent decades, Candida glabrata has emerged as the second most common cause of invasive fungal infection [1, 34, 35] . Azoles such as fluconazole are the firstline drugs for the treatment of fungal infections caused by C. glabrata. However, resistance to azoles can arise rapidly in C. glabrata during treatment of patients with azoles [36] . An increasing body of evidence has implicated ATP-binding cassette transporters (e.g., Cdr1 and Pdr1) and sterol biosynthetic enzymes (e.g., Erg3 and Erg11) in azole resistance in C. glabrata in both clinical and laboratory settings [1, [34] [35] [36] [37] . The expression of these genes in C. glabrata in response to azoles is not completely understood. Therefore, we set out to establish an in vitro model for investigating azole-inducible gene expression in C. glabrata, using RT-qPCR. For reliable gene expression analysis, a compulsory step is the selection of good reference genes for normalization; however, no validated reference genes have been reported for the relative quantification of the mRNA expression profile in C. glabrata following exposure to azoles.

We have been using ACT1 as the internal control for gene expression analysis by RT-qPCR in clinical isolates of C. glabrata in the absence of drug challenges [1] . Other researchers also use ACT1 as the reference gene for azole-inducible gene expression studies by slot blotting in Candida species [34, 37, 38] . However, the suitability of ACT1 in studies of azole-inducible gene expression in C. glabrata has not been validated. In this work, we evaluated 16 reference genes to establish their suitability as control genes for normalization and identified a set of genes that are suitable for quantitative gene expression analysis by RT-qPCR in C. glabrata following fluconazole treatment.

All five C. glabrata strains (Table 1 ) used in the present study were grown in YPD broth (Difco Laboratories, The susceptibility of each C. glabrata strain to fluconazole was determined on YPD agar medium using an E-test (AB Biodisk, Solna, Sweden) according to the manufacturer's instructions (Table 1) .

Total RNA was extracted from C. glabrata logarithmicphase cultures grown in YPD broth, using TRIzol reagent (Invitrogen, Life Technologies, Grand Island, NY, USA) according to the manufacturer's instructions. The concentration and purity of the RNA was determined using a UV spectrophotometer (NanoDrop 2000C; Ther-moFisher Scientific, Waltham, MA, USA) by measuring the absorbance at 230 (OD 230 ), 260 (OD 260 ) and 280 nm (OD 280 ). The OD 260nm /OD 280nm of the samples, reflecting the average purity, ranged from 1.80 to 2.05, and the OD 260nm /OD 230nm was in the range of 2.00-2.60. The integrity of the RNA was further checked in a selected subset of samples by electrophoresis through 1% denaturing and non-denaturing agarose gels.

Reverse transcription (RT) was performed on 1 μg of total RNA using a commercially available kit. Prior to RT, the total RNA samples were treated with DNase for 30 min at 37°C (TURBO DNA-free; Ambion, Life Technologies, Grand Island, NY, USA) according to the manufacturer's instructions. RNA was converted to cDNA using the High Capacity cDNA Reverse Transcription Kit (Applied Biosystems, Life Technologies, Carlsbad, CA, USA). The reaction took place in a thermal cycler (T3 Thermocycler; Biometra, Goettingen, Germany) with a single cycle and incubation periods of 25°C for 10 min, 37°C for 120 min, and 85°C for 5 min. All investigated samples were transcribed with the same reverse transcription reaction conditions. Negative controls, which were run simultaneously, did not contain either RNA (no template control) or no reverse transcriptase (RT negative control), to control for RNA and genomic DNA contamination, respectively.

Primers and probes were designed in our laboratory using the primer analysis software Primer Express 3.0 (Applied Biosystems). TaqMan probes were synthesized by Applied Biosystems, and primers were synthesized by Invitrogen/Life Technologies. The primers and TaqMan probes used in the current study were selected to bind specifically to the cDNAs of Cg84u and other C. glabrata strains ( Table 1 ). The sequences of TaqMan probes and forward and reverse primers, the gene numbers, and the localization for each PCR assay for the 16 reference genes and 10 target genes assessed in this study are listed in Additional file 1, Additional file 2, and Additional file 3.

RT-qPCR for reference gene RNA transcription was performed by SYBR Green chemistry (SYBR Green PCR Master Mix; Applied Biosystems). The increase in fluorescence of the SYBR Green dye was monitored using a 7500 Real-time PCR System (Applied Biosystems). This technique has been successfully used to validate reference gene expression levels in yeast and other cell types [41] . Primers were used at 300 nM each for specific forward and reverse primers and cDNA at 25 ng in 25-μl reactions. Primer sets for the reference genes (Additional file 2) were used to amplify the open reading frame (ORF) region of the genes according to the following conditions: one cycle of 50°C × 2 min, 95°C × 10 min; followed by 40 cycles of 95°C × 15 s, 60°C × 1 min; with dissociation (a melting curve) during the last cycle of 95°C × 15 s, 60°C × 1 min, 95°C × 15 s. The dissociation protocol to determine the melting curve from 60°C to 95°C for each PCR product was added after thermocycling to verify that each primer pair produced only a single product. All samples gave only a single peak, indicating a single pure product and no primer/dimer formation. Real-time PCR efficiencies were acquired by amplification of a standardized dilution series of the template cDNA and were determined for each gene as the slope of a linear regression model. PCR efficiency was determined by measuring the C T to a specific threshold for a serial dilution of cDNA. The corresponding real-time PCR efficiencies were then calculated according to the equation: E ¼ 10 À1=slope À 1 À Á Â 100 . All PCRs displayed efficiencies between 94% and 119%.

To study target gene expression, the amplification was detected in real time using TaqMan chemistry (TaqMan Universal PCR Master Mix; Applied Biosystems) according to the manufacturer's instructions. RT-qPCR was performed in 96-well microtiter plates with a final volume of 25 μl, using a 7500 Real-time PCR System (Applied Biosystems). Primers were used at 300 nM each for specific forward and reverse primers; probes, at 200 nM; and cDNA, at 25 ng in 25-μl reactions. Primer sets and TaqMan probes for the target genes (Additional file 3) were used to amplify the ORF region of the genes under the following conditions: one cycle of 50°C × 2 min, 95°C × 10 min; and then 40 cycles of 95°C × 15 s, 60°L

C × 1 min. The parallel amplification between the reference genes and the target genes was confirmed for each with probe-primer sets. To minimize technical (run-torun) variation between the samples, all samples were analyzed in the same run for both target genes and reference genes.

Evaluation of reference gene expression stability using four different software packages

Non-normalized gene expression levels from our experimental data were analyzed to evaluate the expression stability of potential reference genes, using four different software programs: hkgFinder [42] , geNorm [9, 43] , Best-Keeper [44] , and NormFinder [26] . The hkgFinder software computes the pooled standard deviation (SD) of non-normalized expression data from both phenotypes (i.e., azole-treated and untreated C. glabrata cells), the fold change (FC) values between the two phenotypes, and Student's t-tests of the log 2 fold-change values with Holm-adjusted P-values. The reference genes with the smallest SD and the smallest, non-significant FC are identified as the best potential reference genes (Additional file 4). The geNorm software computes a stability value (M) and a pairwise variation (V), which are used to evaluate each individual reference gene candidate or each combination of reference genes, called a normalization factor (NF). The pairwise variability V of two genes j and k is the standard deviation of all log 2 ratios of a j /a k , while the stability value M of gene j is the mean of all possible pairwise variations V jk . Graphs of the M values help identify the best individual reference genes, and graphs of the V values identify the optimal number of reference genes for an NF. Note that an earlier version called geNorm Excel was produced as an add-in for MS Excel and it required several hand calculations to convert crossing point (CP) values into relative expression values. That version is now unavailable and the new geNorm PLUS from Biogazelle does not require those hand calculations. The BestKeeper software uses pairwise correlation to determine whether potential reference genes should be included in a BestKeeper Index, which is simply the geometric mean of the CP or cycle threshold (C T ) values. The NormFinder software computes a different type of stability value (ρ ig ) based on the intragroup and intergroup variation of the expression data. The software instructions from each package were followed when inputting the RT-qPCR data, fetching the output, and interpreting the analysis results.

Stability of RNA transcription of reference genes in C. glabrata following azole stimulation

In the present study, 16 reference genes were chosen from among commonly used reference genes in published studies with yeast and mammalian cells, paying close attention to selecting genes that belong to different functional classes; their full names, symbols, functions, and gene numbers are listed in Additional file 1. Our aim was to identify reference genes with minimal variability under our experimental conditions. To this end, RT-qPCR was used to measure the RNA transcription levels of 16 reference genes in C. glabrata cells following fluconazole treatment. To compare the different RNA transcription levels after azole exposure, the C T values of the reference genes were directly compared between the drug-treated (t) and untreated (ut) samples using the formula:

The C T is defined as the number of cycles needed for the fluorescence signal to reach a specific threshold level of detection and is inversely correlated with the amount of template cDNA present in the reaction. Thus, a higher value of C T Change indicates lower stability of a reference gene, considering that the expression of a reference gene should not change significantly with azole treatment. As expected, the RNA transcription levels of the reference genes varied ( Table 2 ). The three ribosomal RNA subunits RDN5.8, RDN18, and RDN25 were the most stable reference genes, with C T Change values less than 0.5, while UBC13, PGK1, and UBC7 were relatively stable with C T Change values of only around 1.5. By contrast, the other 10 reference genes showed marked variation in response to fluconazole. Among them, the most prominent variation was found in the RNA transcription levels of SDHA, ACT1, and RPL13A; as seen in Table 2 , the C T Change values of these reference genes were as high as 5.03. To validate the stability of candidate RNA transcription under our experimental conditions, the levels were compared with the RDN5.8 RNA transcription level. We chose to use RDN5.8 as a normalizer because it meets the requirement for both stability and suitability as a reference. First, we calculated the ΔC T between the C T values of reference genes and RDN5.8 from fluconazoletreated (t) and untreated (ut) cells:

In the second step, we subtracted the change in RNA transcription in untreated samples from the change in treated samples to obtain the ΔΔC T (t):

ΔΔC

Thus, ΔΔC T (t) indicates the change in RNA transcription caused by fluconazole treatment after normalization to RNA transcription changes in RDN5.8. A high ΔΔC T (t) value indicates a significant fluconazole-related change in the RNA transcription level of the tested gene. A positive ΔΔC T (t) value indicates down-regulation of transcription, whereas a negative ΔΔC T (t) indicates upregulation of a gene's transcription following azole treatment. We then transformed ΔΔC T (t) into a 2 -ΔΔCT value, which indicates the fold change in RNA transcription of a reference gene in response to fluconazole as compared with the level in untreated cells. The calculated ΔΔC T (t) and 2 -ΔΔCT values of the 16 tested reference genes in drug-treated samples are given in Table 2 .

Following stimulation with fluconazole, the RNA transcription of SDHA, ACT1, and RPL13A was highly regulated in C. glabrata cells, with changes ranging from 13to 23-fold compared with transcription in the untreated cells. There was almost no regulation of RDN5.8, RDN18, and RDN25 RNA transcription, while PGK1, UBC7, and UBC13 RNA transcription were only approximately 2-fold induction in response to drug treatment ( Table 2) .

To choose the best reference genes, the reference gene stability was evaluated using four different software packages: hkgFinder, geNorm, BestKeeper, and Norm-Finder. Each of these software packages uses a slightly different metric to evaluate the candidate reference genes. Our goal was to compare the findings from these four different methods and look for the best-scoring reference genes that might be common to these different methods.

The hkgFinder software identifies the best reference genes by ranking the candidate genes according to their SD and FC values ( Table 3 ). Among the 16 potential reference genes, the SDs ranged from 0.19 to 2.76, and the FCs ranged from 1.2 to 32.7. The best three reference gene candidates were RDN18, RDN25, and RDN5.8. The next three best candidate reference genes, which also had reasonable SD and FC values, were UBC13, PGK1 and UBC7.

The geNorm software evaluates reference genes by their M-stability values and V-pairwise variability values.

Low M values represent more stable expression and thus the most suitable reference genes ( Figure 1 ). The geNorm analysis identified RDN18, RDN25, and RDN5.8 as the three most stable genes; UBC13, PGK1, and UBC7 as relatively stable genes; and RPL13A, ACT1, and SDHA as the three least stable genes under fluconazole treatment in C. glabrata. Interestingly, the ranking of expression stability of the 16 reference genes was identical between the geNorm program and the hkgFinder tool ( Figure 1 ; Tables 3 and 4 ). The geNorm program also estimates the optimal number of reference genes that could be used in combination as an NF value ( Figure 2 ). Each NF was calculated as the geometric mean of the two most stable genes, then the pairwise variability V was computed between NFn and NFn + 1 for n = 2, . . ., 15. Vandesompele et al. [9] proposed 0.15 as a cutoff value for V below which additional reference genes do not need to be added to the NF. Adding the third gene to the most stable two reference genes, RDN18 and The fold change shown in the table represents the difference in reference gene expression between azole-treated and untreated C. glabrata without normalization to an internal control gene. The stability ranking is based on the values of standard deviation (SD), log fold change, and fold change of reference genes in fluconazole-treated C. glabrata cells. The best reference genes will have the smallest SD and smallest fold-change values. 

Determination of the optimal number of reference genes as internal references for normalization using geNorm analysis. The geNorm program calculates a normalization factor from at least two reference genes and the mean pairwise variation (V) between every combination of sequential normalization factors in order to determine the minimum number of reference genes required for accurate normalization in the samples. For example, V5/ 6 represents the comparison of the normalization factors from five and six reference genes, respectively. On the left-most side is the pairwise variation when the number of reference genes is increased from two to three (V2/3). Stepwise inclusion of less stable genes generates the subsequent data points. A decrease in the V value indicates a positive effect and means that the added gene should preferably be included for calculation of a reliable normalization factor. The cutoff value for V, below which the inclusion of an additional reference gene does not result in a significant improvement of normalization, was set at 0.15. It was apparent from the analysis of all studied samples that the combination of the two most stable reference genes is the best option and the combination of the five most stable reference genes is the second-best option for accurate normalization.

The NormFinder program was also used to rate candidate reference gene stability according to a stability value computed from the intragroup and intergroup expression variability. The least reliable reference genes identified by this program were RPL13A, ACT1, and SDHA, which were identical to the worst reference genes identified by geNorm and hkgFinder analyses (Table 4) . However, the ranking order of the most stable genes and the relatively stable genes by the NormFinder program was different from that generated by geNorm and hkgFinder (Figure 1; Tables 3 and 4) . The geNorm and hkgFinder analyses graded RDN18, RDN25, and RDN5.8 as the most stable reference genes, followed by UBC13, PGK1, and UBC7 based on gene expression stability, whereas NormFinder rated UBC7, PGK1, and UBC13 as the most stable reference genes, followed by RDN5.8, RDN18, and RDN25 (Figure 1 ; Tables 3 and 4) .

Finally, the BestKeeper program was used to grade candidate reference gene stability. This approach permits a comparative analysis across reference genes. Ten reference genes analyzed were correlated and were combined into an index. Subsequently, the correlation between each reference gene and the index was calculated. The best correlations between the reference genes and the BestKeeper index were obtained for UBC7, UBC13, and PGK1 (r = 0.983, 0.978, and 0.974, respectively; Table 4 ). The rankings of the top three and the last three reference genes identified by the BestKeeper program were the same as those generated by the NormFinder analysis, although the order of stability of the other reference genes differed slightly between the two programs ( Table 4 ).

Following the identification of the most stable reference genes from the full gene panel of 16 genes, the comparative C T method was used to validate their suitability.

The comparative C T method, also referred to as the ΔΔC T method, is a relative quantitation of gene expression between a specific target gene and a reference gene. For the comparative C T method to be valid, the efficiency of the target amplification and the efficiency of the reference (internal control) amplification must be approximately equal, and this must be determined in a validation experiment.

To this end, we first determined the amplification efficiency of 10 target genes (CDR1, PDH1, PDR1, SNQ2, YOR1, ERG2, ERG3, ERG4, ERG10, and ERG11) and six reference genes (ACT1, PGK1, RDN5.8, RDN18, UBC7, and UBC13). Standard curves were generated by plotting the dilutions of the cDNA of each gene against the C T values. The linear correlation coefficient (r 2 ) for all 10 target genes and the six reference genes ranged from 0.98 to 1.0. Based on these slopes of the standard curves, the amplification efficiencies of the cDNA standards, derived from the formula E = (10 −1/slope −1) × 100, ranged from 94 to 119%. The C T values of all 16 genes in the samples were within the range of the standard curves. Next, the ΔC T (ΔC T = C T target − C T reference ) was calculated using the C T values generated from standard curve mass points (target vs. reference gene). These ΔC T values were then plotted versus log 10 input amount of cDNA to create a semi-log regression line. The slope of the resulting semi-log regression line was used as a general criterion for passing a validation experiment. In a validation experiment that passes, the absolute value of the slope of ΔC T versus log 10 input cDNA would be <0.1, meaning the two C T versus log 10 concentration curves are nearly parallel. As seen in Table 5 , our validation experiments passed for all reference genes analyzed except RDN18, which had an absolute value >0.1 for the slopes of ΔC T versus log 10 input cDNA for all 10 target genes evaluated. Thus, the amplification efficiency of RDN18 was clearly different from that of the ten target genes, whereas the other five reference genes (ACT1, Each value is the slope m of the line (Y = mx + b) of the validation experiment and reflects the correlation of reference gene and target gene amplification efficiencies. An absolute slope value <0.1 is generally used as a criterion for passing a validation experiment, as it indicates that the amplification efficiency is approximately equal between the reference and target genes. The slope m = 0 indicates that the efficiencies of the two PCR reactions are equal.

PGK1, RDN5.8, UBC7, and UBC13) had PCR efficiencies that were similar or relatively equivalent to the target amplification efficiencies ( Table 5) .

To test the effect of azole on the expression of pleiotropic drug resistance genes in C. glabrata, we assessed the fluconazole-induced expression of two ABC genes (CDR1 and PDR1) and one ERG gene (ERG4) in five C. glabrata strains, including the PDR1 mutant strain CgB4 (Table 1) . For comparison, we used both RDN5.8 and ACT1 as references for normalization. As shown in Table 6 , fluconazole markedly induced increases in ERG4 mRNA levels in all C. glabrata strains examined when normalized to RDN5.8. Fluconazole also significantly increased CDR1 and PDR1 mRNA expression in all of the strains except CgB4, consistent with the critical role of PDR1 in azole-induced transactivation of ABC transporter gene expression, but not ergosterol biosynthesis gene expression, in C. glabrata. In contrast, when using ACT1 as the reference gene for quantification, fluconazole appeared to down-regulate the expression of all three target genes in the five C. glabrata strains (Table 6 ). Finally, we compared the fluconazole-inducibled mRNA expression levels of four target genes (CDR1, PDR1, ERG4, and ERG10) in C. glabrata (Cg84u strain) after normalizing to different reference genes (ACT1, PGK1, RDN5.8, and UBC13), individually and in pairs. Differences in quantitation were detected according to the reference genes used. As seen in Table 7 , normalization of the RT-qPCR data against the reference genes suggested as optimal by the four software packages (hkgFinder, geNorm, BestKeeper, and NormFinder) or the 2 -ΔΔCT method, gave comparable relative expression levels of the target genes under fluconazole treatment in C. glabrata. However, normalization against ACT1 resulted in relative expression levels of the targets that were substantially different from those normalized using other reference genes, implying that ACT1 is not a suitable reference gene for these studies. Taken together with the data shown above, these results demonstrate that the relative quantification of azoleinducible gene expression varies largely depending on the reference gene and the number of reference genes used for normalization. This highlights the importance of choosing a suitable reference gene or reference gene pair when using RT-qPCR to determine the level of target gene expression in this model system.

In any gene expression study, the selection of a valid normalization or internal control gene to correct for differences in RNA sampling is critical in order to avoid misinterpretation of results and to obtain reliable conclusions. When choosing a reference gene as the internal endogenous control for gene expression studies by RT-qPCR, two important criteria must be met. The expression of the reference gene must remain stable throughout the given intervention (i.e., stability), and the amplification efficiency of the reference gene should be similar to that of the genes of interest (i.e., suitability). In the present study, we used five different methods to evaluate 16 reference genes for potential use as internal controls and found that the reference genes performed differently in terms of stability and suitability in C. glabrata cells upon exposure to fluconazole. To our Values indicate the fold change in RNA transcription for each target gene in fluconazole-treated C. glabrata, as compared with untreated cells.

ΔC T = C T target − C T reference ; ΔC T of each target gene was calculated by using RDN5.8 or ACT1 as the reference. ΔΔC T (t) = ΔC T (t) − ΔC T (ut); ut, untreated; t, fluconazole-treated. Fold change = 2 -ΔΔCT (t). Fold changes of the target genes were also computed with the hkgFinder tool by using RDN5.8 or ACT1 as the reference. # P > 0.05 and *P < 0.05 for the azole-treated group vs. the untreated group after normalization to RDN5.8 or ACT1. *Values indicate the fold change in RNA transcription for each target gene in fluconazole-treated C. glabrata after normalization to one or more reference genes, as compared with untreated cells.

ΔΔC T (t) = ΔC T (t) − ΔC T (ut); ut, untreated; t, fluconazole-treated. Fold change = 2 -ΔΔCT (t). Fold changes of the target genes were also computed with the hkgFinder tool by using one or more reference genes as the internal control for normalization P < 0.05 vs. untreated cells, for all values shown in the table after normalizing to one or more reference genes as indicated.

knowledge, this is the first report to validate reference genes as RNA internal references in C. glabrata.

The poor performance of ACT1 in C. glabrata cells was surprising, given that this gene has been used frequently as the reference gene in earlier gene expression studies [1, 34, 37, 38] . Our data clearly demonstrate the unsuitability of ACT1 as an internal control for gene expression studies in C. glabrata following fluconazole treatment. The initial results gained from using ACT1 as the internal control suggested that target gene expression was not up-regulated (Tables 6 and 7 ). In fact, the only substantial change caused by azole treatment was a greater increase in ACT1 RNA transcription compared with target gene transcription. While these findings are relevant to our specific study, it appears that numerous other studies have also shown the potential of ACT1 to detrimentally affect the accuracy of results [11, 12, 15, 19, 26, 27, 29, 30, 32, 33, [45] [46] [47] . We have been using ACT1 as the internal control for quantitation of gene expression by RT-qPCR in clinical isolates of C. glabrata, and we find that this gene works well as the reference in cells without azole or other agent stimulation [1] . However, Edlind and colleagues used ACT1 as the reference gene for azole-inducible gene expression studies by slot blotting in Candida species, and their data clearly show the variation of ACT1 expression in response to azoles in their systems [34, 38] . Thus, our data revealing the instability of ACT1 in C. glabrata following antifungal treatment, combined with evidence from mammalian and other fungus studies, add to the growing body of evidence that ACT1 expression is unstable across various cell types and under different experimental conditions [11, 12, 15, 26, 27, 29, 30, 32, 33, [45] [46] [47] .

Although ACT1 gene expression was variable in response to fluconazole, the three ribosomal RNAs (i.e., RDN5.8, RDN18, and RDN25) remained unaffected and showed stable expression in azole-treated C. glabrata. These results indicate that ribosomal RNA expression offers superior consistency compared with the expression of ACT1 and the other reference genes assessed. The stable expression levels of 18S and 28S rRNAs relative to other reference genes under a variety of experimental conditions has previously been described for numerous systems, including both mammalian and yeast cells [10, 11, 15, 18, 23, 30, [48] [49] [50] . The levels of ribosomal RNA, which represents 80% of total RNA, are thought to be less likely to vary under conditions that affect the expression of mRNAs because they are transcribed by a distinct RNA polymerase. As an example, Thellin et al. and other groups have recommended the use of 18S or 28S rRNA as an internal control for mRNA quantification studies because mRNA variations are weak and cannot highly modify the total RNA level [10, 11, 15, 18, 30, 39, [48] [49] [50] [51] [52] . However, our further validation experiments showed that of slopes of ΔC T versus log 10 cDNA were sufficiently parallel between RDN5.8 and the target genes, but not between RDN18 or RDN25 and the target genes. With RDN18 and RDN25, the absolute slope values of the ΔC T versus log 10 input cDNA lines were >0.1 for all target genes. These data indicate that the amplification efficiency of RDN5.8 was similar to the efficiencies of the target genes, whereas the amplification efficiencies of RDN18 and RDN25 were different from the target gene amplification efficiencies. This may be attributable to the much higher abundance of RDN18 and RDN25 than RDN5.8 compared with target mRNA transcripts, making it difficult to accurately subtract the baseline value in RT-qPCR data analysis. Therefore, although all three ribosomal RNA subunits were stable during fluconazole stimulation, only RDN5.8 may offer a more accurate and suitable alternative to ACT1 as an internal control for gene expression studies in C. glabrata.

GAPDH, a glycolytic enzyme, is encoded by a single gene and has the advantage of being highly conserved across different species [53, 54] . Like 18S rRNA and β-actin, GAPDH has been commonly used as an internal control, often without testing. In the present study, GAPDH showed much higher variability than any of the ribosomal RNAs in fluconazole-treated samples. These data demonstrate that GAPDH is not an appropriate control gene for these studies, as has been pointed out in previous examples, and that it may lead to incorrect results under specific experimental conditions [12] [13] [14] [15] 19, [29] [30] [31] 33, 55] . Previous studies have indicated the instability of GAPDH in mammalian systems, and this study broadens the scope of this phenomenon to C. glabrata as well.

PGK1 also plays important roles in the glycolytic pathway, and PGK1 and GAPDH are potentially co-regulated [56] . In our data, however, their potential co-regulation was not significant. PGK1 mRNA levels remained relatively stable, in contrast to the marked variation in GAPDH mRNA levels in C. glabrata cells, following fluconazole challenge. Moreover, our comparative C T calculations showed that the efficiency of PGK1 amplification was approximately equal to the efficiencies of the target gene amplifications. Although PGK1 shows some variation as a reference gene, this may not affect experimental results as long as the intergroup difference being measured is greater than the reference gene variation, that is, a reference gene RNA that has an error of 1 log 2 may not be ideal, but it would be sufficient to measure a 2 log 2 change in a gene of interest. Thus, it is inferred that PGK1 may be a suitable reference gene for the analysis of expression for genes with higher azoleinducible mRNA levels, such as ERG4 and ERG10.

Ubiquitin is a small regulatory protein that has been found in almost all tissues of eukaryotic organisms. The UBC gene codes for a polyubiquitin precursor protein [57] . Due to its ubiquitous existence in different tissues and cells in eukaryotes, there are an increasing number of studies in the literature using the UBC gene as the internal standard for gene expression analysis in different eukaryotic cell systems [26, 29] . Out of curiosity, we validated three UBC genes (UBC4, UBC7, and UBC13) in this study. Interestingly, we found that UBC7 and UBC13 mRNAs (particularly the latter) were relatively stable in C. glabrata during fluconazole treatment. In addition, our validation experiments demonstrated that the amplification efficiencies of these genes were approximately equal to those of the target genes. These findings indicate that like PGK1, UBC13 and UBC7 may also be suitable internal controls for quantifying the expression of specific genes with higher azole-inducible mRNA levels in C. glabrata, such as some ergosterol biosynthesis genes.

To successfully select reference genes for our studies involving azole treatment, we also investigated seven other reference genes, in addition to the reference genes mentioned above, with a diversity of functions. These reference genes can be generally classified into several groups: transcription-related genes (EF1α), structure/ cytoskeleton-related genes (TUB1), protein synthesisrelated genes (RPL2A, RPL10, and RPL13A), and finally, genes that cannot be clearly categorized, including PPIA and SDHA. These potential reference genes such as PPIA, RPL13A, and TUB1 are other examples of commonly used internal controls [28, 29, 32] . For example, PPIA has been used as a reference gene because of its remarkable evolutionary conservation and broad cellular and tissue distribution [58] . Although these seven reference genes have been used as internal standards for normalization in countless studies, all of these genes showed an unacceptable variable expression in our model system, with values ranging from a 4.4-fold induction with TUB1 to a 23-fold induction with SDHA after antifungal treatment. Altogether, these results suggest that the choice of internal controls is highly specific to a particular experimental condition, thus highlighting the importance of validating reference genes for each experimental model before commencement of RT-qPCR studies.

Although it is now widely accepted that normalizing to a single reference gene represents a strategy that is simple to use and can control for every stage of the RT-qPCR, some researchers also advocate the use of two or more reference genes, rather than relying on a single RNA transcript [9, 10, 44, 59] . This is a robust method for providing accurate normalization and is consequently preferable when fine measurements are to be made. According to Vandesompele et al. [9] , the purpose of normalization is to remove the sampling difference (such as RNA quantity and quality) in order to identify real gene-specific variation. They provided evidence that a conventional normalization strategy based on a single gene can lead to erroneous normalization. However, it is not always possible to measure multiple reference genes because of limited sample availability and cost. Furthermore, even when multiple genes are chosen, the resolution of the particular assay remains dependent on the variability of the chosen reference genes. As to our case, the geNorm analysis using the geometric mean of the expression of the 16 candidate cDNAs suggested the use of RDN5.8, RDN18, and RDN25 in combination or the combination of these three ribosomal RNAs plus UBC13 and PGK1 as the reference control in the current study. However, the geNorm assessment is based solely on the variability of reference genes and does not take other factors into account. For example, we found that although RDN18 and RDN25 were quite stable, their amplification efficiencies were not equal to the amplification efficiencies of all target genes tested; thus, they may not be suitable as internal controls in our system. Therefore, when multiple reference genes are necessary, we believe that the combination of RDN5.8 plus UBC13 and/or PGK1 would be a better choice for quantitation of gene expression by RT-qPCR in C. glabrata following azole stimulation.

In this study, we evaluated 16 reference genes for potential use as internal controls for RT-qPCR analysis of gene expression in C. glabrata (Cg84u strain) following 2 hours of exposure to fluconazole at 200 μg/ml. To our knowledge, this is the first identification and validation of RDN5.8, UBC13, and PGK1 as the most suitable and stably expressed reference genes among the 16 reference genes tested. Therefore, we recommend the use of RDN5.8, UBC13, or PGK1 alone or the geometric mean of these genes as standards for normalization when analyzing differences in gene expression levels in C. glabrata during antifungal treatment. More specifically, RDN5.8 may be a more suitable reference gene for the analysis of expression for genes with lower azoleinducible mRNA levels, while UBC13 and PGK1 may be better internal controls for quantifying the expression of genes with higher azole-inducible mRNA levels in C. glabrata. In contrast, we demonstrated that 10 reference genes commonly used in published reports, including ACT1, GAPDH, PPIA, RPL13A, and TUB1, had significant differences in their expression upon azole challenge, and thus were not validated as good endogenous controls in this model. As a main conclusion, this study emphasizes the importance of evaluation studies for the selection of the most appropriate internal controls for each experimental model used for quantitative expression studies.

Microarray and molecular analyses of the azole resistance mechanism in Candida glabrata oropharyngeal isolates

Quantitative real-time RT-PCR-a perspective

The real-time polymerase chain reaction

Twenty-five years of quantitative PCR for gene expression analysis

Real-time PCR for mRNA quantitation

Absolute quantification of mRNA using real-time reverse transcription polymerase chain reaction assays

An overview of real-time quantitative PCR: applications to quantify cytokine gene expression

Guideline to reference gene selection for quantitative real-time PCR

Accurate normalization of real-time quantitative RT-PCR data by geometric averaging of multiple internal control genes

Housekeeping genes as internal standards: use and limits

Betaactin-an unsuitable internal control for RT-PCR

beta-Actin and GAPDH housekeeping gene expression in asthmatic airways is variable and not suitable for normalising mRNA levels

GAPDH as a housekeeping gene: analysis of GAPDH mRNA expression in a panel of 72 human tissues

Control selection for RNA quantitation

Direct comparison of GAPDH, beta-actin, cyclophilin, and 28S rRNA as internal standards for quantifying RNA levels under hypoxia

Mechanical injury of cartilage explants causes specific time-dependent changes in chondrocyte gene expression

Gene expression studies in prostate cancer tissue: which reference gene should be selected for normalization?

Effect of experimental treatment on housekeeping gene expression: validation by real-time, quantitative RT-PCR

Validation of housekeeping genes for normalizing RNA expression in real-time PCR

Comparison of human adult and fetal expression and identification of 535 housekeeping/maintenance genes

Validation of endogenous controls for gene expression analysis in microdissected human renal biopsies

Variation in epidermal housekeeping gene expression in different pathological states

Control genes in quantitative molecular biological techniques: the variability of invariance

Quantitative assessment of gene expression in highly purified hematopoietic cells using real-time reverse transcriptase polymerase chain reaction

Unsuitability of using ribosomal RNA as loading control for Northern blot analyses related to the imbalance between messenger and ribosomal RNA content in rat mammary tumors

Normalization of real-time quantitative reverse transcription-PCR data: a model-based variance estimation approach to identify genes suited for normalization, applied to bladder and colon cancer data sets

Reference gene selection for quantitative real-time PCR in chrysanthemum subjected to biotic and abiotic stress

Selection of reference genes for quantitative gene expression normalization in flax (Linum usitatissimum L.)

Selection of housekeeping genes for use in quantitative reverse transcription PCR assays on the murine cornea

Ribosomal 18 S RNA prevails over glyceraldehyde-3-phosphate dehydrogenase and beta-actin genes as internal standard for quantitative comparison of mRNA levels in invasive and noninvasive human melanoma cell subpopulations

Housekeeping gene expression during fetal brain development in the rat-validation by semi-quantitative RT-PCR

Reference gene selection for quantitative real-time PCR analysis in virus infected cells: SARS corona virus, Yellow fever virus, Human Herpesvirus-6. Camelpox virus and Cytomegalovirus infections

Standardization strategy for quantitative PCR in human seminoma and normal testis

Azole resistance in Candida glabrata: coordinate upregulation of multidrug transporters and evidence for a Pdr1-like transcription factor

Candida glabrata PDR1, a transcriptional regulator of a pleiotropic drug resistance network, mediates azole resistance in clinical isolates and petite mutants

Mechanism of increased fluconazole resistance in Candida glabrata during prophylaxis

Upregulation of ERG genes in Candida species by azoles and other sterol biosynthesis inhibitors

Histone deacetylase inhibitors enhance Candida albicans sensitivity to azoles and related antifungals: correlation with reduction in CDR and ERG upregulation

Collaborative comparison of broth macrodilution and microdilution antifungal susceptibility tests

Function of Candida glabrata ABC transporter gene

Application of a master equation for quantitative mRNA analysis using qRT-PCR

Office of Cyber Infrastructure and Computational Biology (OCICB), National Institute of Allergy and Infectious Disease (NIAID)

qBase relative quantification framework and software for management and automated analysis of real-time quantitative PCR data

Determination of stable housekeeping genes, differentially regulated target genes and sample integrity: BestKeeper-Excel-based tool using pair-wise correlations

Effects of fasting on the expression of gastrin, cholecystokinin, and somatostatin genes and of various housekeeping genes in the pancreas and upper digestive tract of rats

Changes in beta-actin mRNA expression in remodeling canine myocardium

Regulation of hypoxanthine phosphoribosyltransferase, glyceraldehyde-3-phosphate dehydrogenase and beta-actin mRNA expression in porcine immune cells and tissues

Housekeeping gene variability in normal and carcinomatous colorectal and liver tissues: applications in pharmacogenomic gene expression studies

The quantification of gene expression in an animal model of brain ischaemia using TaqMan real-time RT-PCR

Validation of oligonucleotide microarray data using microfluidic low-density arrays: a new statistical method to normalize real-time RT-PCR data

Relative quantitative RT-PCR to study the expression of plant nutrient transporters in arbuscular mycorrhizas

Expression profile analysis of the lowoxygen response in Arabidopsis root cultures

Various rat adult tissues express only one major mRNA species from the glyceraldehyde-3-phosphate-dehydrogenase multigenic family

Comparison of glyceraldehyde-3-phosphate dehydrogenase and 28 S-ribosomal RNA gene expression as RNA loading controls for northern blot analysis of cell lines of varying malignant potential

Validating internal controls for quantitative plant gene expression studies

Enzyme co-localization in pea leaf chloroplasts: glyceraldehyde-3-P dehydrogenase, triose-P isomerase, aldolase and sedoheptulose bisphosphatase

Regulatory mechanisms involved in the control of ubiquitin homeostasis

Cyclophilin: distribution and variant properties in normal and neoplastic tissues

A: Evidence based selection of housekeeping genes

Evaluation of reference genes for real-time quantitative PCR studies in Candida glabrata following azole treatment

We sincerely thank Huei-Fung Tsai, Jason Noble, and Bryan Walker for useful discussions and gratefully acknowledge Jason Noble for helpful technical assistance. This research was supported by the Intramural Research Program of the National Institute of Allergy and Infectious Diseases, the National Institutes of Health.

Additional file 1: Summary of the reference genes evaluated in this study.Additional file 2: Primers and probes for RT-qPCR analyses of reference gene RNA transcription in this study. 

The authors declare that they have no competing interests.

The hkgFinder software is currently available as an R source script and it can be downloaded for free from the NIAID Exon website [http://exon.niaid.nih. gov/hkgFinder/] with a sample data set and complete instructions. It requires installation of R 2.11 or higher on computers using the Microsoft Windows operating system for complete compatibility with its graphic user interface (GUI) elements. Experienced R users should be able to use hkgFinder on an Apple Macintosh or Linux/Unix operating system with some reasonable adjustments. An hkgFinder webtool should be available soon on the NIAID Exon website listed above. The geNorm software is available from Biogazelle [www.biogazelle.com/genormplus/] with a free 15day trial download as a part of the qBase PLUS software system. The geNorm software has now been integrated into the qBase PLUS software, where calculation of relative quantities and geNorm analysis are combined in a single program to speed up analysis. Presently, many manual precalculations are not needed, and cross point (CP) values from RT-qPCR can now be directly used for the gene stability analysis using the qBase PLUS software. The qBase PLUS software currently requires Microsoft Windows XP or above or Apple Mac OS X 10.6 (Snow Leopard) with Java 1.6 or later. Support for Linux is also available, but no requirements are listed on the manufacturer's website. The BestKeeper software [http://gene-quantification. com/bestkeeper.html] is available for free download. Please note that the software requires a password generated by an automatic email response from genequan@wzw.tum.de or password@gene-quantification.info. It requires the Microsoft Windows operating system and Microsoft Excel, but specific versions are not listed on the manufacturer's website. The NormFinder software is available for free download from the manufacturer's website [http://www.mdl.dk/publicationsnormfinder.htm]. It requires Microsoft Windows operating system and Microsoft Excel 2003 or above.Authors' contributions QQL conceived of the project, conducted the studies, performed all the experimental procedures, carried out the analysis and interpretation of data, wrote the manuscript, and is the primary author of this paper. JS developed the hkgFinder software and technically helped with the use of other software packages in the present study. JEB participated in the design and coordination of the study and critically reviewed the manuscript. All authors have read and approved the final manuscript.