First of all, sorry if this is a trivial question, as I am totally new to using DeSeq2. I am attempting to use DESeq2 to look for differential abundance in microbiome data from soil. As such, I have several different factors describing where the samples came from. These are Compartment, with levels "Bulk Soil", "Rhizosphere", and "Endosphere; Rotation, with levels "2-yr", "3-yr" and "4-yr"; and lastly Sampling Time with levels "V1", "V5", and "R2". I am trying to use a grouping variable, which in this case are a combination of the three factors above, to look at contrasts between 2 rotations while keeping everything else the same (e.g. there will be 9 total contrasts, 3 for each of the Compartments comparing the 2-yr vs 4-yr rotation at each sampling time). I put in my model and ran DESeq2 with the following code:
dds = phyloseq_to_deseq2(ps_alpha.trim.fungi, ~Compartment + Rotation + Sampling.Time + Compartment:Rotation + Compartment:Sampling.Time + Rotation:Sampling.Time + Compartment:Rotation:Sampling.Time)
dds$group <-factor(paste0(dds$Compartment, dds$Rotation, dds$Sampling.Time))
design(dds) <- ~ group
dds <- DESeq(dds)
The commands I used to generate the contrasts follow the format below:
bulk2v4_v1 <- results(dds, contrast = c("group", "Bulk.Soil2.yrV1", "Bulk.Soil4.yrV1"), alpha = 0.05)
The problem I am running into is that for one specific comparison (Bulk Soil 2-yr rotation at R2 vs the Bulk Soil 4-yr rotation at R2), I receive an error. The code I run similar code to that above, substituting the other factor levels where necessary as follows:
bulk2v4_R2 <- results(dds, contrast = c("group", "Bulk.Soil2.yrR2", "Bulk.Soil4.yrR2"), alpha = 0.05)
Then I receive the following error:
Error in cleanContrast(object, contrast, expanded = isExpanded, listValues = listValues, : Bulk.Soil2.yrR2 and Bulk.Soil4.yrR2 should be levels of group such that group_Bulk.Soil2.yrR2_vs_Bulk.Soil2.yrR2 and group_Bulk.Soil4.yrR2_vs_Bulk.Soil2.yrR2 are contained in 'resultsNames(object)'
Let me know if there's something obvious I am missing, and thanks in advance.
Kevin,
Thanks for your assistance in this matter. First of all, let me provide you the output of the commands you mentioned in your reply.
If I understand your comment correctly, my default reference level for these comparisons is "Bulk Soil 2-yr R2." My issue is that since no level is really the "control" for the experiment, I don't know what reference level to set. Should I just be picking one arbitrarily, and if so, couldn't I leave it as-is now that I know which is being used?
For your second point regarding fold-change shrinkage, I was unaware that using
contrastdoesn't allow for the default estimator, so thanks for the heads-up! If I wanted to useresults(dds, coef = ...)to run my specific comparisons, would you be willing to provide an example for how I would do this? And is it still possible, without using contrast, to look at all of the comparisons I want with one reference level, or do I need to change it for each comparison usingrelevel?Lastly, would it be possible to use a different estimator for the fold-change shrinkage, and continue to do contrasts as I currently am? Would using
ashrbe appropriate for contrasts, and if so, what is necessary to make this change? Is it worthwhile to switch at this point, or would releveling as above be easier?Sorry for all the questions, and thanks again,
Conard