I have a small toy example. I have a treated and untreated samples from RNA-seq experiment. I calculated fold changes treated/untreated (Assuming, I have only one treated and one untreated sample and I want to look at the genes that are upregulated in the treatment, so fold changes >1). The values are not normal distributed, they are more like from extreme distribution. And now I want to calculate p-values or some kind of significance-value. Zero hypothesis is that there is no different in treated vs untreated, so the fold change will be 1.
Below is a of the distribution of a fold change of treated_sample_1 / untreated_sample_2 (vertical line is at 1). So, not being strict, I would tell that all genes above >1 are up-regulated in treatment vs non-treatment. But I need some kind of a significance value how far are they from 1. Is it possible to calculate a p-value out of it?
I could standardize it, so that the mean=0 and sd=1 but it still would not be normally distributed (plot below). From z-scores I could obtain p-values like 2*pnorm(-abs(z-scores)) but I do not think it is appropriate.
1 answer
You can't compute p-values from fold changes alone.
The edgeR and DESeq2 packages both give advice on what to do with RNA-seq counts in the N=1 situation. Both provide shrunk log fold changes and notional p-values, with appropriate caveats.
Log in to answer this question.
-
Transform log2 fold changes in z-scores
written by gg •
-
How to compare different groups p value and fold change
written by F.eghbalpour.amt •
-
Change in variance Treated-Vehicle
written by bioinfonerd •
-
Delta z score between control and treatment
written by bioinfonerd •
-
considering q- value for up/down regulated genes
written by ahmadvandra •
-
considering q- value for up/down regulated genes
written by ahmadvandra •
-
p-value distribution of edger or deseq2
written by tonja.r
-
high fold changes due to the low counts
written by tonja.r
-
GFOLD posterior distribution
written by tonja.r
-
Error Bars Following Normalisation Of Real Time Pcr Data
written by Dolores Hamilton •
The values are not normally distributed because fold-change is a fraction. Consider the division of two normal bell curves: a lot of your divisors will be nearly zero; and thus very unstable. You should log-transform the values to reign in the extra-large scores.
After log2-transformation I would get a normal distribution, apply standardization (z-scores). Would it be applicable then to calculate p-values as 2*pnorm(-abs(z-scores)?
It's not exactly correct yet, a p-value is from a hypothesis test you never specified, regarding variability you never measured. The values still won't be exactly a bell-curve, but you'll get pretty close for most uses.
The hypothesis is that some genes are expressed in the treated and not in the non-treated. Also, lets say I have two treated samples and two untreated. Fold-changes are calculated only for treated1/untreated1,treated2/untreated2 and log2 transformed. Is there anyway to say that a log2(fold_change) of 4 is more significant than log2(fold_change) of 2? Like, is it possible to assign a significance value?
The program you're looking for is called GFOLD. Note that one should always be hesitant in believe results from N of 1 experiments.
Could you please explain in a bit more details, why approach with a p-value would not work if I have my hypothesis that non-differential genes would have 0 (in log2 scale) as a fold change. So, I would just rank the genes with the p-value?
No it's not possible. Genes with large fold changes but small counts are not necessarily more significant than those with smaller fold changes but larger counts.
If you have two treated samples and two untreated, why not use edgeR or DESeq2 in the usual way (which I know you are familiar with from your other posts)? There's no shortcut.