1
$\begingroup$

I obtained a matrix of RNA-seq count data that has been normalized by DESeq2's median of ratio method.

I know that DESeq2 wants to take in un-normalized counts, but I do not have access to those data.

How do I best proceed here if I want to perform DEG analysis using DESeq2?

I know I can always start from .fastq files, but that would be so much extra work.

I don't think I can un-normalize the data because I don't know the scaling factors.

Any suggestions?

$\endgroup$
1
$\begingroup$

You can read in the normalized count table and don't normalize the data, but my advice here is not to do that. The rounding of the normalized matrix introduces some noise, but I think the larger issue is how sure are you that the table you are working with, is exactly a count table of normalized counts from DESeq2 ?

It really doesn't take much time to align and get count table again, rather than going back and forth wondering whether you did the normalization correctly.

I can show below with an example how it will look like, first we make an example dataset and obtain the normalized counts:

library(DESeq2)
set.seed(111)
sz = runif(6,min=0.5,max=1.5)
x = makeExampleDESeqDataSet(sizeFactors=sz,m=6)
x = estimateSizeFactors(x)
ncounts = counts(x,normalize=TRUE)
x = estimateDispersions(x)

Now suppose we feed this into a new DESeq2 object, but before that you need to do is check most of the size factors are close to 1:

estimateSizeFactorsForMatrix(round(ncounts))
 sample1  sample2  sample3  sample4  sample5  sample6 
1.039179 1.045991 1.048276 1.045872 1.032420 1.052006

If you see something above that suggest a deviation from 1, something is wrong from the normalization and do not even think about using it.

Below you can use the normalized counts and you can see the difference is not huge:

dds = DESeqDataSetFromMatrix(round(ncounts),colData(x),~condition)
sizeFactors(dds) = 1
dds = estimateDispersions(dds)

Compare the new count table:

head(counts(dds,normalize=TRUE))
      sample1 sample2 sample3 sample4 sample5 sample6
gene1       0       2       0       1       0       0
gene2       9       8      16       6       1       2
gene3      23      40       8      38      11      41
gene4       0       0       0       6       8       0
gene5       3       2       0       5       1       2
gene6       0       1       3      19       4       9

And your input:

head(ncounts)
        sample1    sample2   sample3    sample4   sample5   sample6
gene1  0.000000  2.2658464  0.000000  0.9310803  0.000000  0.000000
gene2  8.693480  8.3081033 16.251453  5.5864817  1.121781  2.086619
gene3 23.472395 40.0299524  7.584012 38.1742919 11.217808 40.689076
gene4  0.000000  0.0000000  0.000000  5.5864817  7.852465  0.000000
gene5  3.477392  1.5105642  0.000000  4.6554014  1.121781  2.086619
gene6  0.000000  0.7552821  3.250291 18.6216058  4.487123  9.389787

We can check the dispersion, for lowly expressed genes you see some difference, but on the whole difference

with(elementMetadata(x),plot(log10(baseMean),dispGeneEst,pch=20,cex=0.7,col="#FFC93C80"))
with(elementMetadata(dds),points(log10(baseMean),dispGeneEst,pch=20,cex=0.7,col="#DBF6E980"))

enter image description here

Again, I would say I agree with @ATpoint comments that if you can, start from the fastq again.. don't work with data you have doubts about.

$\endgroup$
1
-1
$\begingroup$

You might be able to cheat by giving DESeq a custom size factor matrix where all the values are 1. Then the software will feel like it's done the normalization step, but the 'normalized' values will be what you want.

$\endgroup$
2
  • $\begingroup$ I'll give this a try. Caveat is that I need to first round the normalized counts to integer, but I guess I'll have to live with that for now. $\endgroup$
    – Nomad420
    Mar 22 at 22:27
  • 1
    $\begingroup$ Nah, don't do that. This is most likely going to mess up things as the DESeq2 model runs on raw counts and everything (normalization, covariates) is done via offsets of the model. The DESeq2 author has been asked this a million times at Bioconductor support page, and if things were that easy then it would be in the vignette, at least as a quick'n'dirty workaround. Start from fastq, you are only adding uncertainty with these custom approaches that in the end cost you time and give inaccurate results. Learn from this, try to backup files and save your R environment incl. the raw data. $\endgroup$
    – ATpoint
    Mar 23 at 11:10

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.