Abstract
Background
Fluorescent data obtained from realtime PCR must be processed by some method of data analysis to obtain the relative quantity of target mRNA. The method chosen for data analysis can strongly influence results of the quantification.
Results
To compare the performance of six techniques which are currently used for analysing fluorescent data in realtime PCR relative quantification, we quantified four cytokine transcripts (IL1β, IL6 TNFα, and GMCSF) in an in vivo model of colonic inflammation. Accuracy of the methods was tested by quantification on samples with known relative amounts of target mRNAs. Reproducibility of the methods was estimated by the determination of the intraassay and interassay variability. Cytokine expression normalized to the expression of three reference genes (ACTB, HPRT, SDHA) was then determined using the six methods for data analysis. The best results were obtained with the relative standard curve method, comparative Ct method and with DARTPCR, LinRegPCR and Liu & Saint exponential methods when average amplification efficiency was used. The use of individual amplification efficiencies in DARTPCR, LinRegPCR and Liu & Saint exponential methods significantly impaired the results. The sigmoid curvefitting (SCF) method produced medium performance; the results indicate that the use of appropriate type of fluorescence data and in some instances manual selection of the number of amplification cycles included in the analysis is necessary when the SCF method is applied. We also compared amplification efficiencies (E) and found that although the E values determined by different methods of analysis were not identical, all the methods were capable to identify two genes whose E values significantly differed from other genes.
Conclusion
Our results show that all the tested methods can provide quantitative values reflecting the amounts of measured mRNA in samples, but they differ in their accuracy and reproducibility. Selection of the appropriate method can also depend on the design of a particular experiment. The advantages and disadvantages of the methods in different applications are discussed.
Background
Reverse transcription (RT) followed by polymerase chain reaction (PCR) is at present the most sensitive method for the detection of specific RNA molecules. Quantification of nucleic acids using the PCR has been significantly simplified by the development of the realtime PCR technique, where the fluorescent signal reflecting the PCR product accumulation is detected in every amplification cycle. In biological applications examining gene expression, it is mostly not necessary to know the absolute amount of the measured mRNA (number of molecules in a sample). Relative mRNA quantification is an approach determining the amount of target mRNA in samples relative each to other. To compensate for differences in the RTPCR input quality and quantity, the target mRNA amount in each sample is normalized to one or more internal controls. Selection of an optimal normalization strategy has been widely discussed [1,2] and is out of the scope of the present study.
Fluorescent data obtained from realtime PCR must be processed by some method of data analysis to obtain the relative quantity of target mRNA. There are several techniques used for realtime PCR data analysis and adequate attention should be paid to the selection of the appropriate method. Skern et al.[3] demonstrated that quantification results can vary dramatically depending on the method chosen for data analysis, and different analytical approaches may even lead to opposing biological conclusions. Methods for analysis of fluorescent realtime PCR data used in relative mRNA quantification can be classified in various ways, depending on the criteria applied. All methods determine the RTPCR template quantity (designated "R_{0}" throughout the present study) from the accumulation of the PCR product during the amplification process. Most techniques utilize exclusively the exponential phase of PCR to determine the amplification efficiency (designated "E" throughout the present study) and the R_{0 }value [4,5]. The methods can be based on the determination of a "crossing point" between the PCR product fluorescence and a chosen benchmark. The benchmark is a point in the amplification curve (a graph of PCR product fluorescence versus amplification cycle number) that represents the same amounts of PCR product in every amplification. The number of amplification cycles needed to reach the benchmark is usually denoted as CP [6]. The most commonly used form of the benchmark is the threshold fluorescence (which can be set manually by the user or automatically by the software of a realtime PCR instrument), and the number of amplification cycles needed for reaching the threshold fluorescence is usually denoted as "Ct". The basic principle of the "thresholdbased" methods is the same – the lower the RTPCR template amount, the more amplification cycles are needed to reach the threshold fluorescence. Apart from the "thresholdbased" methodologies that currently predominate, there are methods which use linear regression analysis of the fluorescent data from the exponential phase of PCR to determine the E and/or R_{0 }values [7,8]. Moreover, a method that utilizes fluorescent data from the whole course of the amplification curve has been developed [911].
To compare the performance of six techniques which are currently used for analyzing fluorescent data in realtime PCR relative quantification, we determined the mRNA levels of four proinflammatory cytokines (IL1β, IL6 TNFα, and GMCSF) in mice with trinitrobenzene sulphonic acid (TNBS) – induced colitis. TNBSinduced colitis is a widely used experimental model for studying gut inflammatory processes such as ulcerative colitis and Crohn's disease [12].
Results
Relative standard curves, determination of amplification efficiency
Relative standard curves of the analyzed targets are shown in Figure 1. The difference in Ct between duplicate reactions of the standard dilutions did not exceed 1% and it was in the range 0.05% – 0.75% for ACTB and in the range 0.15%–1.0% for other genes. One exception was the 800fold dilution of IL6 template where the difference in Ct significantly exceeded 1% (it reached 2.3%) indicating less reliable quantification at this dilution of the IL6 template. In GMCSF (the gene with the highest Ct values), the 800fold template dilution did not produce a measurable fluorescent signal.
Figure 1. Relative standard curves. The standard curves were generated by the Mx3000P software by plotting cycles at threshold fluorescence (Ct) against the logarithmic values of standard RNA amounts. Quantities of standard RNA were expressed as dilution factors of the RNA preparation (1, 0.5, 0.25, 0.125, 0.0125, 0.00125). Correlation coefficients (square values – r^{2}) and amplification efficiencies (E) are shown. SYBRGreen I fluorescences were corrected with passive reference dye (ROX) fluorescences („dRn“).
Amplification of the added luciferase mRNA showed very similar Ct values in all six dilutions of the standard RNA (arithmetical mean: 19.4, coefficient of variation: 1.28) indicating a similar efficiency of reverse transcription in all dilutions. Moreover, high correlation coefficients of the standard curves (Fig. 1) indicate that both, the efficiencies of reverse transcription and the PCR amplification efficiencies are similar in all dilutions.
Amplification efficiency (E) representative for each gene was determined from equation 7 in the relative standard curve method (Table 1). Amplification efficiency for each individual reaction was calculated using equation 10 in the DARTPCR and LinRegPCR methods; in the Liu & Saintexp method, the amplification efficiency of each reaction was calculated from equation 11. Table 1 shows average amplification efficiencies determined for each gene using the DARTPCR, LinRegPCR and Liu&Saintexp methods. On the whole, E values determined with the LinRegPCR and Liu &Saintexp methods were rather lower and the ones determined with the DARTPCR method were rather higher than the E values determined from relative standard curves. In genes with low(ACTB, IL1β) and medium(SDHA, HPRT, TNFα) Ct values, the average amplification efficiencies determined by the four methods were close to 1.0; in genes with highCt values (IL6, GMCSF), the efficiencies were lower. However, the amplification efficiency of the gene with the highest Ct values (GMCSF, average Ct value from dilutions 1× – 80×: 29.79) was higher than the efficiency of the other „highCt“ gene (IL6, average Ct value from dilutions 1× – 80×: 26.69).
Table 1. Amplification efficiency determined by various methods of realtime PCR data analysis. Arithmetical means ± SD and coefficients of variation (in parentheses) are shown.
Determination of target mRNA quantity in known sample dilutions
Linear regression analysis generated by plotting values of target mRNA quantities (R_{0}) against diluting factors of the total RNA preparation was performed to evaluate the accuracy how the R_{0 }values determined by different methods of realtime PCR data analysis reflect the RTPCR template dilutions. The relationship between determined R_{0 }values and sample dilutions was closest in the standard curve method, comparative Ct method, DARTPCR method with average E values, Liu & Saintexp method with average E values and the combined LinRegPCRCt method using average E values, where average Pearson's correlation coefficients (obtained from the coefficients of the seven genes) reached value at least 0.999 (Table 2). In the sigmoid curvefitting method the average Pearson's correlation coefficient was 0.9953; manual selection of the cutoff cycle (see Methods) was necessary in SDHA, improving the value of the Pearson's correlation coefficient from 0.9986 to 0.9999. In the remaining three methods (LinRegPCR, DARTPCR and Liu & Saintexp methods with individual E values) the average Pearson's correlation coefficients reached values between 0.9577 and 0.9733 (Table 2). In most cases, the Pearson's correlation coefficients in the gene with the highest Ct values (GMCSF) were lower than in other genes.
Table 2. Pearson's correlation coefficients obtained from the linear regression plotting R_{0 }values against diluting factors of the total RNA (RTPCR template). The R_{0 }values were obtained by transformation of fluorescence data using the following methods for realtime PCR data analysis: St cur, relative standard curve; Comp, comparative Ct; SCF, sigmoid curvefitting; DART ind E, DARTPCR with individual E values; DART av E, DARTPCR with average E values; Liu&S ind E, Liu & Saintexp with individual E values, Liu&S av E, Liu & Saintexp with average E values; LinReg ind E, LinRegPCR (using individual E values); LinRegCt av E, LinRegPCR combined with Ct (using average E values).
We normalized the quantity of cytokines (IL1β, IL6, TNFα, GMCSF) with the quantity of one reference gene (ACTB, a "lowCt" gene or HPRT, a "highCt" gene) or with the normalization factor to test the influence of the normalization on the accuracy of the evaluated methods. Normalization factors were determined in all dilutions of the RTPCR template using the geNorm software. The parameter "M", calculated by the software, is a measure of the gene expression stability; the lower is the M value the higher is the gene expression stability. The parameter "V2/V3" (i.e., the third gene added to the two genes) represents pairwise variations determining the optimal number of reference genes which should be included for calculation of a reliable normalization factor [13]. Parameters of the geNorm software calculated for the three reference genes (used for the normalization factor calculation) are shown in Table 3. Values of the M parameter for ACTB, SDHA, and HPRT were similar within all methods of analysis. However, comparison of the M values between individual methods of analysis showed that the M values were 2 – 10 fold higher in the three methods utilizing individual E values (DARTPCR, Liu & Saintexp and LinRegPCR with individual E) than in other methods. Similarly, values of the V2/3 parameter were higher in the three methods utilizing individual E values than in other methods indicating a requirement of higher amount of reference genes to calculate a reliable normalization factor.
Table 3. Parameters of the geNorm software for ACTB, SDHA, and HPRT. The parameters were determined from quantities (R_{0 }values) of the three genes in the RTPCR template dilutions. M – measure of the gene expression stability, V2/V3 – pairwise variations determining the optimal number of reference genes. Indicated methods for realtime PCR data analysis were used for the transformation of fluorescence data to the R_{0 }values. Designation of the methods is the same as in Table 2.
Since the determined quantity of cytokines and reference genes should reflect the template dilutions by the same manner, the same normalized quantity of cytokines should be determined in each dilution. Variability of the normalized quantity between the known template dilutions can then be used as a measure of the quantification accuracy – the lower the variability between dilutions the higher the accuracy of quantification. Figure 2 shows coefficients of variation (CV) for normalized quantity of IL1β in six RTPCR template dilutions (similar results were obtained for IL6, TNFα, and GMCSF; data not shown). In most methods, similar CV values were found after normalization to the normalization factor, to the „low Ct“gene or to the „high Ct“gene (which is caused by similar expression stability of the three reference genes as shown in Table 3); differences detected by the three methods utilizing individual E values were not consistent nor among the three methods neither among the four measured cytokines. On the other hand, comparison of CV values between individual methods of analysis showed markedly higher coefficients of variation in the three methods utilizing individual E values (DARTPCR, Liu & Saintexp and LinRegPCR with individual E) than in other methods (Fig. 2).
Figure 2. Coefficients of variation for normalized quantity of IL1β in serial dilutions of the RTPCR template. Quantities of IL1β, HPRT, SDHA and ACTB mRNAs (the R_{0 }values) were determined in each dilution of the RTPCR template (six dilutions of total RNA) using all the tested methods of realtime PCR data analysis. The amount of IL1β mRNA in each dilution was then divided by the relative amount of HPRT, ACTB, or by the normalization factor (NF, geometric mean of HPRT, SDHA and ACTB amounts) of the dilution. Arithmetical mean, standard deviation and coefficient of variation of the normalized IL1β quantity for each type of normalization (HPRT, ACTB, NF) and each method of realtime PCR data analysis was then calculated. Coefficients of variation (CV) are shown: the first columns are CV values after normalization with the normalization factor, the second columns are CV values after normalization with ACTB, and the third columns are CV values after normalization with HPRT. Methods for realtime PCR data analysis: St. C., relative standard curve; COM, comparative Ct; SCF, sigmoid curvefitting; DAR iE, DARTPCR with individual E values; DAR aE, DARTPCR with average E values; L&S iE, Liu & Saintexp with individual E values, L&S aE, Liu & Saintexp with average E values; LR iE, LinRegPCR (using individual E values); LRCt, LinRegPCR combined with Ct (using average E values)
Intraassay and interassay variability
Reproducibility of the tested methods for realtime PCR data analysis was evaluated by performing experiments where the intra and interassay variability was measured. Two genes were selected for the measurements – IL1β, representing a target with low Ct values and IL6, representing a target with high Ct values. Results from the experiment measuring the intraassay variability are shown in Table 4. In both genes, the lowest coefficients of variation were obtained when the following methods of realtime PCR data analysis were used for the determination of target mRNA quantity (R_{0}): standard curve method, comparative Ct method, DARTPCR method with average E values, Liu & Saintexp method with average E values and the combined LinRegPCRCt method (using average E values). Using the sigmoidal curvefitting method, the coefficients of variation were markedly higher. The highest coefficients of variation were obtained using the LinRegPCR method, DARTPCR method and Liu & Saintexp method with individual E values. Similar results (except of the SCF method) were obtained in the experiment measuring the interassay variability (Table 5). Coefficients of variation calculated from Ct values of 15 replicate PCR reactions in the intraassay experiment were 0.52 for IL1β and 0.63 for IL6; in the interassay experiment, the coefficients were 1.1 for IL1β and 1.18 for IL6.
Table 4. Intraassay variability. Coefficients of variation calculated from R_{0 }values of 15 replicate PCR reactions. Indicated methods for realtime PCR data analysis were used for the transformation of fluorescence data to the R_{0 }values. Designation of the methods is the same as in Table 2.
Table 5. Interassay variability. Coefficients of variation calculated from R_{0 }values of 8 PCR reactions. Indicated methods for realtime PCR data analysis were used for the transformation of fluorescence data to the R_{0 }values. Designation of the methods is the same as in Table 2.
In the standard curve method, comparative Ct method, and the three methods using average E values (DARTPCR, Liu & Saintexp, LinRegPCRCt – with average E), the intraassay variability was higher in IL6 than in IL1β; a similar trend was found for the interassay variability (except of two methods – DARTPCR and Liu & Saintexp, which showed similar values for both genes). The three methods (DARTPCR, Liu & Saintexp, LinRegPCRCt) utilizing individual E values (and showing high CV values) gave inconsistent results and in some cases, the determined intraassay variability was even higher than the interassay variability. The CV values obtained in the intra and interassay experiment with the SCF method were comparable, and they were higher in IL6 than in IL1β which is in accordance with the five wellworking methods (Table 4, Table 5).
Normalized cytokine expression in the experimental samples
Figure 3 and Figure 4 show normalized expression of IL1β and IL6 as determined by all the tested methods of realtime PCR data analysis in the experimental samples. Estimation of IL1β and IL6 expression using the KruskalWallis test showed that differences between the experimental groups were highly significant (P ≤ 0.01 using the LinRegPCR method with individual E values and P ≤ 0.001 with all other methods). Estimation of the difference in IL1β and IL6 expression between the group Un (untreated colitic animals) and the other groups of animals (with the MannWhitney test) showed significantly increased expression of these cytokines in the untreated colitic animals, using all the tested methods of realtime PCR data analysis. All methods also detected a significant decrease of the two cytokines expression after treatment B (with lower level of significance using LinRegPCR or DARTPCR methods with individual E values). In the case of IL1β, some of the methods also detected significant decrease in the cytokine expression after treatments A, or C.
Figure 3. IL1β mRNA expression. Quantities of IL1β, HPRT, SDHA and ACTB mRNAs (the R_{0 }values) were determined in each sample using all the tested methods of realtime PCR data analysis. The amount of IL1β mRNA in each sample was then divided by the normalization factor (geometric mean of HPRT, SDHA and ACTB amounts) of the sample. Values are arithmetical means + SEM, n = 5–8. Statistical significance of the differences between the group of untreated colitic animals (Un) and other groups of animals was assessed with the MannWhitney test: * P ≤ 0.05, ** P ≤ 0.01. Methods for realtime PCR data analysis: Stand curv, relative standard curve; Comparat, comparative Ct; SCF, sigmoid curvefitting; DART ind E, DARTPCR with individual E values; DART av E, DARTPCR with average E values; Liu&S ind E, Liu & Saintexp with individual E values, Liu&S av E, Liu & Saintexp with average E values; LinReg ind E, LinRegPCR (using individual E values); LinRegCt av E, LinRegPCR combined with Ct (using average E values) Designation of the animal groups: the first columns, control sham animals (Sh); the second columns, untreated colitic animals (Un); the third columns, colitic animals with the treatment A; the fourth columns, colitic animals with the treatment B; the fifth columns, colitic animals with the treatment C
Figure 4. IL6 mRNA expression. Quantities of IL6, HPRT, SDHA and ACTB mRNAs (the R_{0 }values) were determined in each sample using all the tested methods of realtime PCR data analysis. The amount of IL6 mRNA in each sample was then divided by the normalization factor (geometric mean of HPRT, SDHA and ACTB amounts) of the sample. Values are arithmetical means + SEM, n = 5–8. Statistical significance of the differences between the group of untreated colitic animals (Un) and other groups of animals was assessed with the MannWhitney test: * P ≤ 0.05, ** P ≤ 0.01. Designation of the methods and animal groups is the same as in Figure 3.
In TNFα and GMCSF, significant differences were found only between the group of untreated colitic animals (Un) and the group of sham control animals (Sh). In TNFα, all the tested methods of realtime PCR data analysis showed significantly higher amount of the cytokine in the group Un than in the group Sh. In GMCSF, the difference between the group Un and Sh was detected as significant only with the standard curve method, comparative Ct method, and the three methods utilizing average E values (DARTPCR, Liu & Saintexp, LinRegPCRCt – with average E); the methods which used individual E values (DARTPCR, Liu & Saintexp, LinRegPCRCt – with individual E) did not detect the difference as statistically significant (data not shown).
Discussion
To compare the methods currently used for analyzing fluorescence data in realtime PCR relative quantification, we determined mRNA levels of four proinflammatory cytokines in the in vivo model of colonic inflammation, using six different techniques. Three housekeeping genes were also quantified to serve as reference genes for the normalization of cytokine mRNA quantity. The compared methods differ principally in the mathematical function used for modeling the PCR process, in the necessity to create a dilution series of the RTPCR template, in the necessity and means of amplification efficiency (E) determination, in the way of calculation of the target mRNA quantity (R_{0}).
Most of the methods utilize fluorescence data only from the exponential phase of PCR amplification (the exponential model of PCR) and require setting a threshold fluorescence (common for all compared samples) to determine the Ct value ("thresholdbased" methodologies). The relative standard curve method determines target mRNA relative quantities in samples from known relative quantities of standard RNA or cDNA. To obtain correct results, amplification efficiency in the dilutions of the standard preparation and in samples must be similar. The E value determined from the slope of the standard curve (or from a dilution series of a representative sample) can be utilized for the calculation of the target mRNA quantity (R_{0}) in techniques of relative quantification, which can be denoted as "comparative Ct methods". Equation 5 can be used for the transformation of Ct values to the R_{0 }values; the only parameter which differs among the compared samples is Ct. The term "comparative Ct method" is most frequently used for the " 2^{ΔΔCt }" method introduced by Livak [14]; in this method E is assumed to be equal to 1 and the formula for the R_{0 }calculation can be modified to R_{0 }= 2^{Ct}. Pfaffl's model of the comparative Ct method [15] incorporates a correction for amplification efficiency differing from the optimal value 1. Most published versions of the comparative Ct method use normalization of the target gene quantity to a reference gene quantity, very often referring to one of the samples (a "calibrator" sample), and the formula for calculation of the target gene relative quantity is more complicated [1419]. In any case, formulas used in comparative Ct methods for calculation of the normalized relative amount of target gene quantity can be derived from equation 5. Pfaffl et al.[20] developed the software named REST (Relative Expression Software Tool) which can perform the comparative Ct quantification in two experimental groups (with or without the E value correction) followed by a statistical test. In the present study, we used E values determined for each gene from the relative standard curves and calculated the R_{0 }value for each gene in each sample using the equation 5. The cytokine genes quantity was then normalized with the factor obtained from the quantities of three reference genes.
In the method proposed by Liu and Saint [21], the amplification efficiency for each sample is determined from the amount of fluorescence and the number of cycles at two arbitrary fluorescence thresholds along the exponential phase of the PCR amplification (Eq 11). The fluorescence thresholds can be set individually for each amplification, or a threshold level common for all reactions can be set. In the first case equation 4 is used for the R_{0 }calculation, and in the second case equation 4 or equation 5 can be used. We compared the quantification results obtained with this method using either individual E values (determined for each amplification) or using an average E value.
Ramakers et al.[7] developed a computer program entitled LinRegPCR, which determines the target mRNA quantity (R_{0}) and amplification efficiency (E) by linear regression analysis (Eq 9, Eq 10) of fluorescence data obtained from realtime PCR. Like the abovementioned methods, linear regression exploits only the exponential phase of PCR amplification, but the method is not "thresholdbased" – no benchmark (threshold fluorescence) is needed for the calculations. We analyzed our fluorescence data with the LinRegPCR software and with a technique which combines linear regression analysis and the thresholdbased methodology. The combined technique (we have designated it "LinRegPCRCt") utilizes E values determined by the LinRegPCR software (for calculation of an average E) and Ct values determined by the Mx3000P realtime PCR instrument. A similar strategy was applied by Karlen et al.[22] and Schefe et al.[23]. The DARTPCR program developed by Peirson et al.[8] provides an automated analysis of realtime PCR fluorescence data utilizing the combined approach (linear regression for E determination and threshold fluorescence for Ct determination). Similarly as for LinRegPCR, we compared the quantification results obtained with DARTPCR using individual or average E values for the R_{0 }calculation.
The method in which a sigmoid mathematical model that fits the kinetics of the whole realtime PCR process is applied [9,10], represents an approach completely different from the abovementioned techniques. The sigmoid curvefitting (SCF) method utilizes all fluorescence data recorded during the amplification process (not only the data from the exponential phase) for determination of the R_{0 }value. Moreover, the method can carry out quantification without the knowledge of amplification efficiency and without determination of Ct. Rutledge [11] found that amplification cycles within the plateau phase of PCR deviate from that predicted by sigmoid curvefitting, and their exclusion from the curvefitting process is necessary. He proposed the selection of a cutoff cycle beyond which further cycles are excluded from the fitting of the amplification curve. The criterion used for the selection of the cutoff cycle was based on repetitive curvefitting in which the last cycle was sequentially excluded and the R_{0 }value was calculated at each individual curvefitting. Plotting the calculated R_{0 }values against the cutoff cycle revealed a highly regular pattern in which the calculated R_{0 }value decreased with subsequent cycle removal, and after reaching a minimum a small increase in the R_{0 }value appeared; the minimumcalculated R_{0 }value was selected as the resulting R_{0}. In our work, we found that the shape of the graph of R_{0 }dependency on the cutoff cycle can be influenced by fluorescence data provided by the realtime PCR system. Using backgroundsubtracted fluorescence data from the Mx 3000P system (Stratagene), we were not able to identify the regular trend of the curve R_{0 }vs cutoff cycle as described by Rutledge [11]; backgroundsubtracted data from Opticon2 DNA Engine (MJ research Inc.) were used in Rutledge's study. On the other hand, raw fluorescence data from Mx 3000P provided a regular pattern of the R_{0 }value as a function of the cutoff cycle (R_{0 }decrease followed by a small increase after reaching a minimum) in most curves. In some instances the cutoff cycle with the minimal R_{0 }value was determined in the region of the amplification curve which contained an insufficient amount of fluorescent data, and manual selection was necessary. Karlen et al.[22] testing the performance of several methods for realtime PCR data analysis found the sigmoid curvefitting method (together with the method fitting PCR amplification to the exponential function) as the least suitable for quantitative PCR analysis. On the contrary, our results indicate that the SCF method can provide reasonable results. Similarly, Qiu et al.[24] obtained comparable results using the SCF method and a classic thresholdbased method. The differences in the SCF method performance were probably caused by different number of amplification cycles included into the fitting process. An appropriate selection of the optimal cycle number (exclusion of late cycles) is probably the key factor for obtaining satisfactory performance of the SCF method, but the choice of a suitable criterion for determination of the "cutoff cycle" can be difficult (as discussed above).
In our study, we determined amplification efficiencies (E) of seven genes using four methods of realtime PCR data analysis, and found some differences in the determined E values. On the other hand, all the methods were capable to identify the two genes whose E values significantly differed from the others. Interestingly, the amplification efficiencies determined by the two methods which employ linear regression for the calculation (DARTPCR and LinRegPCR) were less close each to other than to efficiencies determined by the relative standard curve or Liu & Saintexp method. This can be caused by differences in the way the two methods determine the exponential phase of amplification. Lower amplification efficiency found by all tested methods in two genes with high Ct values (IL6 and GMCSF) could suggest some influence of Ct value on the determination of the E value. However, the comparison of E and Ct values in the two „highCt genes“ indicates that a higher Ct value is not necessarily leading to obtaining of a lower E value. This finding is in accordance with results of Karlen et al.[22] who did not find a dependency of E value on Ct value; the authors defined the amplicon and primer sequences as the main factor influencing the efficiency of amplification.
To compare the accuracy of relative quantification conducted with application of different methods for realtime PCR data analysis, we determined quantities of seven target mRNAs in serially diluted preparation of total RNA. We found that R_{0 }values determined with the relative standard curve method, comparative Ct method and with the three methods using an average E value for the calculations (DARTPCR, Liu & Saintexp, LinRegPCRCt – with average E) most accurately reflected the RTPCR template dilutions. Less effective was the SCF method, and the worst results were obtained with the three methods which used individual E values (DARTPCR, Liu & Saintexp, LinRegPCR – with individual E). Normalization of cytokines mRNA quantity (to a single reference gene or to a normalization factor calculated from reference genes with comparable expression stability) did not influence the differences in the performance of the methods for realtime PCR analysis tested in our study. On the other hand, our results indicate that targets with high Ct values can be quantified with a lower accuracy than targets with medium and low Ct values. Reproducibility of the tested methods was estimated by determination of the intra and interassay variability and showed the same result as the accuracy test. The highest reproducibility was found in the relative standard curve method, comparative Ct method and in the three methods using an average E value for calculations. The SCF method was less precise, and the worst results were obtained with the three methods which used individual E values. Our results also showed a negative effect of higher Ct values on the reproducibility of the tested methods.
Comparing results of normalized expression of IL1β and IL6 in the experimental samples, all the tested methods were able to detect significant changes between the control animals, untreated colitic animals, and animals undergoing treatment B. The effects of treatments A or C on the IL1β mRNA level were detected as statistically significant only by some of the tested methods. In TNFα and GMCSF, all the tested methods showed higher amount of the cytokines in untreated colitic animals than in control animals, hovewer the three methods which utilized individual E values (DARTPCR, Liu & Saintexp, LinRegPCRCt – with individual E) were not able to detect the difference in GMCSF expression as statistically significant.
Application of corrections for individual sample efficiency should theoretically improve the accuracy and reproducibility of the quantification. But our results showed that, on the contrary, the use of individual E values for the R_{0 }calculation impaired the quantification. Similar findings were presented in studies where linear regression was utilized for the calculation of individual amplification efficiencies [8,22,23]. We used two approaches for determination of the individual amplification efficiencies – linear regression (in DARTPCR method and LinReg PCR method which differ in the way of exponential phase determination) and the method which utilize setting of two fluorescence thresholds along the exponential phase (Eq 11). Independently on the way used for the calculation of E values, all the three methods showed the negative effect of the individual E values on the quantification which indicates that this is probably a feature connected with the limited precision of individual data and not with the mathematical approach used for E value determination. The individual E values (determined with linear regression or with the Eq 11) are derived from individual sample kinetics represented by the fluorescence values obtained in particular amplification cycles. Above mentioned results suggest that fluorescence values detected by realtime systems do not reflect the reaction kinetics with the precision which would be sufficient for reliable E value determination from the exponential phase of the individual amplification. The average E value obtained from the group of amplifications eliminates individual imprecisions enabling to find a reliable value of the amplification efficiency.
Conclusion
Choosing the appropriate method for realtime PCR data analysis can depend on conditions in a particular application. The relative standard curve method is widely used and can provide reliable results, especially in the case when the same quality of RTPCR template is ensured in standards and samples. Sometimes this is not possible, for instance if the quantification is performed on different tissue types, or if the amount of tissue is limited (e.g. tissue biopsies, preimplantation embryos). Similar restrictions apply also for comparative Ct methods. These methods require serial dilutions of a representative sample to determine the E value which is usually done in a validation experiment preceding a series of measurements. This approach enables to measure more samples in a PCR run (no standard curve in the run), but identical conditions for all measurements must be ensured. The other three methods tested in the present study (Liu & Saintexp, DARTPCR, LinRegPCR) do not require serial dilutions of the RTPCR template for E value determination. We found that these methods provide more reliable results when an average E value and not individual E values (determined for each amplification) are used for the R_{0}calculation. The Liu & Saintexp method is simple – the formula for E value calculation can be easily implemented into spreadsheet programs such as Microsoft Excel and the selection of the exponential phase of amplification is not difficult (the human eye is good enough for distinguishing a straight line). Our results with the Liu & Saintexp quantification showed that although the criterion for selection of the exponential phase is subjective, reasonable results can be obtained if the same criterion is applied for all compared samples. The LinRegPCR and DARTPCR methods use a more complicated calculation (based on linear regression) for the E value determination. The DARTPCR software combining the linear regression analysis with thresholdbased methodology enables R_{0 }values to be calculated using an average E value; the LinRegPCR software do not enable the automated use of an average E value, but manual combination with a thresholdbased technique is possible. The SCF method differs from all the other methods in the mathematical model used for the calculations – it is not necessary to look for the exponential phase of amplification and to determine the E value. For reliable quantification with this method, fluorescence data including at least the beginning of the plateau phase are needed, which can be a disadvantage when genes with low expression are quantified or when low sample amounts are available. In summary, our results show that all the tested methods for realtime PCR data analysis can provide quantitative values reflecting the amounts of measured mRNA in samples, but they differ in their accuracy and reproducibility. Although selection of the appropriate method can be limited by the design of a particular experiment (e.g. tissue type, gene abundance, number of experimental groups) the use of more than one analytical method is recommended for validation of results.
Methods
Animal experiment, sample preparation, and realtime PCR
All the fluorescent data used in this study were obtained from a previous study examining the effects of plant essential oils thyme and oregano on trinitrobenzene sulphonic acid (TNBS)induced colitis in mice [25]. Briefly, colitis was induced in male Balb/c mice by administration of TNBS and the animals were treated with three increasing concentrations of the plant oils (treatments A, B and C; the group of untreated colitic animals was designated as Un and sham control group as Sh in the present study). Total RNA was isolated from the strips of colonic tissue with TRIzol Reagent (Invitrogen Life Technologies, Karlsruhe, Germany). The RNA preparations were then cleaned and DNase I treated with an RNeasy Micro kit (Qiagen, Hilden, Germany). Complementary DNA (cDNA) was then synthesized from the RNA (0.75 μg from each sample) using Superscript™ II Rnase HReverse Transcriptase (Invitrogen Life Technologies; for more details see [25])
Serial dilutions of total RNA (1× – i.e. no dilution, 2×, 4×, 8×, 80×, 800×) were prepared from the pool of colon RNA obtained by combining aliquots of samples from the colitic animals. Complementary DNA was then synthesized from each dilution as described above. To compensate for different RNA amounts in the reverse transcription reactions yeast total RNA was added in appropriate amounts to the colon RNA dilutions (all oligonucleotide primers were checked so as not to create any PCR product on the yeast cDNA template). The cDNAs then served as a PCR template for construction of relative standard curves. Amplifications performed on the standard cDNAs were also utilized for comparison of the accuracy of quantification by the tested methods for realtime PCR data analysis (determination of known sample dilutions); in the relative standard curve method two sets of total RNA dilutions (and cDNA preparations) were used – one set was utilized for the construction of standard curves and the other set served as known sample dilutions.
To test the efficiency of reverse transcription in the serial dilutions of standard RNA, identical amount (10 pg) of luciferase mRNA (Promega, Madison, WI) was added into the six RNA dilutios (1×, 2×, 4×, 8×, 80×, 800×) and then cDNA was synthesized from each dilution (as described above).
PCR reactions were carried out in duplicates using SYBRGreen I as a fluorescent detection dye and they were performed in the realtime PCR system Mx 3000P (Stratagene, La Jolla, CA). Backgroundsubtracted fluorescences were used for data analysis by all tested methods except for the sigmoid curvefitting method, where raw fluorescences were used (see below). Specific oligonucleotide primers for amplification of mouse interleukin 1 beta (IL1β), interleukin 6 (IL6), tumor necrosis factor alpha (TNFα), granulocyte macrophagemacrophage colony stimulating factor (GMCSF), beta actin (ACTB), hypoxanthine guanine phosphoribosyl transferase 1 (HPRT) and succinate dehydrogenase complex subunit A (SDHA) were used (for sequences of the primers and PCR reaction conditions see [25]).
Amplification of luciferase mRNA: 1 μl of each cDNA (six cDNAs was prepared, see above) was amplified in 25 μl PCR containing 1× SYBR Green/ROX PCR Master Mix (PA012, SuperArray Bioscience Corp., Frederick, MD), and 0.4 μM primers; 40 amplification cycles at 95°C for 20 s, 60°C for 60 s, and 82°C for 20 s (fluorescence acquiring) was used. Primers specific for the luciferase sequence (5'GCTTACTGGGACGAAGACGAAC3', 5'CTTGACTGGCGACGTAATCCAC3') amplified a 247 bp PCR product.
Normalization of cytokine mRNAs quantity
To ensure correctness of the quantification we normalized cytokine (IL1β, IL6 TNFα, and GMCSF) expression to the three reference genes (ACTB, SDHA, and HPRT) whose expression was found to be stable in our previous work [25]. After determining quantities of mRNAs (R_{0 }values) of the three reference genes in each sample (using tested methods of realtime PCR data analysis) the sample normalization factor was calculated as a geometric mean of the three R_{0 }values. GeNorm software was utilized for the calculation [13]. The amount of cytokine mRNA in each sample (cytokine R_{0 }value) was then divided by the normalization factor of the sample.
Data analysis
The mathematical equations used in most methods for analyzing data obtained from realtime PCR are derived from the basic formula describing the PCR amplification in the exponential phase of the reaction:
where X_{n }is the amount of PCR product at cycle n, X_{0 }is the starting amount of PCR template (which we are interested in) and E is the amplification efficiency which can have a value between 0 (no amplification) and 1 (doubling of the PCR product in each amplification cycle). In fluorescent realtime PCR it is assumed that accumulation of reporter dye fluorescence (R, fluorescence readings after background subtraction) is proportional to the accumulation of PCR amplification product, and equation 1 can then be written as:
and starting fluorescence R_{0 }can be then calculated as:
where R_{n }is the intensity of reporter dye fluorescence (proportional to the amount of PCR product) at cycle n, and R_{0 }is the theoretical starting fluorescence which is proportional to the amount of starting PCR template. Thus, the R_{0 }value represents the target quantity expressed in arbitrary fluorescence units. There are several techniques based on this equation which are used for calculating the R_{0 }value. In thresholdbased techniques where "Ct" (the number of amplification cycles needed to reach the fluorescence threshold) is measured, "n" in Eq 3 can be replaced by "Ct":
"R_{Ct}" then represents the threshold fluorescence which can be set for each of the compared amplifications individually, or a threshold value (R_{Ct}) common for all compared amplifications can be used. In the latter case, the numerical value of R_{Ct }in Eq 4 can be ignored (replaced for example by 1.0) and the R_{0 }value can be calculated:
Other possibilities to obtain the R_{0 }and E values are techniques which use relative standard curve, linear regression (utilizing equations derived from the basic formula – Eq 1) or fitting the PCR process to the sigmoid function (see below). The quantities determined from the relative standard curve as well as the quantities determined with the other methods were designated as „R_{0}“ throughout the manuscript.
Relative standard curve method
For each gene, standard cDNAs (see above) were amplified along with sample cDNAs in the same PCR run. Standard curves were generated by Mx 3000P 2.0 software (Stratagene). The threshold fluorescence common for all compared samples was set into the exponential phase of the amplifications by the Mx 3000P system. The target mRNA quantity in each sample (R_{0}) was determined from the relative standard curve (using sample Ct values) and expressed in arbitrary units corresponding to the dilution factors of the standard RNA preparation. Amplification efficiency (E) representative for each gene was determined using equation of the standard curve:
Comparative Ct method
Amplification efficiency (E) for each gene was determined from the relative standard curve (see above). The Ct value for each reaction was determined by the realtime PCR system Mx3000P setting the threshold fluorescence (common for all compared samples) into the exponential phase of the amplifications). The target mRNA quantity in each sample (R_{0}) was then calculated from equation 5.
LinRegPCR method
A computer program entitled LinRegPCR developed by Ramakers et al.[7] utilizes linear regression analysis of fluorescence data from the exponential phase of PCR amplification to determine the target mRNA quantity (R_{0}) as well as the amplification efficiency (E). Following equations were used:
LinRegPCR software utilizes an iterative algorithm (considering the number of data points, regression coefficient and slope of the regression line) for the selection of the exponential phase in each PCR amplification. We analyzed our fluorescence data using this software and obtained the R_{0 }and E values for each reaction. Since ANOVA detected no significant differences in amplification efficiencies between the sample groups (Sh, Un, A, B, C) we also applied a combined analysis using an average E value (arithmetical mean of E values of all samples). In the combined analysis ("LinRegPCRCt method"), the R_{0 }was calculated from equation 4 using the average E value and threshold fluorescence values (R) with corresponding Ct values (determined by the realtime PCR system Mx3000P).
DARTPCR method
DARTPCR (Data Analysis for RealTime PCR) Excel workbook developed by Peirson et al.[8] determines the E value of each individual reaction using linear regression analysis of the fluorescence data from the exponential phase of each amplification (Eq 10). For the selection of the exponential phase a midpoint (M) for each PCR amplification is calculated, using maximal and minimal fluorescence levels. DARTPCR determines the Ct value for each reaction, offering the possibility of using individual threshold fluorescences or a threshold fluorescence common for all compared reactions (the second possibility was used in this study). Target mRNA quantity (R_{0}) is then determined using equation 4. Oneway analysis of variance (included as a component of DARTPCR workbook) detected no significant differences in amplification efficiencies between the sample groups (Sh, Un, A, B, C), so we also used the average E value (arithmetical mean of E values of all samples) for the R_{0 }calculation in each sample. Since the DARTPCR requires inputting of fluorescence data in triplicates and our PCR reactions were carried out in duplicates, we created the "third replicate" as the mean value of our duplicate fluorescence readings.
Liu and Saint method using exponential model of PCR ("Liu & Saintexp method")
Two arbitrary fluorescence thresholds (lower R_{1 }and higher R_{2}, common for all compared samples) were manually set in the exponential phase of amplification curves and corresponding Ct values (provided by the Mx 3000P) were recorded. For the thresholds selection the semilogarithmic graph of the amplification curves (log of fluorescences – log R, against cycle number – n) was used because the exponential phase of PCR amplification can be simply identified on the graph (fluorescent data acquired in the exponential phase of the reaction produce a straight line on the semilogarithmic graph). Amplification efficiency (E) of each reaction was calculated using equation:
The target mRNA quantity (R_{0}) was calculated using equation 4. Oneway analysis of variance (ANOVA) detected no significant differences in amplification efficiencies between the sample groups (Sh, Un, A, B, C), so the average E value (arithmetical mean of E values of all samples) was also used for the R_{0 }calculation in each sample (using equation 4).
Sigmoid curvefitting (SCF) method
Raw fluorescences (i.e. fluorescences without background subtraction) from Mx 3000P were fitted to the fourparametric sigmoid function using the nonlinear regression function of SigmaPlot (Version 10, Systat Software, Richmond, CA, USA). Following equations were used:
where R_{na }is the aggregate reaction fluorescence at cycle n, R_{b }is the background reaction fluorescence, R_{n }is the fluorescence generated by the PCR product at cycle n, R_{max }is the maximal fluorescence generated by the PCR product, n_{1/2 }is the cycle number at which fluorescence reaches half of the R_{max}, k describes the slope of the sigmoid curve.
Repetitive regression analyses with sequential removal of the last amplification cycle (the cutoff cycle) were performed until the curve fitting failed, due to insufficient data. The R_{0 }value was calculated at each curvefitting using equation 13. For all amplification curves, the graph of dependence of the calculated R_{0 }values on the cutoffcycle was constructed and the minimum R_{0 }value was selected as the resulting R_{0 }value – the target mRNA relative quantity [11]. For the data treatment a macro for SigmaPlot provided by Qiu et al.[24] was utilized. In some instances the cutoff cycle with the minimal R_{0 }value was located in the linear region of the amplification curve. The resulting R_{0 }value was then selected manually by shifting the cutoff cycle into the region located between end of the linear amplification and entry of the reaction into the plateau phase (i.e. into "upper arc" of the amplification curve).
Intraassay and interassay variability
In the intraassay variability experiment, fifteen replicate PCR reactions amplifying IL1β or IL6 were ran. In the interassay variability experiment, eight PCR reactions amplifying IL1β or IL6 were ran in eight separate PCR (Mx 3000P) runs. The fluorescent data obtained from each reaction replicate were then transformed into mRNA quantity (the R_{0 }value) by the tested methods for realtime PCR data analysis. Arithmetical means, standard deviations (SD) and coefficients of variation (CV) were then calculated from the determined R_{0 }values (15 and 8 replicates for each gene, respectively).
Statistics
All statistics were performed using Statistica (StatSoft, Tulsa, OK). To enable better comparison of the results obtained from the different methods of realtime PCR data analysis the values of target mRNA quantities (R_{0}) were transformed percentually (setting the median of R_{0 }values of all samples as 100%). Differences in amplification efficiencies between the sample groups were assessed using oneway analysis of variance (ANOVA). The KruskalWallis test was used for the comparison of differences in normalized cytokine expression between groups, and the MannWhitney U test was used to compare differences between the group of untreated colitic animals (Un) and other groups of animals. Values of P < 0.05 were considered as significant.
Authors' contributions
SC conceived and designed the study, analyzed and interpreted the data, and drafted the manuscript. AB participated in the fluorescent data analysis and interpretation, and in the manuscript preparation. JK participated in the data interpretation and manuscript preparation. All authors read and approved the final manuscript.
Acknowledgements
This work was supported by the Slovak Research and Development Agency under Contract No. APVT51015404.
References

