DESeq2 group comparison
1
0
Entering edit mode
rbronste ▴ 60
@rbronste-12189
Last seen 4.4 years ago

Quick question about running some 4-way comparisons, have the following set up:

table(sampleSex, sampleTreatment)

dds<-DESeqDataSet(se, design= ~batch + sex + treatment)

dds$group <- factor(paste0(dds$sex, dds$treatment))

design(dds) <- ~0 + group

dds <- dds[ rowSums(counts(dds)) > 1, ]

dds <- DESeq(dds, betaPrior = FALSE)

res1<-results(dds, format = c("GRanges"), contrast=c(1, -1/3, -1/3, -1/3)

res2 <- subset(res1, padj<.001)

res2_sorted = res2[order(res2$padj), ]

In doing it this way I get results that show the first member of the group contrast (1) as having the smallest number of counts in the matrix I am importing as the basis for DESeq2 - however I would like to have the opposite with the first member of the group having the most counts and the remaining members having the fewest. Not sure how to change the design to reflect that, thanks!

deseq2 deseq differential binding analysis • 1.0k views
ADD COMMENT
0
Entering edit mode
@mikelove
Last seen 3 hours ago
United States

Sorting by (adjusted) pvalue has no preference for the sign of the estimated coefficient.

The fact that you observe more or one sign or the other at the top of the list indicates that these just have smaller pvalue. 

If you want a one sided test, see lfcThreshold in the vignette or in ?results.

ADD COMMENT
0
Entering edit mode

Thanks for the info. I am fine sorting by the (adjusted) pvalue, Im moreso interested in matching the way I do the 4-way analysis here to the following which is what I am doing for straightforward pairwise comparisons (where I don't do the onesided lfcThreshold):

table(sampleTreatment)

dds<-DESeqDataSet(se, design= ~batch + treatment)

dds$group <- factor(paste0(dds$treatment))

dds$treatment<-relevel(dds$treatment, ref="vehicle")

design(dds) <- ~batch + group

dds <- dds[ rowSums(counts(dds)) > "1", ]

dds <- DESeq(dds, betaPrior = FALSE)


Me_pw<-results(dds, format = c("GRanges"), contrast=c("group", "vehicle", "EB"))

Me_pw_padj.001 <- subset(Me_pw, padj<.001)

Me_pw_padj.001_sorted = Me_pw_padj.001[order(Me_pw_padj.001$padj), ]

Me_pw_padj.001_sorted
ADD REPLY
0
Entering edit mode

Is there a remaining question? I think I answered why you see one or the other sign of LFC at the top of the table sorting by pvalue. If you are only interested in one sign of coefficient, use a one sided test. If you want to look at the genes where the LFC is positive or negative, subset the table.

ADD REPLY

Login before adding your answer.

Traffic: 813 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6