Specifying unequal replicates in deseq2 with htseq-count input

Question

[kousi31]

I have two controls and four treatments. I calculated gene-counts using htseq.

directory <- "~/d/new_de/gene_counts/"
sampleFiles <- c("P-1.deseqfile",
                 "P-2.deseqfile",
                 "P-3.deseqfile",
                 "P-4.deseqfile",
                 "NP-3.deseqfile",
                 "NP-4.deseqfile"
)

sampleNames <- c("P-1", "P-2", "P-3", "P-4","NP-3", "NP-4")

sampleCondition <- c("Treatment", "Treatment", "Treatment", "Treatment","Control", "Control")
sampleTable <- data.frame(sampleName = sampleNames, fileName = sampleFiles, condition = sampleCondition)
View(sampleTable)
ddsHTSeq <- DESeqDataSetFromHTSeqCount(sampleTable = sampleTable, directory = directory, design = ~ condition)
dds<- DESeq(ddsHTSeq)
res<- results(dds)

When I switch the sampleFiles and sampleNames I again get the same treatment vs control results instead of control vs. treatment. I think I am wrong in specifying the design/treatment conditions properly. Kindly pin-point the error. Thank you!

[benformatics]

I'm not 100% but it could just be taking them in the order you provide them and doing the first obvious comparison.

Try just running the following.

sampleFiles <-  rev(sampleFiles)
sampleNames <-  rev(sampleNames)
sampleCondition <-  rev(sampleCondition)

Before you create your sampleTable data.frame object and see if that fixes it.

Alternatively what is the problem with treatment vs. control ? Just reverse the signs and it should be the same as control vs. treatment?

[kousi31]

Specifying res <- results(dds, contrast=c("condition","B","A")) worked. Thank you.

Answer 1

[benformatics]

If my above comment doesn't work then reading the manual should solve your problem.

http://bioconductor.org/packages/devel/bioc/vignettes/DESeq2/inst/doc/DESeq2.html#factorlevels

Note on factor levels By default, R will choose a reference level for factors based on alphabetical order. Then, if you never tell the DESeq2 functions which level you want to compare against (e.g. which level represents the control group), the comparisons will be based on the alphabetical order of the levels. There are two solutions: you can either explicitly tell results which comparison to make using the contrast argument (this will be shown later), or you can explicitly set the factors levels. In order to see the change of reference levels reflected in the results names, you need to either run DESeq or nbinomWaldTest/nbinomLRT after the re-leveling operation. Setting the factor levels can be done in two ways, either using factor:

dds$condition <- factor(dds$condition, levels = c("untreated","treated"))
  

…or using relevel, just specifying the reference level:

dds$condition <- relevel(dds$condition, ref = "untreated")
  

[swbarnes2]

Yes. Swapping the input order won't change them alphabetically. You have to relevel or specify the contrast. Personally, I always specify the contrast, because it's easier to see in 6 months what you did.