Bustin SA: Quantification of mRNA using realtime reverse transcription PCR (RTPCR): trends and problems.
J Mol Endocrinol 2002, 29:2339. PubMed Abstract  Publisher Full Text

Huggett J, Dheda K, Bustin SA, Zumla A: Realtime RTPCR normalisation; strategies and considerations.
Genes Immun 2005, 6:279284. PubMed Abstract  Publisher Full Text

Skern R, Frost P, Nilsen F: Relative transcript quantification by quantitative PCR: Roughly right or precisely wrong?
BMC Mol Biol 2005, 6:10.
DOI 10.1186/14712199610
PubMed Abstract  BioMed Central Full Text  PubMed Central Full Text 
Pfaffl MW: Quantification strategies in realtime PCR. In AZ of Quantitative PCR. Edited by Bustin SA. International University Line (IUL), La Jolla; 2004:86120.

Wong ML, Medrano JF: Realtime PCR for mRNA quantitation.
BioTechniques 2005, 39:7585. PubMed Abstract

Rasmussen R: Quantification on the LightCycler. In Rapid cycle realtime PCR Methods and applications. Edited by Meuer S, Wittwer C, Nakagawara K. SpringerVerlag, Heidelberg; 2001:2134.

Ramakers Ch, Ruijter JM, Deprez RHL, Moorman AFM: Assumptionfree analysis of quantitative realtime polymerase chain reaction (PCR) data.
Neuroscience Letters 2003, 339:6266. PubMed Abstract  Publisher Full Text

