Unable to run DESeq command in DESeq2 without transforming data
1
0
Entering edit mode
watsons2 • 0
@watsons2-22383
Last seen 4.4 years ago

I've used DESeq2 previously on my 16S microbiota data, where I ran the phyloseqtodeseq function then went straight on to the deseq function, e.g.

MERGEDnonsinglet.age.dds <- phyloseqtodeseq2(MERGEDnonsinglet.age, ~ Ageclass) MERGEDnonsinglet.age.dds2 = DESeq(MERGED_nonsinglet.age.dds1, test="Wald", fitType="parametric")

...where 'MERGED_nonsinglet.age' is my phyloseq object with singletons removed (but not rarefied).

However, with recent updates to the package (V1.22.2 I think) I can no longer do this as I get the following error:

estimating size factors Error in estimateSizeFactorsForMatrix(counts(object), locfunc = locfunc, : every gene contains at least one zero, cannot compute log geometric means

I've managed to get the script running but only if I include the geomeans function:

MERGEDnonsinglet.age.dds <- phyloseqtodeseq2(MERGEDnonsinglet.age, ~ Age_class)

gmmean = function(x, na.rm=FALSE){ exp(sum(log(x[x > 0]), na.rm=na.rm) / length(x))} geoMeans = apply(counts(MERGEDnonsinglet.age.dds), 1, gm_mean)

MERGEDnonsinglet.age.dds1 = estimateSizeFactors(MERGEDnonsinglet.age.dds, geoMeans = geoMeans) MERGEDnonsinglet.age.dds2 = DESeq(MERGEDnonsinglet.age.dds1, test="Wald", fitType="parametric")

MERGEDnonsinglet.age.dds.res = results(MERGEDnonsinglet.age.dds2, cooksCutoff = F) View(MERGED_nonsinglet.age.dds.res)

...and I've read conflicting information on whether or not you should transform your data. Also the output results are a bit weird. If I look at differences between, for example, males and females, then do the same for adults vs juveniles...3 ASVs appear significant with sex and 2 ASVs with age class, but those 2 ASVs for age class are the same as those for sex - if that makes sense? Additionally, those ASVs are not what I would expect to come out - I can see some ASVs that clearly have a hugely different abundance that aren't picked up.

After I looked in to this more, I've seen that there is some debate over whether or not you should transform the data using geomeans so I'm unsure how to proceed. Is it this geomeans step that could be giving these weird results? If yes, what should I be doing differently to make it work?

Thanks for any advice you may have, Sophie

deseq2 • 590 views
ADD COMMENT
0
Entering edit mode
@mikelove
Last seen 1 hour ago
United States

I’m not convinced DESeq2 / RNAseq Methods is the right tools for this type of data.

That said you aren’t using the poscounts method for size factor estimation above, which was created for this handoff. So that will make the code work but I’m not sure it’s the most appropriate or best model.

ADD COMMENT

Login before adding your answer.

Traffic: 735 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