Question

how to make a contrast of uniquely DE genes

0

Entering edit mode

decppced • 0

@7b74a5ea

Last seen 3.6 years ago

United States

Hello all thanks for the help in advance I have 4 groups, two different cell lines (w and b) which were non-treated (nt) and treated (t). I wanna find genes that are uniquely DE in the w-treated group compared to the other groups.

Code should be placed in three backticks as shown below

contr.matrix <- makeContrasts(                                  
  unique_w_treated = w_t- (w_nt + b_nt + b_t)/3,
  levels = colnames(design)

so, when I do this I get 164 DE genes, but doesn't make much sense to me because when I compare w_t vs w_nt alone, I only have 13 DE genes? Am I doing something wrong? Thanks

sessionInfo( ) R version 4.0.2 (2020-06-22) Platform: x86_64-apple-darwin17.0 (64-bit) Running under: macOS Mojave 10.14.6

Matrix products: default BLAS: /System/Library/Frameworks/Accelerate.framework/Versions/A/Frameworks/vecLib.framework/Versions/A/libBLAS.dylib LAPACK: /Library/Frameworks/R.framework/Versions/4.0/Resources/lib/libRlapack.dylib

locale: [1] en_US.UTF-8/en_US.UTF-8/en_US.UTF-8/C/en_US.UTF-8/en_US.UTF-8

attached base packages: [1] grid stats4 parallel stats graphics grDevices utils datasets methods base

edgeR limma RNA-seq DEGseq • 1.7k views

ADD COMMENT • link 3.8 years ago decppced • 0

score 1 · Accepted Answer · 2021-03-17

1

Entering edit mode

James W. MacDonald 67k

@james-w-macdonald-5106

Last seen 9 hours ago

United States

You are comparing the average of three groups to the w_t group. There are many ways that the average of those three groups could be much different, even if there aren't many differences between w_t and w_nt.

Analogously, if you were to compare heights between two gymnastics teams, they might not be that different. But if you compared the heights of a gymnastics team to the average of a gymnastics team and a basketball team you shouldn't be surprised if that comparison is significant.

ADD COMMENT • link 3.8 years ago James W. MacDonald 67k

0

Entering edit mode

That makes sense. It is pretty obvious! lol. I just didn't know how to interpret the code completely thanks

ADD REPLY • link 3.8 years ago decppced • 0

0

Entering edit mode

however, that still doesn't answer the main question. Since this is the wrong contrast for my question, which one is it? thanks

ADD REPLY • link 3.8 years ago decppced • 0

0

Entering edit mode

I don't know. What does 'uniquely DE' mean in this context? Are you looking for genes that are DE in all three comparisons, and have the same sign? Or just DE in all three comparisons? Something else?

ADD REPLY • link 3.8 years ago James W. MacDonald 67k

0

Entering edit mode

thanks for answering. I am trying to find genes that are DE-ed in only one but not the other three conditions when I compare them individually. The condition I am interested. Is w_t. Is that possible? Thanks

ADD REPLY • link 3.8 years ago decppced • 0

0

Entering edit mode

DE means differentially expressed, which by definition is a comparison between groups (e.g., you can restate that as 'expressed differently'). If you want to know the genes that are DE in w_t as compared to each of the other groups, then you can make all three comparisons and then do a Venn diagram to see which genes are significantly different in all three comparisons.

But that brings me back to my original question. Say Gene X is up-regulated when you compare w_t to w_nt and w_t vs b_t (e.g., the expression is higher in w_t than those other two groups), but is down-regulated in w_t vs b_nt (lower in w_t vs b_nt). Is that gene 'unique to w_t'? I mean it's different between w_t and each of the other three groups, but not in the same way.

So what I think you want is to make all the comparisons, run decideTests and then do something like

dt <- decideTests(fit) 
## where fit is an MArrayLM from running eBayes
## if you don't care which direction
ind <- rowSums(abs(dt)) == 3L
## if you do care about directionality
ind <- abs(rowSums(dt)) == 3L

Then what you do is dependent on what you want. If you just want the gene names, you can use that indicator variable to extract things from the genes item in your MArrayLM object. Or you could get the genes and statistics using topTable with sort.by = "none" and subset using the indicator.