Hi,
I have RT-qPCR data that I am analysing currently, but am unsure of how to correctly express the statistical analysis of it.
The experiments (two independent replicates) contain the following;
Control sample, test sample
Endogenous control gene, 6x target genes
I want to express the data as fold change, but need to add error bars and statistically analyse for significance. I have thus far calculated the delta delta CT value for each set of samples for each target gene relative to the endogenous control, and have calculated fold change (=2^-ddCT) of the test samples relative to the control samples.
Do I use the Std Dev of the ddCT values as the error bars, and how do I analyse for significance? (my thought was using test versus control samples with independent samples T test, but with what values?)
Thanks in advance
Daniel:
My understanding is that you should not present qPCR data using CT or dCT values because of the scale, thus for presentation of the data - I express the results as fold change, with error bars to show the SD.
NO! - you should present dCT values or ddCT values! If you exponentiate these values (to get "relative expression" [2^dCT] or "fold.change" [2^ddCT]), then you should use a diagram with a log-Y axis.
Regarding your question:
DO NOT DO STATISTICS ON "FOLD VALUES".*
DO ALL STATISTICS ON (D)DCT VALUES.
You can average the dCT values and calculate SDs and CIs as you wish.If you want (despite all warnigs) go for the relative expressions (2^dCT), then you will need the average and the upper and lower limits of the CI for the dCT values, let's call them M, U, and L, and these have to be exponentiates to get 2^M, 2^U and 2^L as the respective statistics on the linear scale that you can subsequently use to plot the fold-changes with respective error-bars.
If you want to show fold-changes (2^ddCT) then you have to get the mean and the CI for the ddCT value and transform them to the linear scale as for the dCT values.
* means "relative expression values" like 2^dCT or "fold-changes" like 2^ddCT.
Consult Jochen Wilhelm on this forum and all Dr. Jo Vandesompele and Dr. Michael Pfaffl publications on this topic. You will find that the appropriate/responsible algorithms for this have been tucked inside of Biogazelle and REST software programs. The questions you ask here are not simple ones. First off, you should use efficiency correction and not assume 100% efficiency for all reactions (as the 2-ddCt method does). Also, Applied Biosystems Bulletin #2 is a good thing to read before you embark on this journey. This is near the start of the problem/question/answer. Ignoring the error associated with estimations of amplification efficiencies is the common caveat many fall into. ABI literature (original machine manuals) instruct when and where to pay attention to error and where to ignore it. This is still a major point of concern in all of qPCR. Hopefully you are dealing with much greater than 2-fold changes to begin with, and hopefully your reference gene Cq values are indeed unaffected by treatment conditions.
Given your design is good, the expression of your reference gene is stable, and your PCRs run according to the required assumption of the ddCt method, then (!! see Jack's post !!) the solution to your problem is quite simple:
You have two groups: control and test.
For each group and each target gene you get several dct values (dct as the difference of ct[ref.gene] and ct[target gene] for each sample).
The null hypothesis of equal relative expression translates to the null hypothesis of equal dct values, and this hypothesis can be tested with the two-sample t-test (using the dct values).
An extra note: When you have several different target genes, I'd suggest to use the pooled standard error, that's more robust. Practically an elegant way to get the results is from a linear model for the dct values on the interaction between "target gene" and "group".
Thank you both for your input.
My understanding is that you should not present qPCR data using CT or dCT values because of the scale, thus for presentation of the data - I express the results as fold change, with error bars to show the SD.
Would it be correct to combine all control samples into an average and use the multiple test samples as individual samples to calculate individual fold changes, and then calculate the SD of the fold change for these to use as error bars, and if so, how are these statistically analysed?
Daniel:
My understanding is that you should not present qPCR data using CT or dCT values because of the scale, thus for presentation of the data - I express the results as fold change, with error bars to show the SD.
NO! - you should present dCT values or ddCT values! If you exponentiate these values (to get "relative expression" [2^dCT] or "fold.change" [2^ddCT]), then you should use a diagram with a log-Y axis.
Regarding your question:
DO NOT DO STATISTICS ON "FOLD VALUES".*
DO ALL STATISTICS ON (D)DCT VALUES.
You can average the dCT values and calculate SDs and CIs as you wish.If you want (despite all warnigs) go for the relative expressions (2^dCT), then you will need the average and the upper and lower limits of the CI for the dCT values, let's call them M, U, and L, and these have to be exponentiates to get 2^M, 2^U and 2^L as the respective statistics on the linear scale that you can subsequently use to plot the fold-changes with respective error-bars.
If you want to show fold-changes (2^ddCT) then you have to get the mean and the CI for the ddCT value and transform them to the linear scale as for the dCT values.
* means "relative expression values" like 2^dCT or "fold-changes" like 2^ddCT.
You have two groups: control and test.
For each group and each target gene you get several dct values (dct as the difference of ct[ref.gene] and ct[target gene] for each sample).
@Jochen Wilhelm I strongly suggest you to do something like a Nature methods because I find RTqPCR data analysis done in a wrong way practically in each lab I have visited.
Actually, the p-value (if calculated correctly) refers to the log fold-change. The tested hypothesis is that the expected log fold-change is zero (according with the assumption that individual log fold-change values scatter symmetrically around the expected value -> "normal distribution").
The dCt values are measures that are proportional to a log expression, and a t-test using two groups of dCt values is a test about the hypothesis that the expected difference of the dCt values between the groups is zero. The expected difference of dCt-values is an expected difference of log expressions, what is just the expected log ratio of the expressions (i.e., then log fold-change).
p-values + SD + sample size allow to recompute means, here log-fold change (without the sign, however), but that's not a very practical way to express results. p-values alone are of very little interest. The best way would be like for any comparison of means between two samples: « mean DCT changed from 24±2 to 28±4 between A and B samples, and this was statistically significant (+4, 95% CI [+1 ; +7], p = 0.021 [T-test for 2 independent samples with equal variances]) » or variations around this. If A and B sample sizes were given somewhere, typically in the methods, and if there was no problems of missing or below the LOQ values, all information needed is available (or else, the full dataset is needed...)
@Emmanuel,
although you are correct in principal, the effect measure (and CI) you provide is atypical for quantitative real-time PCRs. The interesting measure is simply the (expected) difference in DCT values (in your example: 24 - 28 = -4), whereas you give the relative difference (4.95%). The DCT values are proportional to log-concentrations, and their difference is a log-ratio of concentrations. The ratio of DCT values is something strange that is hard to interpret, practically.
EDIT: This comment is wrong. Emmanuel correctly wrote that the ddCt (B-A) is 4, with a 95%-CI ranging from 1 to 7. See Emmanuels response below.
@Jochen: I intended to give « +4 » and « 95% CI : » were CI was for confidence interval, with the comma as a separator — but I agree it is easily misread as « +4.95 % », especially for readers whose native language use the comma as a decimal separator (like me!). Would have been less error-prone if CT would have a unit... Or if I used a semi colon, « +4 ; 95 % CI ... »
Hello!
As far as I know, it is not false to show the data in log fold changes, even if you were analysing the dCt values. It can be more easily interpreted, and there is a well defined mathematical operation between dCt and log fold change. If you use any transformation on your data to analyse it, you can still use the original data set to show. Does it really matter if you announce the winner of a 100 m running competition via her time spent, or average velocity = 100/t ?
I am still confused on how to get a p value for the fold change? Can we give fold change along with p value calculated from delta ct compsrison between the cases and groups or therd is another way to calculate the p value for fold change? Foe exam if I have fold change of 0.5 in cases and the independent sample t test calculates the p value of 0.00 for the same set of delta ct of cases and controls , can I present this as “ The vegf gene expression was significantly higher in cases (0.5 fold, p value 0.00).
If you calculate the p value on dCt values, that is good. Then you transform dCt into fold change and draw a figure. dCt to log fold is a simple mathematical operation, just like any kind of data transformation. Make sure you state statistics calculated for dCt values, no matter what form of the data you draw on the figure. The thing you are interested in is the difference between the samples, not in log fold changes. You just give log fold change because it is easier to understand, and statistics is better on dCt.
I do not recommend that, because you evaluate the control value with an uncertainty. If you take it for one, you lose the information on uncertainty. Leave the log fold changes, dCt values are better for statistical analysis, regardless what you draw on the figures. Just state what you really did with the data. If somebody asks why you analysed the dCt, tell them because log fold changes are from skewed distribution, which you would have to transform into something simiar like dCt anyway, so why bother.
The t-test is not the right tool to test fold-changes, because fold-changes are not normally distributed. However, dCt values are normally distributed. So the easy and correct way is to use a t-test to demonstrate a mean difference between dCt values. Such a mean difference of dCt values is an estimate of the log fold-change. Thus, a t-test on the dCt values is testing the log fold-change.
Although to my opinion it is absolutely sufficient to report the mean difference of the dCt values (what is the ddCt value), preferrably together with its 95% confidence interval, you can transform this (log fold-change) by exponentiation to express it as fold-change. That the p-value tests the hypothesis using the dCt values should be described in the methods section, so it is perfectly clear to the reader if you write, for instance, that "gene X was induced 3.4-fold (1.2; 7.1), p = 0.023" (the values in the parentheses are the limits [in "folds"] of the 95% confidence interval; note that the estimate [3.4] is not in the center of this interval).
Any half-way good stats software will give you the confidence interval.
I think this will be a useful thread for many investigators, with particular thanks to @Jochem!
Rather then create a new thread, it might be useful to have the answer here for continuity - I wanted to ask what format is recommend for plotting, e.g., heatmaps or ordination plots?
Assuming one has several primary cell lines, under different conditions, and testing several genes, a data matrix can be generated with the genes in row 1 and the samples in column 1. So would the data matrix be filled with the values from "dCT" (HK gene - target gene), or "ddCT" (dCT test - dCT control), or other such as "fold change"?
An issue with ddCT is there is no control sample per se.
Input appreciated!
Chris
@Christopher, I think the best here is to show the dCt values in a heatmap. However, it's difficult to normalize the expression between different cell lines. It is crucial to have a set of reference genes that is expressed at similar levels between all the cell lines - and that may be hard to find.
What if you have more than one treatments? For example, I have one control condition and three treatments. Could you recommend the best statistical test to use in this situation?
You would always compare dCt values between two groups using a t-test. The difference is in how you calculate the t-statistic. If you have data from more than two groups, you can use a pooled variance estimate to get more robust t-statistic. This is done in multivariable models.
There is also the option to test if treating in general is better than not treating, or if a subgroup of treatments is better (than the control, or than another subgroup). This requires approriate contrasting within multivariable models.
If you screen among treatments, you will have to take care of the multiplicity problem and you should consider to adjust the p-values for multiple testing. This is done, for inststance, by Tukey's HSD tests (which are modified t-tests, using the pooled variance and controlling for multiple testing for all-pariwise comparisons). Other options are Bonferroni/Holm, Dunnett's MCP (multiple-to-one) and a myriad of other methods, each with its particular advantages and disadvantages in paricular situations. If you screen among a very large number of treatments, it might not make sense to interpret individual p-values at all and use them only to rank the treatments by their observed (statistical) signal-to-noise ratio.
@ Melissa: if you want to compare all treatments "in turn" to the control, but not the treatments between them, the best choice would be Dunnett's method, assuming usual conditions are met (that is, normality, independance and equality of variances).
If you have several genes to test, do not forgot to also correct for multiple testing from that point of view. Or to replace the dCt method by the method in
https://academic.oup.com/bioinformatics/advance-article/doi/10.1093/bioinformatics/bty629/5053322
that automatically handles the problem of the multiplicity correction for several genes (not for three treatments to one control) and does not rely on the invariance of "reference" genes expression.
I have RT-qPCR data that I am about to analyze, but am unsure of how to correctly express the statistical analysis of it.
The experiments (two replicates) of a cross over study 4 weeks treatment-4 weeks washout-4 weeks treatment (I don't know who had the drug vs placebo at this time) contain the following:
reference gene, the genes I am testing ex. GAPDH, with two time points after 4 weeks of "treatment (placebo or drug)" and then after 12 weeks of treatment.
here is what an example of my numbers look like
Patient A week 4 A week12
ACTA1 13.85 14.16
Reference Gene 14.21 14.14
AVG 14.03 14.15
Patient A Week 4 AWeek 12
GAPDH 17.55 18.02
17.71 18.06
AVG 17.63 18.04
ACTA1 14.03 14.15
minusACTA1 3.60 3.89
DD CT 0.28
Fold Change 2 -DD CT 0.82
so do I analyze 3.6 and 3.89 and do t - test there and also get the rest of the parameters or how do you recommend I get this done? Also could someone please add a table of how this is reported in an article without having to report all of these values.
I don't think my numbers above come out as nice as when I typed them in therefore I have a word document with just a sample of my data.
Thank you in advance.
If you have a cross-over design, you should analyze your data respecting this design, that is using a model with a random-effect patient factor, a period factor, a sequence factor and a treatment factor. Then interpret the adequate coefficient, or sum of square (depending on your usage... and model coding), in corresponding models. That's more complicate than just T-tests, especially if it happens that there is a period or a sequence effect... So the DDCt method is definitely inadequate, it just can't handle such a design.
You may do this on the DCt, or with alternative methods (see the previous post). But probably should also consult someone that has enough statistical background to help you in analyzing a crossover...
Hello Everyone!
I'm no good with statistics but I'm trying to analyse my qPCR data.
What I get from this old discussion (mostly Jochen Wilhelm ) is that you should use FoldChange for graphic representation and use dCt for statistics (sorry if I got it wrong).
But Jack M Gallup said " First off, you should use efficiency correction and not assume 100% efficiency for all reactions (as the 2-ddCt method does)."
Then your FolfChange take into account your efficiency and your dCt doesn't. Isn't it a problem?
Another question: I draw the standard deviation of the FoldChange of my different biological samples (per condition) on the plot of the mean fold change per condition. Does it seem right to you?
Sorry to start again this topic.
Hi Morgane,
I never recommend using FoldChange anywhere. Also not for graphical representation. It's common, but it's bad.
Jack it right, if the efficiency is way off from being optimal. However, I'd doubt results from an assay based on obviousely suboptimum efficiency anyway (that's a different topic, though). But to your question: you can (you should, if suboptimum efficiency is a concern) efficiency-correct the Ct-values.
The standard deviation of a FoldChange simply is the square-root of the variance of the FoldChange ... apart from this it is hard (to impossible) to understand what it tells you. FoldChanges don't have a symmetric distribution, the variance depends on the mean, and there are other nasty properties of FoldChanges. It makes pretty little sense to add this "information" to a plot, like as error-bars. In the most friendly look, they are simply worthless, but I'd rather say this is misleading (misguiding the interpretation of the data) and therefore it's actually bad science to do things like that. But I repeat myself. If you have read my posts you should know my recommendations regarding how to best plot qPCR data and why I consider this better than barplots with error-bars.
Hi Morgane,
For the first question, you should consider that in a first step you manage to treat your raw data in a way that gives correct « Ct » (or Cq or...) or « Fold change », but only one of the two. And in the second step, when you consider this value as the best possible representation of the molar (or mass) fraction of RNA in your sample, you define CT = log( fold change ) (or fold change = exp( CT ) — the basis of the logarithm/exponential has really no importance, as far as you use the same for all conditions), so that they are absolutely equivalent, and you analyse the CTs and represent the one you prefer. So either both of them or none of them correct for different efficiency.
For the second question: FoldChange being exponential of the CT and by constraint positive, standard deviation quite badly represent the dispersion of your data because the distribution is asymetric. A better way would be to compute a pseudo-prediction interval on the CTs (mean CT ± standard deviation of the CT), and then exponentiate this interval and draw it on the FoldChange scale. But beside this drawbacks, it's not wrong per se (but is it suited to your problematic strongly depends on your problematic).
Jochen Wilhelm sorry I misread a comment were you quoted someone else about plotting FoldChange.
I never thought about plotting anything else than FoldChange because I was taught this way (Just did a rapid check _ in the small cephalopods/molluscs qRTPCR world _ it seems using using Log2 of the Fold Change is the rule).
Emmanuel Curis Ok I get your point for my 1st question, , so I guess as I already corrected my data when calculating fold change, I can do my statistics on the Log of the FoldChange? (Or I find a way to correct my Cts).
Thanks a lot
Morgane Bonadè: yes, if the fold changes are the final results of your first step, then do all the analysis on their logarithm (and call this logs Cq for instance, it seems it's the suggested notation now to distinguish these from Ct).
Don't forget that this analysis strongly depends on the experimental design and the question you ask, DDCt (DDCq) beeing usable in a very very few contexts on this aspect, and that most often you cannot conclude on single genes but only on ratios of genes, including ratio with reference genes.
I saw most articles use Student T test on the ddCt values or logFoldChange, is it correct? My standard deviation are quite close to each other and I guess they should follow a normal law if the dCt do follow one... but I'm not 100% sure:
- I compare three different situation between them (t; t+2; t+3). I'm just afraid the first comparison (control to t) is different as all the control are set to 0 than the other comparison (all the ddCt values take control as a 0).
On the other hand, I tried doing the T test on the dCt and I reach significativity way more easily which makes me feel quite uncomfortable. (I'm not a statistician, it might be perfectly normal)
Hi @Morgane Bonadè
Student's T-test is valid if you are comparing 2 groups. In your case, you have more than 2 - you have three groups which are (i) group t (ii) group t+2 and (iii) group t+3. Analysis of Variance (ANOVA) is used when comparisons are to be made between "more than" 2 groups which is the case in the current discussion. I would not want to bore you with too many technical jargon on statistics. Nonetheless, please do not be tempted to conduct multiple T-tests in the form of;
GROUP 1 vs GROUP 2
GROUP 1 vs GROUP 3
GROUP 2 vs GROUP 3
it would be a very bad idea.
Kindly spare a few minutes to look though this clip to get a feeling of how ANOVA works >>> https://www.youtube.com/watch?v=-yQb_ZJnFXw
On a side node, subjecting a dataset containing only 2 groups to ANOVA should return a T-test based results to you, but T-test will not return an ANOVA-based results as long as you are comparing more than 2 groups.
This link might also might be useful to you as well >>> http://blog.mcbryan.co.uk/2013/06/qpcr-normalisation.html?m=1
If you check back to see what Emmanuel Curis and Jochen Wilhelm said one more time, about estimating the spread of your qPCR data, you have to take special care in computing your standard deviation. It may be the reason why you have observed what you mentioned.
ddCt is already a comparison. It is a log fold-change, and "change" means from one condition (group) to another. There is nothing like a fold-change within a group. If the groups are defined by time-points (A:t, B:t+2, C:t+3), then there are thrre possible comparisons:
B vs. A, C vs A, and C vs. B.
For each of them you get a (single, mean) ddCt. This ddCt can be positive (indicating a higer expression in the first group as compared to the second group) or negative (indicating a lower expression in the first group as compared to the second group). If you have the mean ddCt, its standard error (SE), and the number of (dCt-)values used to calculate the ddCt (n), then you can calculate the t-value as ddCt/SE, which is t-distributed with n-2 degrees of freedom. You can thus use the t-distribution to find the p-value for your ddCt-value under the hypotheses that the expected value of the ddCt is 0. However, the difficult part ist to calculate the SE of the ddCt value (google "pooled standard error").
You may not be familiar with this way to calculate a p-value. Usually, it's simpler and more straight-foreward to compare two (independent) samples of dCt values with the usual two-sample t-test (in each group you have several dCt values). This is also about the mean difference (in dCt-values) between the groups, what is just the ddCt value (ha!). So this test is equivalent to what I explained for the ddCt values, and the tricky calculation of the SE of the ddCt is done automatically and under the hood.
There is no need to be concerned with the assumption of normal distribution. Experience shows that the distribution of dCt values in fact always approximate the normal distribution model quite well. There is also a paper demonstrating that the expression levels must follow a log-normal distribution, and sinde dCt values are proportional to log expression, their distribution must be normal. So there are theoretical and empirical arguments to believe that this assumption is apt.
@ Morgane: to complement what Jochen and Damilola wrote:
1) standard deviation never follows a normal law. I think you wanted to say that your data/your residuals/your errors/your variations around the mean followed a normal law.
2) Having t, t+2, t+3 suggests potential paired (longitudinal) data, where the same subject is measure at different times. If indeed so, that changes completely the analysis, that must take this pairing into account.
3) if you have all 0 for one group/gene, there is something wrong. That is not possible. So the analysis method is simply not suited. It may well explain you unexpected, all significant tests.
You should give a precise description of your experimental design and of your question to have sensible advices on how to analyze your data. A first hint that is always applicable is anyway :
1) first think what kind of analysis would be suited if you could analyse data for a given gene alone (or, if you prefer, if your data were, let's say, patients weights instead of CT / FC)
2) then apply this analysis method to the (log of) ratio of gene quantities [or to differences of CT between two genes, in the CT=log scale; this one of the delta in the delta-delta-Ct method]
You may either consider ratios to a predefined gene (control gene; the most usual, typically the DDCt method) or all possible ratios (paper in Bioinformatics, 2018) ; in all cases, beware of multiple comparisons issues and unverifiable additional hypothesis (like « reference gene quantities do not change between conditions ») in the interpretation of the results.
In order to be clearer Emmanuel Curis :
I'm looking at 4 different biological stages of an embryo. So I chose stage 1 as control. And I want to see two things for my several genes:
1. If the expression is modified from stage to stage? (null hypothesis stage n and n+1 have the same level of expression)
2. Does all/some of my gene have the same pattern of expression?
For the second question I know I need an ANOVA (Damilola Akinyemi)
For the first one I just had trouble to have some stats about the compared expression from stage 1 to 2 as my stage 1 is also used as the control. But I guess whether I just used dCts (and find a way to correct them with efficiency or just assumed my efficiency aren't so different).
Jochen Wilhelm I understand your point about SE being easier to calculate for dCt values. For now I've just been using the "Student T test function" on Excell on my different ddCt values so I guess it was just calculating the SE wrong. I'm invastigating the "pooled standard error".
In case you were wondering most my efficiencies are ranging from 1.91 to 2.01. Sadly one is 1.86 and another 2.14 but I'm currently trying to get better efficiencies for these two genes by trying another set of primers.
@ Morgane: I'm not sure to understand why different efficiencies matter, once you know them: just use them to have better estimates of the initial quantities of RNA in your middle; then use log of these quantities as your "new CT" and do all analysis with them.
The fact that stage 1 is control does not matter neither; you're stuck with the delta-delta Ct method which is irrelevant here.
The right question is do you expect a monotonic, eventually linear, increase or decrease of your expressions with stage ? If so, the natural model would be to make a linear regression of the RNA [log] quantities on the stage number / time of development. If not, the most general model would be to do an ANOVA of RNA [log] quantities on the stage category.
Now, since raw interpretation of CT/RNA quantities in qRT-PCR is irrelevant because of the way they are obtained [efficiency changes being one reason, but not any more after correction for that], just do this analysis on the [log] ratio of your gene and another gene.
Typical ratio would be of your gene and one reference, control, housekeeping gene (or the average of several such genes).
For the second question, ANOVA does not help at all; ANOVA is what you need for your first question. You may perform a multivariate ANOVA, and suited post-hoc tests. But basically, use of a control gene as a denominator prevents comparing genes between them, you can just say if your genes behave in a similar pattern than the control gene or not.
For your second question, I would suggest the use of all pairwise ratios and the method described in the following paper
https://academic.oup.com/bioinformatics/article-abstract/35/2/258/5053322
that is just meant for that [self citation…]. This method should also limit the effects of different efficiencies between genes (as far as a given gene as the same efficiency at all stages of development — that is, efficiency is gene dependent and not experiment dependent).
Emmanuel Curis I guess estimating "real" Ct is trivial but I've never learnt. Should I just do a cross product?
Ct[real]=Ct*efficiency/2
For my first question what I was doing for now was a T-Test compairing stage 1 to stage 2 then one compairing stage 2 to stage 3 and a last one comparaing stage 3 to stage 4. Based on the ddCt values I get when I compare all the both stage to stage 23. (I don't if it's clear: I calculate my ddCts by compairing everyone to stage 23, then I run the Student T test on this value for 2 following stages)
I guess an ANOVA might be better but it wouldn't give significativity values that I could put on the graph if I plot my ddCts (or my dCts if I do the stats on them)?
I hope it was clear enough.
If you mean the efficiency-corrected Ct-value, this is calculated differently.
From some starting amount N0, the amount of PCR product increases exponentially with efficiency a (1 ≤ a ≤ 2). After c amplification cyles, the amount is
N(c) = N0*ac
The measured fluorescence signal F is proportional to this amount. With p being the proportionality factor, the signal is
F(c) = p*N(c) = p*N0*ac
The Ct value is that cycle, after which the signal crosses a given threshold Ft:
Ft = p*N0*aCt
To get the efficiency-corrected Ct value, we ask what the Ct would have been, had a been 2.0. I denote this (corrected) Ct as Ct':
Ft = p*N0*2Ct'
This refers to the same threhold, so both equations are equal:
p*N0*2Ct' = p*N0*aCt
The factors p and N0 cancel out, so we have 2Ct' = aCt. This can be solved for Ct'. Take the log2 on both sides we get
log2 2Ct' = log2 aCt
Ct' = Ct * log2 a
That's it.
----------------------------------------------------------------------------
Example: a=1.7, Ct = 28: Ctcorr=28*log2(1.7) = 21.4 (6.6 cycles earlier)
Example: a=1.7, Ct = 12: Ctcorr=12*log2(1.7) = 9.8 (2.2 cycles earlier)
----------------------------------------------------------------------------
The correction can change the Ct quite a bit. If the efficiency is not determined precisely, the corrected Ct values can have a huge error.
What if efficiency a is determined with an error of ±h:
The following equation gives the number of cycles the corrected Ct differs between assuming a+h and a-h:
Ct*log2(a-h) - Ct*log2(a+h) = Ct*(log2(a-h) - log2(a+h))
This can be interpreted as the possible bias introduced by the correction. It can also be interpreted as the minimum uncertainty in the corrected Ct value.
----------------------------------------------------------------------------
Example: a=1.7, Ct = 28, h = 0.1: 28*(log2(1.8) - log2(1.6)) = 4.8 cycles!
Example: a=1.7, Ct = 28, h = 0.05: 28*(log2(1.75) - log2(1.65)) = 2.4 cycles, still!
Example: a=1.7, Ct = 12, h = 0.1: 12*(log2(1.8) - log2(1.6)) = 2.0 cycles!
----------------------------------------------------------------------------
It's thus extremely important that the efficiency is known quite precisely. The best way to ensure this is to have a PCR that runs with (almost) perfect efficiency, so that a must be very close to 2.0. If the -estimated!- value is clearly below that, you take a huge risk to get grossly wrong results, no matter whether you correct or not.
So am I correct in assuming, Jochen Wilhelm that you would recommend using an efficiency corrected Cq, [D]DCq method, rather than something like Pfaffl method? Because I have always used Pfaffl, but most of the examples used in this thread are impossible using Pfaffl. Indeed I have always assumed that -
"the difference of ct[ref.gene] and ct[target gene] for each sample"
- is a meaningless calculation because these primers have different efficiencies and so it is comparing apples and oranges. But if I apply an efficiency correction directly to my Cq values (which is anyway just a different way of correcting just as Pfaffl does), then it makes all my subsequent analysis much easier?
Hi
FGRs in the ICTs of the research @hefri Yodiansyah-
Gread authors,
Hefri Yodiansyah
Alun Parsons, one can efficiency-correct Cq values, what makes it possible to use the simple analysis (dCq as response in general linear models). However, I am not a fan of efficiency-correction. If the efficiency is estimated wrongly or with considerable uncertainty, the results will easily be severely biased or have such a large uncertainty that they are useless (at best). My alarm bell rings when I see that I am not able to get an ideal efficiency during the early cycles of the qPCR. Under the conditions of the qPCR (very short amplicons, vers small sample volumes, rapid temperature changes, great overshoot of enzyme, primers, dNTPs, presence of everything that makes the enzyme feel comfortable (Mg++, BSA, betaine, ...)) it is really strange that not every molecule should be amplified in each sample (unless there is a serious contamination with some inhibitors, what is a problem of sample preparation that can usually be solved prior to qPCR). Mostly, artefacts of the analysis of the amplification profiles (e.g. suboptimal background correction, chosing a suboptimal threshold) and/or of standard curves (not realizing the linear range).
Jochen Wilhelm OK, but I'm not sure I am getting you, it isn't only a failure of primers to amplify every target transcript that can affect efficiency, right? Some of your primers are going to have transient interactions and therefore transient amplification of unintended targets, which means you're always introducing at least some dsDNA in every round that isn't your target DNA, and your SYBR green is going to bind that DNA, affecting your efficiency. Although we would all love to have primer pairs that all share identical annealing temperatures, and which never have transient interactions with unintended targets, the more target genes you have in your qPCR, the greater the variation in primer efficiency you're going to get. I mean the more you lower your annealing temperature, the more you're going to get amplification of unintended targets, but the higher you increase your annealing temperature the more you're going to get failure of annealing / amplification for some transcripts. Excuse my scepticism, one can't get an efficiency of 2 with every primer pair in every qPCR, PCR simply never works like that. I see it even is my standard PCRs, if I introduce DMSO I can get a nice sharp band because I'm inhibiting annealing, but I'm going to get less product because my efficiency is compromised, but if I omit DMSO then often I'm going to get unintended amplicons even if I get a lot more of my intended product too. There simply is no way to square this circle, you are always compromising between low-efficiency priming, and non-specific priming, when you set an annealing temperature in a PCR reaction that contains heterogeneous groups of primer pairs on the same plate.
Maybe I'm misunderstanding you, but you seem to be implying that we must always get an efficiency of almost exactly two for every primer pair in every qPCR we run, or else our qPCR has failed, and from my experience that is impossible.
I am saying that primers are not the only possible problem. You must distinguish the apparent efficiency (the value measured in your assays) and the actual efficiency (what you can not know, only infer from the apparent efficiency estimate). There are several things that lead to an apparent efficiency that seems to be suboptimal. Or that seems to be super-optimal: there are reported cases wher this apparent efficiency is 2.2 or so! How would one explain this? In such cases, authors usually claim that the efficiency is 2.0. In cases where the apparent efficiency is 1.8, authors start becoming sceptical. My point is that there are many possible reasons for getting such values (2.2 or 1.8), even if the "true" efficiency is just 2.0. One has to carefully look at this, and I know too many examples where people just made quite a mess wehn determining the efficiency (extreme example: just 5 dilution steps, 2 of them giving Cq values around 35, anyway a line is fitted through all the points, and the slope is taken as a measure of the efficiency, the uncertainty in this estimate is nowhere considered).
I see your point that primers can be the problem and I agree that they are often a problem. To my experience, it is often (not always, though) possible to design primers that do work well. I usually suggest ordering a set of 2-3 forward and reverse primers, all close together (maybe simply shifted by 1-5 nucleotides) and then test all possible combinations of foreward and reverse primers. This often leads to a pair that has not been identified as "best" by some design-algorithm, but that works batter than the original pair that has been suggested.
Unspecific amplification is only a problem when it really competes for the resources before the Cq of the target gene is reached. It does not matter if unspecific product accumulates without competition or if if accumulates in later cycles. I have seen cases where "primer-dimers" heavily accumulated in late cycles, especially in reactions with a low initial concentration of the sepcific template. However, the amplification of the specific template in earlier cycles was not compromised.
But I agree, if you have only such bad primers with heavy non-specific and inefficient priming, you have a problem. If you can't find better primers, you have to cope with this, and efficiency-correction is a possible solution. If you are careful that the efficiency you determine is not biased, it should be ok. What you should take into account then is the uncertainty of the efficiency estimate. There is no simple and proper way how this can be done (there are some papers dealing with error-propagation, but that's not always applicable), but you can at least use the "worst-case" scenarios (using the upper and lower bound of the credible range of efficiencies) and check if they lead to similar conclusions.
If you are really concerned with efficiency, you may also use methods that do not rely on a precise value of the efficiency to obtain a pseudo-Ct, that is a value that is the log of a quantity proportional to the initial amount of the target RNA.
Such methods are, for instance, based on modelling the whole fluorescence curve instead of searching for a threshold or a linear part. See for instance the paper by Carr & Moore in PLoS one, may 2012, vol. 7 n° 5, e37640, that assumes a possible change of efficiency from a cycle to another.
I'd like to place a waring here. Using the whole amplification curve is tricky and mostly leads to very wrong results. The data prior to the Ct are by definition not informative, as they represent only noise (that is, in these early cycles, orders of magnitude larger than any specific signal from the amplicons), and at later cycles the shape of the curve also provides limited to no information about the actual amplification process, as this is the result of several processes like re-annealing, depletion of resources, fatique of the enzyme, accumulation of inhibitory side products (e.g. PPi), and the limitation of probes. Using these parts of the amplification profiles for a quantification is quite likee reading tea leaves.
Roche tried a similar method in their LightCycler software. The method is based on the analysis of the shape of the amplification curve, to get rid of selecting a threshold. They took the point of the largest curvature of the curves (this is where the 2nd derivative takes its maximum, hence the name "2nd derivative maximum method"). The results were never really satisfying in practice. To rescue this somehow, they introduced a "plausibility check" of the results to at least identify those results that were grossly wrong (Roche tried to sell this as a "nice feature", however...).
Hello, Using Z ~ 1.96 assumes you perfectly know the standard deviation, which is irealistic in most situations, and especially here for such a small sample.
You can see in the Graph Pad Prism that the correct value to use is 6.14 instead of 1.96, that makes a huge difference...
And the notation ±0.28 to give a confidence interval is definitely a non standard and very bad one... Give [0.53, 3.02] as Graph Pad or any notation that is really an interval... There is so many situations when the interval is not centred on the estimation...
Note that the two values given by GraphPad are the two ends of the CI, NOT 2 CIs as your answer suggest.
I think you'd better follow a basic statistics curse before interpreting GraphPad Prism (or any statistical software) results.
Hello Dr. Daniel Morse
i suggest you read these:-
1- Basic of quantitative real time PCR
Presentation Basic of quantitative real time PCR
2- Analysis of Relative Gene Expression Data Using RealTime Quantitative PCR and the 22DDCT Method
http://www.gene-quantification.net/livak-2001.pdf
3- Analyzing real-time PCR data by the comparative CT method
https://www.nature.com/articles/nprot.2008.73
4- A New Method for Quantitative Real-Time Polymerase Chain Reaction Data Analysis
Article A New Method for Quantitative Real-Time Polymerase Chain Rea...
Hello! Please show me how to do the statistic analysis of my qPCR data.
I would like to compare the relative expression of a GOI between the wild type strain and the mutant strain of a bacteria using 1 reference gene (HKG). Samples were collected at different time points. For example, the result of samples at 6 hours of incubation: 1. Cq of GOI of Wild type strain in triplicate ( 19.92; 19.92; 21.91), and of Mutant strain (26.73; 25.37; 26.3). 2. Cq of Reference gene (HKG) of Wild type in triplicate (19.95; 19.59; 20.72), and of Mutant strain (24.01; 22.29; 23.59) 3. The efficiency of GOI and HKG are 2.01 and 1.9 respectively. Using the Pfaffl method I got the result that gene expression ratio (mutan/ wild type) is 0.18 My question is how to statistic test to conclude that the expression of GOI of the wild type strain is significantly higher than that of mutant strain?
Hi all, and thanks for such a useful discussion.
@Jochen Wilhelm may I ask for a clarification? It's about @Zaima Ali's question: "But if we are applying the t test for delta Ct values, how do we get limits for fold change for 95 percent CI as you gave in the above example , won’t we get limits for delta Ct. or there is some other test for fold change?". I understood it's best to run normality tests and corresponding statistics on the Ct dataset (which will give us CI, significance, etc) and that the delta (delta) Ct method is just an easier/more descriptive representation method. However, how can we describe CI limits for the transformed (log fold-change) dataset? Is it something we derive from CIs of the original dataset or must we submit log-transformed values to additional tests? I believe this applies to any fold change data representation (i.e. not only qPCR data).
Thanks!
Andrea Miccoli,
It's best do do all statistics on the delta-Ct values. Simply consider that these delta-Ct values are your measurements, your read-outs. Each biological sample gives you a delta-Ct value.
It is not really important, what a delta-Ct value is exactly (the exact meaning is difficult to describe). It is important to understand that the delta-Ct value is some monotone transformation of the expression of the target gene. So there is a monotone relationship between delta-Ct and expression, but the precise relationship is not known.* It is further important to understand that the measured variable (the delta Ct value) has an approximate normal distribution. That's all.
You measure something - let's call it Y -, and you know that Y is monotonuously related to the gene expression and that it is approximately normal distributed.
It's simple now to estimate and test the expected difference in Y between groups/treatments. The estimate (the mean difference) will either be positive or negative, implying a difference in the gene expression. What the sign means depend on how you calculated the delta-Ct values. If you calculate them as dCt = Ct[ref] - Ct[goi] , then a positive sign of the difference mean(dCt[treat]) - mean(dCt[ctrl]) = ddCt indicates a higher expression in treat vs. ctrl (and a negative sign indicated a lower expression, respectively). The t-test of the dCt values would tell you if your data are sufficient to interpret its sign (if can you confidently conclude the direction of the regulation). You are actually done here.
I wonder why people want to know fold-changes. There is not more insight you will get from knowing the fold-changes. The relevant information will still be if the gene is induced or repressed, and this question is already answered by looking at the mean difference of the dCt values (ddCt). It would be helpful to know the fold changes if we could interpret these values biologically. But we can't (at least usually). We don't have a clue what expression change is really biologically relevant in the system under study. A 0.9 or 1.1 fold change may or may not be very relevant, nobody knows that. We only have the vague feeling that fold changes below 0.5 or above 2 should usually be relevant for the cell, and we feel that folds below 0.1 or above 10 are quite surely relevant. But that's not different to think that |ddCt| > 1 should usually be relevant and |ddCt| > 3 are quite surely relevant.
As we know (*), and as the dCt values of both groups/treatments were obtained with the same assay, the unknown proportionality factor will cancel out in the ddCt value, so that we can say that the ddCt is a log2 fold-change. The exponentiation of the confidence interval of ddCt (to base 2) will have the same coverage probability over a fold change of 1 than the confidence interval of ddCt has over 0.
One academic problem is that the ddCt as a mean value does not represent the log of the mean fold change.** This means: 2^ddCt is not the mean fold-change. It is rather the geometric mean of the fold-change. If you want to model the mean fold-change, you must really transform the dCt values (2^dCt) and then fit a generalized linear model of the gamma-family with log-link. Nobody does that, as it is not important to estimate the mean fold-change. And even if, this is still an approximative solution. Ct values can not be arbitrarily large or small (they are theoretically bound in (0; 37 [approx.]) and usually in (15; 35). The distribution of dCt values is therefore only approximately normal. One would have to model the Ct values as a function of a Poisson variable - but don't ask me how...
----
*Actually, the relationship is that the delta Ct is proportional to the logarithm of the normalized concentration of the target sequence; the proportionality factor is assay-specific and (typically) unknown.
** E(logb(X)) ≠ logb(E(X)) --> b^E(logb(X)) ≠ E(X)
@ Jochen: assuming a normal distribution of dCT, if one would really want to have the mean fold change (and not the median/geometric mean fold change), an easiest way to have it is to compute it directly using the log-normal distribution properties, typically E(X) = exp(µ + 2 * sigma^2) where µ and sigma are the mean and standard deviation in the log scale (that is, CT scale). That's indeed very academic however.
I would also add that since there is no way to prove that so-called reference genes do not experience expression change, because of the compositional nature of RT-qPCR data (imposed by experimental design), in reality what we test is that the expression ratio of the two genes changes or not between the two conditions, not that the gene of interest expression changes. But that's another problem (but not so academic!).
Yes, thank you for noting that E(X) can be obtained from E(log(X)). The "problem" here is to get a good estimate of sigma². I would think that the error (uncertainty) of exp(m+2s²) (where m and s² are the sample estimates of the log(x) values) is quite large.
Do you know how one would get the confidence interval for this mean?
The point you made that the validity of the results depends on the validity of the constant reference gene expression is surely very important. There are enough examples of badly designed experiments. Chosing good rerefence genes is problematic when very different samples are compared (e.g. different cell types, including tumor cells, or cells in different metabolic states [e.g. growing vs. senescent], or different cell populations [e.g. healthy tissue vs. tissue massively infiltrated by lymphocytes]). This is done very often, and only very rarely the authors carefully think about a valid reference.
A good practice is to use a whole buch of (putative) reference genes, hoping that differences in their expression between groups will cancel out (at least partly).
Sometimes one can standardize the protocols (using a known number of cells, doing parallel RTs) so that it may be ok to not using any reference gene. If RT efficiency is a concern, one can use spike-in RNA controls as reference.
Well, the estimate of sigma² would be exactly the one used when making T-tests in the CT-scale, so the "best" possible for given the sample size.
Obtaining a confidence interval is more problematic, however. « Delta-method » or « bootstrap » would be the tentative answers, however I'm not sure that usual sample sizes will make these methods very efficient (but that would be the same problem for a generalized linear model, whose confidence intervals rely on asymptoric results...). Exact derivation assuming a Gaussian distribution and the fact that under this assumption, sigma² and µ estimators are independent should be possible, since their distribution is known, hence the distribution of their sum could be obtained by convolution. I'm not sure that would add any value to the use of the median/geometric mean whose confidence interval is much straightforward to obtain, as you explained.
For the reference gene problem, I'm not convinced that averaging a lot of genes is a good method, since it increases noise and remains sensitive to a single gene having large variations. Of course, reference gene selection methods try to avoid this, but it's not clear how their use control false detection rate and so on. I strongly prefer approaches assuming the limitation of compositional data, not "strongly" relying on reference gene assumptions.
Spiking in RNA controls is an alternative. Another alternative would be a self-citation...
Jochen Wilhelm You have mentioned in your previous answer the following "What the sign means depend on how you calculated the delta-Ct values. If you calculate them as dCt = Ct[ref] - Ct[goi] , then a positive sign of the difference mean(dCt[treat]) - mean(dCt[ctrl]) = ddCt indicates a higher expression in treat vs. ctrl (and a negative sign indicated a lower expression, respectively)".
I am confused. Isn't the opposite (positive sign indicating lower expression and negative sign indicating higher expression in treat vs ctrl)?
dCt[treat] - dCt[ctrl] = ( Ct[ref,treat] - Ct[goi,treat] ) - ( Ct[ref,ctrl] - Ct[goi,ctrl] ) = (assuming Ct[ref, treat] = Ct[ref, ctrl] as would be expected for a nice reference gene) Ct[goi,ctrl] - Ct[goi,treat]. So having ddCt > 0 means Ct[goi] is lower in treated group, hence that less cycles are needed to detect the goi, hence that it has a higher expression.
Jochen is right, I guess...