edgeR - making 12 contrast but only get 5 contrasts as a result
1
0
Entering edit mode
nooshin ▴ 300
@nooshin-5239
Last seen 6.0 years ago

Hi,

I'm using edgeR and I aim to find differentially expressed genes for 12 contrasts, like following:

fit <- glmFit(dge,design)
cont <- c('oneA-oneB',
           'twoA-twoB',
           'threeA-threeB',
           'oneA-twoA',
           'oneA-threeA',
           'twoA-threeA',
           'oneB-twoB',
           'oneB-threeB',
           'twoB-threeB',
           '(oneA-oneB)-(twoA-twoB)',
           '(oneA-oneB)-(threeA-threeB)',
           '(twoA-twoB)-(threeA-threeB)')
c.fit <- glmLRT(fit, contrast = mc)
res <- topTags(c.fit,adjust.method = "fdr",n = Inf)$table
> dim(res)
[1] 19859     9

Which includes the first 5 contrasts and the following statistics:

logCPM       LR        PValue           FDR

I can do the contrast one by one like

c.fit <- glmLRT(fit, contrast = mc[,1])

But I'm interested to do it all at once in one run.

Could anybody help and show where I do a mistake or what I have to do to get all my contrasts at once?

Thanks in advance and looking forward.

 

 

r edger de contrast design and contrast matrix • 1.6k views
ADD COMMENT
0
Entering edit mode
> sessionInfo()
R version 3.3.0 (2016-05-03)
Platform: x86_64-w64-mingw32/x64 (64-bit)
Running under: Windows 7 x64 (build 7601) Service Pack 1

locale:
[1] LC_COLLATE=English_United States.1252  LC_CTYPE=English_United States.1252   
[3] LC_MONETARY=English_United States.1252 LC_NUMERIC=C                          
[5] LC_TIME=English_United States.1252    

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

other attached packages:
[1] edgeR_3.14.0 limma_3.28.2

loaded via a namespace (and not attached):
[1] tools_3.3.0     splines_3.3.0   grid_3.3.0      locfit_1.5-9.1  statmod_1.4.24  lattice_0.20-33
ADD REPLY
2
Entering edit mode
Aaron Lun ★ 28k
@alun
Last seen 5 hours ago
The city by the bay

Running all the contrasts together in glmLRT does not do the same thing as running them separately. Each contrast defines a null hypothesis, and when you run each contrast separately, you're testing whether the lone corresponding hypothesis is true for each gene. However, when you run them together, you're setting up for an ANOVA-like test where the joint null hypothesis is that all of the individual null hypotheses have to be true. The output of the ANOVA-like test will tell you which genes have some DE somewhere in one or more of the contrasts, but will not specify which contrast contains the DE. Whether this is desirable or not depends on what you want to do, but in general, running pairwise tests to get inferences for specific contrasts is more useful.

Now, assuming that you do indeed want to run an ANOVA-like test, recall that each of your contrasts defines a null hypothesis. For example, 'oneA-oneB' would be saying that oneA - oneB = 0. Then:

oneA - oneB = 0 # and
twoA - twoB = 0 # then obviously,
(oneA - oneB) - (twoA - twoB) = 0

In other words, the last contrast doesn't need to be specified, as it is implied by the first two if they are present in your contrast matrix. edgeR is smart enough to know this and throws away the last one to avoid wasting time on having to do redundant things. Losing the third contrast doesn't affect the p-value calculation, though it is a bit annoying if you're interested in the corresponding log-fold changes; the solution is just to go through the missing contrasts with glmLRT and compute the log-fold changes directly.

ADD COMMENT
0
Entering edit mode

Thanks you very much for very nice and comprehensive explanation.

Now I understand it much better :)

 

ADD REPLY
0
Entering edit mode

Just as a very short question, is this concept in edgeR then different from limma?

Since in limma, I get all specific contrasts separately when I run:

c.fit=contrasts.fit(fit,mc)

 

ADD REPLY
0
Entering edit mode

Yes, the behaviour of the code is different (though the concept of the ANOVA is similar).

ADD REPLY
0
Entering edit mode

Thanks a lot. It was my mistake to use edgeR similar as limma. Now all is clear, tnx

ADD REPLY

Login before adding your answer.

Traffic: 394 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6