Question: How to interpret DESeq2 result?
gravatar for helen
5 weeks ago by
helen30 wrote:

Hi all,

I have 6 subjects 1-6, and would like to identify the differentially expressed genes in condition A vs.B and consider the gender effect.

I arrange my coldata as follow to mimic the example in this section of DESeq2 manual: , which gender can be approximate to grp, subject to ind, condition to cnd

gender  subject  condition
     M        1          A
     M        1          B
     M        2          A
     M        2          B
     M        3          A
     M        3          B
     F        4          A
     F        4          B
     F        5          A
     F        5          B
     F        6          A
     F        6          B

Add a column subject.n

coldata$subject.n <- factor(rep(rep(1:3,each=2),2))

so the new coldata would be:

gender  subject  condition  subject.n
     M        1          A          1
     M        1          B          1
     M        2          A          2
     M        2          B          2
     M        3          A          3
     M        3          B          3
     F        4          A          1
     F        4          B          1
     F        5          A          2
     F        5          B          2
     F        6          A          3
     F        6          B          3

and the columns of my model.matrix(~ gender + gender:subject.n + gender:condition, coldata) are

[1] "(Intercept)"             "gendermale"                  
[3] "genderfemale:subject.n2" "gendermale:subject.n2"    
[5] "genderfemale:subject.n3" "gendermale:subject.n3"    
[7] "genderfemale:conditionB" "gendermale:conditionB"

so my design formula would be

dds <- DESeqDataSetFromMatrix(countData = cts,
                              colData = coldata,
                              design = ~ gender + gender:subject.n + gender:condition)

then I filter out low read counts and perform differential expression analysis

keep <- rowSums(counts(dds)) >= 10
dds <- dds[keep,]
dds <- DESeq(dds)

and then extract the result

res <- results(dds, contrast=list("gendermale.conditionB","genderfemale.conditionB"))


out of 18461 with nonzero total read count
adjusted p-value < 0.1
LFC > 0 (up)       : 13, 0.07%
LFC < 0 (down)     : 0, 0%
outliers [1]       : 0, 0%
low counts [2]     : 11095, 60%
(mean count < 565)

only 13 up-regulated genes are found (I was expecting hundreds to thousands of DE genes)

Does anyone know how to interpret this result?

Does this mean under the gender-specific condition effect, comparing conditions A to B would render 13 up-regulated genes? Does this method reach my aim? And why there are so few DE genes?

rna-seq deseq2 • 214 views
ADD COMMENTlink modified 16 days ago • written 5 weeks ago by helen30
gravatar for Asaf
5 weeks ago by
Asaf8.3k wrote:
  1. You mean Sex, not Gender
  2. Why are you expecting a lot of genes? Are M:B and F:B so different? More than M:A and F:A? Try to plot a PCA of the data and see the structure of it.
  3. Getting only 13 genes is usually great news and it means you pinpointed the exact effect, getting a lot of genes means your conditions are too different to learn something useful.
ADD COMMENTlink written 5 weeks ago by Asaf8.3k
  1. Yes, Sex would be more appropriate. Thanks for pointing it out. For 2. & 3. I think I misunderstood what the result is trying to tell. I thought this design is to compare condition A vs. B under specific sex, but in fact this interpretation is not correct. If I understand correctly, the fewer DEGs indicates the sex difference effect is very minimal in this experiment. 13 genes now seems an acceptable number.
ADD REPLYlink written 5 weeks ago by helen30

Yes, Sex:Condition interaction will give you the differences between A and B which are sex specific, in your case there are 13 genes that are differentially expressed between A and B in males only.

ADD REPLYlink written 5 weeks ago by Asaf8.3k
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.3.0
Traffic: 966 users visited in the last hour