[Michael Love]

The estimator for the amount of shrinkage changed (this and all changes are noted in NEWS). Note that rlog is just one of the possible transformations, which we have found useful in many situations, but certainly it can shrink too much or too little given the diversity of datasets.

I'd try out the varianceStabilizingTranformation, or even log2(counts(dds, normalized=TRUE) + pc) for varying pseudocount.

[TomHarrop]

Changes in r-log transformation?

Hi,

Sorry to dredge up this thread again but I have a question about the r-log transformation in the current version of DESeq2 (1.6.1).

I have a set of 8 libraries. Two of them were produced from degraded RNA. Using the previous version(s) of DESeq2, when I made density plots of r-log expression values, the two libraries that were produced from "bad" RNA had distributions that were obviously different from the libraries made from "good" RNA. This was accompanied by clear separation between the "bad" libraries and their matching biological replicates on a PCA plot made from r-log expression values.

With the current version, in density plots produced with the exact same code, all the distributions look the same as each other (i.e. the "bad" libraries now look the same as the "good" ones). The differentiation on the PCA plot is still there, though.

I am wondering if this observation can be explained by changes in the r-log transformation in the most recent version(s) of DESeq2? Unfortunately I can't say exactly which version I was using before, but it would have been the current Bioconductor release around the start of September.

Please let me know if any more information is required.

Thanks a lot for reading!

Tom

[dpryan]

Unique counts. If you want to use the effective counts then put them through voom() and use limma.

[chris_s]

I'm new in DESeq2 and i'm having problems to create my dds.

I have a matrix with eXpress results (read counts) with 2 conditions (SP, EB) and 2 replicates/condition.

This is what i'm doing:

Cele_SPvsEB = read.csv (file.choose(), header=TRUE, row.names=1)
CeleDesign <- data.frame(
  row.names = colnames(Cele_SPvsEB),
  condition = factor(c("SP", "SP", "EB", "EB")))
dds <- DESeqDataSetFromMatrix(countData = Cele_SPvsEB,
                              colData = CeleDesign,
                              design = ~ condition)
dds

And i'm getting this error message:

Error in get(name, envir = asNamespace(pkg), inherits = FALSE) : 
  object '.setDummyField' not found

Any help?
Thanks.

Hi,

I am in a fairly similar situation. I have used eXpress to generate read counts and they specifically recommend using the effective counts output for differential expression analysis. However, when I try to create the DESeqDataSet from the count matrix I get an error saying that :

Error in DESeqDataSet(se, design = design, ignoreRank) : 
  some values in assay are not integers

I know that DESeq2 normally takes in raw read counts, which probably explains why it is looking for integers, but since eXpress advises to use effective read counts which are not integers (see below): "The estimated number of fragments generated from this target in the sequencing experiment, adjusted for fragment and length biases. In other words, this is the expected number of reads from the experiment if these biases did not exist. This is the value recommended for input to count-biased differential expression tools."

I am keen to use the effective reads as this is applies a correction to an important potential bias.

Here is a subset from my countTable:

          A1_effcounts A3_effcounts A10_effcounts
evgtrinCL81615Contig1_79790      2.583273     0.000000      0.000000
evgtrinCL11213Contig1_47858     14.839541     9.894940     19.304148
evgtrinCL34839Contig1_73815     18.764835     0.000000      0.000000
evgtrinCL90968Contig1_58906      0.000000     0.000000      6.883716
evgtrinCL83431Contig1_102122     1.358742     1.354169      1.150552
evgtrinCL34506Contig1_139871     4.719385     2.416909      0.000000

and my design file:

   treatment  time
A1_effcounts        SED   24h
A3_effcounts        SED Week2
A10_effcounts       CON     0
A11_effcounts       CON   24h
A13_effcounts       CON Week2
H1_effcounts        SED   24h
H3_effcounts        SED Week2
H10_effcounts       CON     0
H11_effcounts       CON   24h

So do you use the total counts, the unique counts or is there a work around to use the effective counts from the eXpress output (e.g. rounding them to the nearest integer)?

Thanks for the help!

Cheers,

Chris

[IsBeth]

I'm reading your paper now. I will post my MA plot and ask more as soon I've read it. Thank you for your reply and for the link of your paper! =)

[Simon Anders]

If you don't show us the MA plot, we can't tell you whether it looks right. That's why I suggested you post it. Also tell us about the biology of your experiment.

EdgeR's log fold-changes are not shrunken, DESeq2's are. They are a new feature in DESeq2, and explained in our paper.

Hence, sorry if the following question sounds a bit direct but: By "I don't understand the concept of shrunken fold changes", do you mean that you have read our preprint but found the explanation unclear (in which case we should improve it) or that you have not yet looked at our paper (in which case we would be disappointed that nobody reads what we wrote)? ;-) You don't need to read the whole thing; just the section on shrunken fold changes.

[IsBeth]

In my MAplot are differential expressed genes (in red). Yes, there are also negative values. I think I don't understand the concept of shrunken fold changes, still a newbie

It just seemed strange to me as I compared the lfc values after using edgeR and DESeq2

[Simon Anders]

Maybe there is not much happening in your data?

Can you post an MA plot?

Also, make sure you understand the conspt of shrunken fold changes. (See our paper for details.)


Edit: Values between 0 and 1? No negative values? Sounds like you are looking at the p values, not the log fold changes.

[IsBeth]

Hello,

I'm using DESeq2 and so far it seems to work well. But why are the log fold changes so low? They have values between 0 and 1, even when the genes are differentially expressed.

[kmcarr]

Wasn't Michael's reply (the next post) in regard to that?

In the words of Homer Simpson, "DOH!!"

Thank Devon, and apologies to everyone.

[dpryan]

I'm bumping up KHubbard's question because I have the same one and it seems it didn't get answered the first time around. I have looked at the structure of the DESeqDataSet object (best as an R neophyte can) and don't see them in there. Is there any way to get these values in DESeq2?

Thanks.

Wasn't Michael's reply (the next post) in regard to that?

[kmcarr]

Hello,

In DESeq, it was possible to view/export the baseMean values for each condition (baseMeanA and baseMeanB) in addition to the overall baseMean. In DESeq2, I have only been able to view/export the overall baseMean. Is it possible to acquire the baseMean values for each condition as well? Thanks!

I'm bumping up KHubbard's question because I have the same one and it seems it didn't get answered the first time around. I have looked at the structure of the DESeqDataSet object (best as an R neophyte can) and don't see them in there. Is there any way to get these values in DESeq2?

Thanks.

[Michael Love]

hi,

Sorry, this was a simple bug for when number of samples equaled number of coefficients, which also came up on the Bioconductor mailing list. I pushed a fix to the release branch.

http://comments.gmane.org/gmane.science.biology.informatics.conductor/51865

Mike

[mbahin]

Hi all,

I just found the answser, I was trying to use DESeq2 with only two count files (waiting for the other one to come).
Now that I am using it with 6 files it's ok. It looks like it was the issue.

Cheers

[mbahin]

Hi all,

I am trying to use DESeq2 just after HTSeq-count.

I am following the tuto in the DESeq2 vignette but when I try to do the:

res <- results(dds)

Although I think I just followed the tuto properly, I get this error message

Error in is.numeric(cooksCutoff) : 'cooksCutoff' is missing

I can't find any answer about this error...

Cheers