Search
Question: Treat Method - number of DE genes
0
10 months ago by
Sanches0
Sanches0 wrote:

Hello everybody,

I have been running my microarray data using limma and I'd like to ask a question. What is the best way to reduce the number of differentially expressed genes? As you can see below, I have used the treat method to reduce it and no differentially expressed genes were listed in the Contrast(3).  But,  when I create a plot to represent of differential expression results ( Figure Contrast(3)), it is possible to identify significant genes.  So, Why does it happen? Sorry for my basic question? Thanks in advance.

fit2 <- contrasts.fit(fit, cont.matrix)

efit2 <- eBayes(fit2)

summary(decideTests(fit2))

summary(decideTests(fit2))
Contrast(1)                    Contrast (2)           Contrast(3)
-1       2658                           6566                   942
0       56790                          51236                 59649
1        3528                           5174                   2385

tfit <- treat(fit2, lfc=1)

dt <- decideTests(tfit)

summary(dt)

summary(dt)

Contrast(1)                 Contrast(2)         Contrast(3)
-1         57                            783                     0
0       62672                          61728                 62976
1         247                            465                     0

plotMD(tfit, column=3, status=dt[,1],main=colnames(tfit)[3], xlim=c(5,18))

I do not know how to insert A figure here.

So, let me try to explain the result from the figure. So, when looking at the figure (Log-fold change x Average log-Expression), I can see a few genes highlighted in green and red. These genes are between the range of  -2<log-Fold change <2. My question is: If I got 0 up- and down-regulated genes in the contrast (3), why are these genes highlighted in the contrast(3) figure?

modified 10 months ago by Gordon Smyth34k • written 10 months ago by Sanches0

It isn't clear from your post what the problem is. Evidently the software is working correctly. What are you trying to do that you haven't been able to do?

You've chosen a very high value for 'lfc', so getting no DE genes from treat() isn't surprising.

Hi Gordon,

The software is working correctly.

You just answered my question. But If I set lfc=log2(1.5), I would increase false discovery rate, right? So, in your opinion, what is the best lfc cutoff?

1

treat() always controls the FDR correctly whatever lfc you choose. Relative to that lfc. The whole idea is that you choose a lfc cutoff that is meaningful to you and appropriate for your dataset.

If there was one "right" lfc that was always best then we would set it for you and it wouldn't be an argument for the treat() function.

Dear Gordon,

I just edited this post. Is it possible you help me to understand my further question?

1
10 months ago by
Gordon Smyth34k
Walter and Eliza Hall Institute of Medical Research, Melbourne, Australia
Gordon Smyth34k wrote:

treat() gives more emphasis to genes with larger fold changes, and you can choose whatever lfc seems appropriate depending on how large the changes are in your experiment. In practice, you can choose the lfc cutoff to get the number of DE genes that you want to analyse. Personally though, I would very seldom choose lfc larger than about log2(1.5).

You have chosen an unnecessarily high threshold that prevents you from getting any significantly DE genes. That does seem counter-productive.

5 Days Later

You have now edited your question to ask something different. The problem now is that you used plotMD to plot fold changes from contrast 3 but you highlighted genes according to contrast 1. There are significantly DE genes for contrast 1, so these genes are highlighted on the plot. I'm not sure why you set status for a different contrast than was being plotted, but the plot is just doing what you asked it to do.

Dear Gordon,

I just edited this post. Is it possible you help me to understand my further question?

My mistake. You are right. Thanks again.

0
10 months ago by
Denali
Steve Lianoglou12k wrote:

I'm not really sure what the issue is that you are trying to chase down is. I would recommend, however, that you don't keep updating your question and changing it slightly after interacting with Gordon, but maybe ask a different question that is specific to the new problem you think you are facing.

That having been said, you are interacting with the results of your analysis from something of "a distance" via the plotMD and decideTests functions. You should rather get the full table statistics for your contrast of interest via topTable(tfit, "Contrast(3)", n=Inf) and explore the logFCs, pvalues, t-tstats, and adjust pvalues over your genes directly.