Question

Error bars on DESeq2 visuals

0

Entering edit mode

lauren.czaplicki • 0

@laurenczaplicki-13775

Last seen 6.7 years ago

Hi Bioconductor community!

I'm using DESeq2 to visualize which taxa are significantly enriched in one treatment over a control and realized that the vignette doesn't result in plots with error bars. Is this because what the error bars represent (standard deviations) get wrapped up in DESeq2's p-values or another statistic?

First I get the results from a DESeqDataSet object "DDSfbps1spooled" grouped by "Rep", between "ctrl0" and "treatment" with

res.ctrl.treatment<-results(DDSfbps1spooled, contrast=c("Rep","treatment","ctrl0"), cooksCutoff = FALSE)

res.ctrl.treatment = res.ctrl.treatment[order(res.ctrl.treatment$padj, na.last=NA), ]

alpha = 0.05 #0.01 by default, in my field, env. eng, we define significance at alpha=0.05

sigtab.res.ctrl.treatment = res.ctrl.treatment[(res.ctrl.treatment$padj < alpha), ]

sigtab.res.ctrl.treatment = cbind(as(sigtab.res.ctrl.treatment, "data.frame"), as(tax_table(fbps1spooled)[rownames(sigtab.res.ctrl.treatment), ], "matrix"))

Then I plot it using ggplot2's layers with

treatmenteffect<-ggplot(sigtab.res.ctrl.treatment, aes(x=log2FoldChange, y=Order, fill=Phylum, shape=Kingdom))

treatmenteffect$layers <-treatmenteffect$layers[-1]

treatmenteffect+ scale_x_continuous(limits=c(-50,40))+ ylab("")+xlab("") +

  theme(legend.position="none",text=element_text(size=16, family="serif"),axis.text.x = element_text(angle = 0, hjust = NULL, vjust=0.5))+

  scale_shape_manual(values=c(21,23))+

  scale_fill_manual(values=c("#41ab5d","#dd3497"))+geom_point(size=5, alpha=0.7)

Thanks!!!

~Lauren

deseq2 metagenomics phyloseq ggplot2 • 1.9k views

ADD COMMENT • link 6.7 years ago lauren.czaplicki • 0

0

Entering edit mode

Thanks Mike! A few follow up questions:

Say I were to put the lfcSE bars on the chart without trying to build confidence intervals, what would my N be? In other words, would the lfcSE reflect the error between the starting number of communities I used for the comparison (N=6 if triplicate reactors), or the number of comparisons made (N=9 if all against all)?
I wouldn't be able to compare two different DESeq plots to each other without this CI correction, correct?

ADD REPLY • link 6.7 years ago lauren.czaplicki • 0

0

Entering edit mode

The SE on the LFC are calculated using all the samples in the DESeqDataSet. It doesn't have an 'n' like you would put if you were estimating the mean for a population, and then saying that you observed 'n' from the population.

I guess the multiple test correction may be getting too complicated. I would just mention if you draw the SE that it is the estimated standard error for the log fold change from the model (DESeq2).

ADD REPLY • link 6.7 years ago Michael Love 41k

score 1 · Answer 1 · 2017-08-17

hi Lauren,

For each LFC, we do return a column in the results table called lfcSE which is the standard error on the estimated fold change. Be careful with the standard error and using that to build confidence intervals. Those CI will not be corrected for multiple testing, so you would need to note that if you show it in a figure (it's equivalent to reporting an uncorrected p-value). There are methods for doing something like correction of CI for multiple testing, but we don't provide this functionality in DESeq2.