Question: expression matrix with varianceStabilizingTransformation
gravatar for lirongrossmann
11 months ago by
lirongrossmann10 wrote:

Hey all,

I have an rna-seq expression matrix and I used Deseq2 to compare gene expression between two groups of in my dataset. I then wanted to to sort the genes I got from Deseq2 according to their levels of expression, and I used  the varianceStabilizingTransformation (vsd) to get the normalized expression data. 

The problem I had is the following:

When I changed one of the groups (omitted few samples), and applied again the vsd function, I saw that the expression levels of some genes actually changed for samples that were not removed from the dataset. That is, the expression values for the samples that remained in the dataset were changed just by omitting other samples from the dataset.

Is there a way to get a normalized matrix with expression values for each sample that does not depend on other samples? 

This is my code:

ep<-read.table("expression.txt",header = TRUE, row.names = 1)
dds <-DESeqDataSetFromMatrix(countData = ep,colData = cp,design =~Group)
dds <- estimateSizeFactors(dds)
vsd <- varianceStabilizingTransformation(dds)
ADD COMMENTlink modified 11 months ago by Michael Love19k • written 11 months ago by lirongrossmann10
gravatar for Wolfgang Huber
11 months ago by
EMBL European Molecular Biology Laboratory
Wolfgang Huber13k wrote:

Dear lirongrossmann

The transformation parameters depend on the statistical distribution of the data, so it is to be expected that the transformation changes (a bit) if data are added or removed, especially if these make replicate variances look higher or lower.

If your goal is sorting genes by overall expression, you can do something like

mcr = matrixStats::rowMedians(counts(dds, normalized = TRUE))

and sort by mcr.

Kind regards


ADD COMMENTlink written 11 months ago by Wolfgang Huber13k

Thank you !


ADD REPLYlink written 11 months ago by lirongrossmann10

 Thanks Dear Wolfgang,

One clarification please: if I want to use the genes I got for my Deseq2 to build a machine learning model using a training and validation set, would you use the counts (using the counts data you recommend) or the transformed version of of the data (using the variance stabilizing function) as an input to the learning algorithm?


ADD REPLYlink modified 11 months ago by Michael Love19k • written 11 months ago by lirongrossmann10

We recommend using variance stabilized, transformed data for downstream methods that benefit from homoskedasticity (same scale of variance across the dynamic range)

ADD REPLYlink written 11 months ago by Michael Love19k
Please log in to add an answer.


Use of this site constitutes acceptance of our User Agreement and Privacy Policy.
Powered by Biostar version 2.2.0
Traffic: 289 users visited in the last hour