Peirson SN, Butler JN, Foster RG: Experimental validation of novel and conventional approaches to quantitative realtime PCR data analysis.
Nucleic Acids Res 2003, 31:e73.
DOI 101093/nar/gng073.
PubMed Abstract  Publisher Full Text  PubMed Central Full Text 
Liu W, Saint DA: Validation of a quantitative method for real time PCR kinetics.
Bioch Biophys Res Comm 2002, 294:347353. Publisher Full Text

Tichopad A, Dzidic A, Pfaffl MW: Improving quantitative realtime RTPCR reproducibility by boosting primerlinked amplification efficiency.
Biotech Lett 2002, 24:20532056. Publisher Full Text

Rutlege RG: Sigmoidal curvefitting redefines quantitative realtime PCR with the prospective of developing automated highthroughput applications.
Nucleic Acids Res 2004, 32:e178.
DOI 101093/nar/gnh177.
PubMed Abstract  Publisher Full Text  PubMed Central Full Text 
Panes J: Inflammatory bowel disease: pathogenesis and targets for therapeutic interventions.
Acta Physiol Scand 2001, 173:159165. PubMed Abstract  Publisher Full Text

Vandesompele J, De Peter K, Pattyn F, Poppe B, Van Roy N, De Paepe A, Speleman F: Accurate normalization of realtime quantitative RTPCR data by geometric averaging of multiple internal control genes.
Genome Biol 2002, 3:RESEARCH0034. PubMed Abstract  BioMed Central Full Text  PubMed Central Full Text

