[BioC] LIMMA decideTests result zero from contrast matrix

[James W. MacDonald]

Hi Ekta,

On 5/8/2012 1:49 AM, Ekta Jain wrote:
> Hi Steve,
> Definitely not annoy people here on. I am being given help, which is nice. Apologies for the same.
>
> I wouldn't completely rely on Wikipedia, ever :)
>
> I am a warm person.
>
> Yes, (x,...) is generated from decideTests(fit). The summary(decideTests, lfc=0, method="separate",adjust.method="BH",p.value=0.05) gives 0 genes for upregulation and downregulation. I can only set the threshold for lfc, p-value and not for *adjusted p-values* - http://127.0.0.1:21780/library/limma/html/decideTests.html

It's ironic that you post a (busted) link to a help page that you have 
evidently not read yourself. From that man page:

Usage:

      
decideTests(object,method="separate",adjust.method="BH",p.value=0.05,lfc=0)

and

adjust.method: character string specifying p-value adjustment method.
           Possible values are  "none" ,  "BH" ,  "fdr"  (equivalent to
            "BH" ),  "BY"  and  "holm" . See  p.adjust  for details.

So by _default_ decideTests() uses a BH adjusted p-value to filter genes.

Both Steve and I have pointed out numerous times that you are doing two 
different things and expecting the same result. Yet you persist in 
ignoring what we have been telling you. You may have got the same 
results before with this code, but that was a coincidence. The reason  
you are "clueless why it wouldn't work" is because you seem bent on 
ignoring what people tell you. So let me try once more. Please read the 
following very closely...

When you are using topTable you are not filtering on p-values. In fact 
you aren't filtering at all and stating that there are lots of genes 
with a |fold| > 2. You then use decideTests(), which filters on adjusted 
p-values. Since your smallest adjusted p-value is something like 0.999, 
everything gets filtered out. So you don't filter using topTable(), but 
you do filter using decideTests(). That is why you get different results!

Best,

Jim




>
> Working on the LIMMA source code is next on the list.
>
> Thank you much for your help and time.
>
> Best,
> Ekta
>
>
>
>
> -----Original Message-----
> From: Steve Lianoglou [mailto:mailinglist.honeypot at gmail.com]
> Sent: 08 May 2012 10:44
> To: Ekta Jain
> Cc: bioconductor at r-project.org
> Subject: Re: [BioC] LIMMA decideTests result zero from contrast matrix
>
> Hi,
>
> On Tue, May 8, 2012 at 12:42 AM, Ekta Jain<Ekta_Jain at jubilantbiosys.com>  wrote:
>> Hi Steve,
>> I have been using the lists for a long time now and was never told off before for emailing to both R (r-help at r-project.org) and Bioconductor(bioconductor at r-project.org) at the same time. I think, it saves time and is the quickest way to get your query across people since not everyone can be on both the lists.
> It's also a quick way to annoy people ... I'm not trying to reprimand
> you or anything, I'm just pointing out what "the rules of the road"
> are on most mailing lists so that you can learn and stick to them.
>
>> I can most certainly try and post to the relevant list but 95% of the time I am then told to post to the other list since the topic is not relevant here.
> It seems as if you simply have to flip the sign of your rhelp vs. bioc
> classifier and all will be well ;-)
>
>> I am afraid this process at times takes up an extra day before some good soul out there reads my post and offers help. I am sorry but I did not mean to stress anyone.
>>
>> I could have googled 'cross-posting' but I meant to ask in this context and not as a general rule. Ideally, I wouldn't rely on Wikipedia for all information :-).
> I feel like wikipedia gets a worse wrap than it deserves --  but all
> the same, if you search the bioconductor or r-help mailing list for
> the relevant keywords, you'll get many results, eg:
>
> http://search.gmane.org/?query=please+don%27t+cross+post&author=&group=gmane.science.biology.informatics.conductor&sort=relevance&DEFAULTOP=and&xP=cross%09post&xFILTERS=Gscience.biology.informatics.conductor---A
>
>> I've tried changing things and no results. I reckon, it's probably down to the data I have. It works fine on another dataset with decent p-values. I am clueless why it wouldn't work.
> You're getting warmer, here: -
>
>> You asked me for adjusted p-values, why?
> Because it seems as if the default way you are calling vennDiagram
> results in reporting counts for genes that have *adjusted* pvalues
> below a certain threshold. You've said that your analysis results in 0
> such genes, so this would explain why your vennDiagram is turning up
> empty. -
>
> If you read the help in ?vennDiagram, the fact that "x is .. Usually
> created by `decideTests`". Now, looking back at your code, the
> `results` you are passing into your vennDiagram call are coming from a
> call to decideTests which will, by default (which you're doing), only
> pick the genes that pass a certain threshold for their *adjusted*
> pvalues.
>
> To convince yourself of this fact, it wouldn't be a bad idea to
> minimally read the help pages of the functions you're calling. I'd
> also recommend downloading the source code for limma and poke around
> these functions to see what's going on under the hood ... you can, of
> course, accomplish this w/o d/ling the code, but I find that browsing
> through the source code of packages is a good way to learn, as well.
>
> -steve
>

-- 
James W. MacDonald, M.S.
Biostatistician
University of Washington
Environmental and Occupational Health Sciences
4225 Roosevelt Way NE, # 100
Seattle WA 98105-6099