[BioC] Variance stabilizing transformation seems not working on NGS (bacterial 16srRNA based) data

[Shucong Li]

Hi Michael,

I have rerun the test with 'blind' set as TRUE.  The graph didn't change.  It is interesting that  rlogTransformation seems to work very well. 

Please allow me ask one more question about differential analysis with DESeq().   I also run the differential analysis with the same set of data. Original counts was used as suggested.  I'm wondering if any variance stabilizing transformation is involved in DESeq() procedure?

Shucong

On Nov 27, 2013, at 3:31 PM, Shucong Li <Shucong.Li at umanitoba.ca> wrote:

> Hi Michael,
> 
> The supplemental file 3 in the manuscript is not available yet. I was told it will be available in a week or so.     I will try to rerun the test with blind set as TRUE.
> 
> 
> Shucong 
> On Nov 27, 2013, at 3:13 PM, Michael Love <michaelisaiahlove at gmail.com> wrote:
> 
>> hi Shucong,
>> 
>> I'm replying and including the Bioconductor list. 
>> 
>> It depends on what you mean by 'not working', I see that the plots of standard deviation over mean are not flat using the two transformations from DESeq2. For both transformations, you are specifying blind=FALSE, which means that the group means are taken into consideration when estimating the dispersion-mean relationship. Then for the plot, you are calculating the standard deviation over all the samples, regardless of groups. I would say that the plot does not imply that the transformation did not work, as the dispersion-mean relationship was not calculated with respect to inter-group differences because 'blind' was set to FALSE. I can't tell if the manuscript you reference used blind dispersion estimation or not when they were using the VST (I couldn't find a link to supplemental file 3).
>> 
>> hope this helps,
>> 
>> Mike
>> 
>> On Wed, Nov 27, 2013 at 3:21 PM, Shucong Li <Shucong.Li at umanitoba.ca> wrote:
>> Hello Dr. Love,
>> 
>> I tried to used DESeq2 to do VST transformation with my NGS (bacterial 16srRNA based) data. I came to this approach after I read the manuscript written by Profs McMurdie and Holmes titled "Waste Not, Want Not: Why Rarefying Microbiome Data Is Inadmissible".  However, the VST transformation did work as what suggested in the Vignettes. I started wondering this VST from DESeq2 may not work to my data.  I emailed Prof. Holmes this morning to get her opinion. By the mean time, it would be great i can have your opinion regarding this issue.
>> 
>> I also attached the R codes and data (R workspace file).
>> 
>> Thank you in advance.
>> 
>> Shucong
>> 
>>  
>> 
>> 
>> Shucong Li   PhD Research Associate
>> -----------------------------------
>> Department of Animal Science
>> University of Manitoba
>> 12 Dafoe Rd, Winnipeg, Manitoba
>> R3T 2N2, Canada
>>  
>> Tel: (204) 474 9841
>> Fax: +1(204) 474 7628
>> Email: shucong.li at ad.umanitoba.ca
>> -----------------------------------
>> 
>> 
>> 
>> dds <- DESeq(dds, fitType="local”)
>> 
>> 
>> 
>> # ddsClean <- replaceOutliersWithTrimmedMean(dds)
>> 
>> # ddsClean <- DESeq(ddsClean,fitType="local") 
>> 
>> 
>> 
>> plotDispEsts(dds)
>> 
>> 
>> 
>> rld <- rlogTransformation(dds, blind=FALSE)
>> 
>> vsd <- varianceStabilizingTransformation(dds, blind=FALSE)
>> 
>> 
>> 
>> 
>> 
>> library("vsn")
>> 
>> par(mfrow=c(1,3))
>> 
>> notAllZero <- (rowSums(counts(dds))>0)
>> 
>> meanSdPlot(log2(counts(dds,normalized=TRUE)[notAllZero,] + 1),
>> 
>>            ylim = c(0,2.5))
>> 
>> meanSdPlot(assay(rld[notAllZero,]), ylim = c(0,2.5))
>> 
>> meanSdPlot(assay(vsd[notAllZero,]), ylim = c(0,2.5))
>> 
>> 
>> 
>> 
>> 
>> par(mfrow=c(1,2)) 
>> 
>> dist <- dist(t(assay(rld)))
>> 
>> fit <- hclust(dist, method="average")
>> 
>> plot(fit, labels=dds at colData$treat)
>> 
>> 
>> 
>> dists <- dist(t(assay(vsd)))
>> 
>> fit <- hclust(dists, method="average")
>> 
>> plot(fit, labels=dds at colData$treat) 
>