Livak KJ: ABI Prism 7700 Sequence Detection System. User Bulletin no. 2. [http://docs.appliedbiosystems.com/pebiodocs/04303859.pdf] webcite
PE Applied Biosystems 1997, AB website, bulletin reference: 4303859B 777802002

Pfaffl MW: A new mathematical model for relative quantification in realtime RTPCR.
Nucleic Acids Res 2001, 29:e45. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Livak KJ, Schmittgen TD: Analysis of relative gene expression data using realtime quantitative PCR and the 2^{ΔΔc}_{T }method.
Methods 2001, 25:402408.
DOI 10.1006/meth.2001.1262
PubMed Abstract  Publisher Full Text 
Zimmermann AK, Simon P, Seeburger J, Hoffmann J, Ziemer G, Aebert H, Wendel HP: Cytokine gene expression in monocytes of patients undergoing cardiopulmonary bypass surgery evaluated by realtime PCR.
J Cel Mol Med 2003, 7:146156. Publisher Full Text

Meijerink J, Mandigers C, van de Locht L, Tönnissen E, Goodsaid F, Raemaekers J: A novel method to compensate for different amplification efficiencies between patient DNA samples in quantitative realtime PCR.
J Mol Diagn 2001, 3(2):5561. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Kamphuis W, Schneemann A, van Beek LM, Smit AB, Hoyng PF, Koya E: Prostanoid receptor gene expression profile in human trabecular meshwork: a quantitative realtime PCR approach.
Invest Ophthalmol Vis Sci 2001, 42(13):32093215. PubMed Abstract  Publisher Full Text

Pfaffl MW, Horgan GW, Dempfle L: Relative expression software tool (REST) for groupwise comparison and statistical analysis of relative expression results in realtime PCR.
Nucleic Acids Res 2002, 30:e36. PubMed Abstract  Publisher Full Text  PubMed Central Full Text

Liu W, Saint DA: A new quantitave method of real time reverse transcription polymerase chain reaction assay based on simulation of polymerase chain reaction kinetics.
Anal Biochem 2002, 302:5259.
DOI 10.1006/abio.2001.5530.
PubMed Abstract  Publisher Full Text 
Karlen Y, McNair A, Perseguers S, Mazza Ch, Mermod N: Statistical significance of quantitative PCR.
BMC Bioinformatics 2007, 8:131.
DOI 10.1186/14721058131.
PubMed Abstract  BioMed Central Full Text  PubMed Central Full Text 
Schefe JH, Lehmann KE, Buschmann IR, Unger T, FunkeKaiser H: Quantitave realtimeRTPCR data analysis: current concepts and the novel "gene expression's C_{T }difference" formula.
J Mol Med 2006, 84:901910.
DOI 0900600976.
PubMed Abstract  Publisher Full Text 
Qiu H, Durand K, RabinovitchChable H, Rigaud M, Gazaille V, Clavere P, Sturtz FG: Gene expression of HIF1α and XRCC4 measured in human samples by realtime RTPCR using the sigmoidal curvefitting method.
BioTechniques 2007, 42:355362.
DOI 10.2144/000112331.
PubMed Abstract 
Bukovská A, Čikoš Š, Juhás Š, Il'ková G, Rehák P, Koppel J: Effects of a combination of thyme and oregano essential oils on TNBSinduced colitis in mice.