How basemean/basemeanA/basemeanB and log2foldchange are counted in DESeq2 ?

Question

[IvyChan]

I am new in RNAseq analysis , I was wondering why DESeq2 doesn't show the basemeanA and basemeanB in results just like DESeq ? And how the basemean/basemeanA/basemeanB and log2foldchange are counted ? Thank you.

Answer 1

[Martombo]

As already pointed out, the normal DESeq2 log2FC are shrunken. For this reason, the program doesn't show the two basemeans, since the default log2FC wouldn't match. If you want to know these values, you'll need to get the MLE log2FC, by setting the option addMLE = TRUE when calling the function results(). You can access all intermediate quantities computed by DESeq2 with the function mcols(dds) on the DESeqDataSet object. There you'll find MLE_Intercept column and MLE_contrast_B_vs_A. In case you have a simple single variable design, your first condition basemean will simply be 2^MLE_Intercept and the second condition basemean is 2^(MLE_Intercept + MLE_contrast_B_vs_A). The general basemean, reported in the results, would be the average of the two.

[IvyChan]

Thank you so much . Your answer is really help and I will try it again .

[sutturka]

This is really very helpful information. Can you please suggest how to calculate the log2fc or fc using the values from mcols(dds) output. I calculated first and second basemean as suggested and log2fc by log(basemean2/basemean1). However, it does not match with DeSeq2 output value.

1. MLE_Intercept = 2.470582585
2. MLE_trt_RJn_vs_LJn = 11.96050899
3. basemean1 = 2^2.470582585 = 5.54267565
4. basemean2 = 2^(2.470582585 + 11.96050899) = 22089.78225
5. basemean2/basemean1 = 3985.400489
6. log2(3985.400489) = 11.960

However, in DESeq2 output, it is given as log2FoldChange = 11.77419978

I am not sure why the difference?

Answer 2

[Devon Ryan]

Since DESeq2 shrinks fold-changes I'm not sure how well basemeanB would match what you're expecting. Anyway, "basemean" is essentially the intercept in the GLM, with the caveat that an "extended model matrix" might get used in which case it's more like what you'd get with a ~1 design. The log2foldchanges are then the equivalent coefficients from the GLM. Note that these are shrunken by default, which tends to make them more reliable (so you're not going to calculate them by hand).

[IvyChan]

Thank you for your help.I just curious about what these value mean and how they come from.Thank you , your answer is really help me a lot !