Re: A question about LEfSe

[Curtis Huttenhower]

Hi Ali - likewise, thanks for getting in touch, and I hope the tools have generally been behaving themselves.  It sounds like you're taking the right approach to getting all possible outputs from LEfSe by changing those parameters, but there may be other filtering criteria that are obscuring the rest, depending on how sparse they are.  I've forwarded your email along to the lefse-users list, which is a good place for questions like this - LEfSe folks, any ideas on what might cause this, or would it help to take a look at a snippet of specific data?

Many thanks -

Curtis

On Thu, Feb 13, 2014 at 1:39 PM, Ali May <a....@vu.nl> wrote:

Hi!

Thanks a lot for developing humann, LEfSe and all the other goodies! I have a question about LEfSe.
I downloaded the LEfSe source code and I'm able to run it and reproduce the results I obtain from the web server. What I really would like to achieve with LEfSe is to obtain logarithmic LDA scores and p-values not only for features that are discriminative, but basically for all of them. For this, I tried to max the alpha parameters (to 1.0) and floor the threshold on the absolute value of the logarithmic LDA score (to 0.0), but it did not help getting an LDA score and a p-val for each feature in the input dataset. With default values (alpha = 0.05 for anova & wilcoxon, 2.0 for LDA score) I think I could get around LDA scores and p-vals of ~100 features in the output, and when I used the extreme values mentioned above, this number climbed up to ~270. It'd be awesome if I could get the same information for the rest of the ~4000 features.

Is this any possible? If you could point me places at the source code I can try to make the necessary changes. Any kind of input at this stage is much appreciated.

Best regards,
Ali May
PhD cand. VU University Amsterdam

[Ali May]

Hi all,

Hoping that it will help you help me, I am attaching the raw and processed (by LEfSe's format_input.py) input files, along with the output file I get from LEfSe. Here are the commands I used:

/nsegata-lefse-9e29b430108e/format_input.py lefse_raw_input.csv lefse_input.in -c 1 -s 2 -o 1000000
./nsegata-lefse-9e29b430108e/run_lefse.py lefse_input.in lefse_output.tab -a 1.0 -w 1.0 -l 0.0 --wilc 1 -r lda -e 1 --min_c 1 -y 0
Number of significantly discriminative features: 273 ( 3952 ) before internal wilcoxon
Number of discriminative features with abs LDA score > 0.0 : 273

Let me briefly describe the dataset. There are two conditions (and therefore two classes), and in each condition there are 4 biological replicates (which I denote as subclasses using -s 2 in format_input.py). SUBJ1 in CLASS1 is the same subject with SUBJ1 in CLASS2. If you look at the features you will realize that the dataset is indeed quite sparse. This is an unfortunate consequence of very shallow sequencing, and I'm afraid is the main (maybe the only?) reason for me not being able to get p-values for all my features. After a careful look at the LEfSe code and the verbose output I realized that

test_rep_wilcoxon_r function in lefse.py

returns false for the huge majority of features, even when the alpha parameters for the anova and wilcoxon are set to 1.0 and the LDA threshold is set to 0.0. I guess this means that for some reason these features either failed to 'take' the test or..?

Does this mean that there's no hope for being able to obtain a p-value and a logarithmic LDA score for these features?

Thanks in advance!

Best,
Ali

[Ali May]

Reading the earlier posts in the group made me realize that I've made a simple, though fundamental mistake while denoting the properties of my subjects. As I noted in my previous post, I specified the biological samples as subclasses... Why on earth this sounded like a good idea to me in the first place, I have no clue. Now I changed

./nsegata-lefse-9e29b430108e/format_input.py lefse_raw_input.csv lefse_input.in -c 1 -s 2 -o 1000000

to

./nsegata-lefse-9e29b430108e/format_input.py lefse_raw_input.csv lefse_correct_input.in -c 1 -u 2 -o 1000000

and then

./nsegata-lefse-9e29b430108e/run_lefse.py lefse_correct_input.in lefse_correct_output.tab -a 1.0 -w 1.0 -l 0.0 --wilc 1 -r lda -e 1 --min_c 1 -y 0
Number of significantly discriminative features: 2858 ( 3952 ) before internal wilcoxon

And it's been running for almost an hour.

Though you can still help by explaining from which test the value in the last column of the output comes from.

e.g.
FEAT67    2.66413506676    CLASS2    2.50978535187    0.321121070762

Thanks!

[Nicola Segata]

Hi Ali,

 thanks for getting in touch. Yes, as you found out, specifying the replicates as subjects seems to be the best option in your setting. The p-value is the one from the KW test. As far as computing the LDA score for many many features, I confirm that the computation is rather slow and that the effect size estimation works best when only significant features are considering because of the confounding effect of the other features.

many thanks

Nicola

[Allison H]

Hi,

I have a similar problem - I have an OTU table with species whose abundance I've found to be significantly different between two classes using a Mann-Whitney test, but when I try to run the dataset in LEfSe (with no subclasses), they do not come up as significant. I've also tried to relax the alpha values and LDA score threshold all the way to 1 and 0.1, respectively, which, I would think, would make just about everything come up as significant, but still, there are only 1-2 significantly differential features. 

I've used LEfSe before, and I don't think it's a formatting issue, so I was just wondering if I am missing some step (in the case of no subclasses) beyond the KW test (which is really just a MW test if there are only two classes) and the LDA analysis.

Thanks.

-Allison

[Nicola Segata]

Hi Allison,

 thanks for getting in touch. What is the p-value you get with the Mann-Whitney test for the significant OTU? If it is close to 0.05 it may be a numerical approximation in one of the tests. I'd be happy to take a look at the dataset if you like.

many thanks

Nicola

[Allison H]

Thanks, Nicola. I've attached one of the taxa summaries from the dataset. When I run a simple Mann-Whitney comparing case and control samples, Ureaplasma, for example, comes up as significantly different with a p-value of 0.006. However, when I run the dataset in LEfSe, Ureaplasma doesn't come up as a significantly differential feature no matter how much I relax the alpha value/LDA scores. Thanks again for looking at it.

On Thursday, February 13, 2014 2:34:14 PM UTC-5, Curtis Huttenhower wrote:

[Nicola Segata]

Hi Ali,

 LEfSe has a check looking whether the median are different when very few observations are available in two classes. For a bug, this check was maintained also when many examples are available (as in your case). In practice this occurs when the medians of the two classes are both zero (as in the Ureaplasma case). I just fixed the bug. The bitbucket code is updated, but if you are using LEfSe on galaxy you may need to wait few days for the update.

Also, it seems you have two entries for Ureaplasma (one with an additional "|" at the end) and you may want to fix this in your input file.

thanks a lot for your help

Nicola

[Allison H]

Thanks so much, Nicola. That seems to be working much better. 

-Allison

[Allison H]

Hi Nicola, 

I've been using LEfSe successfully on several datasets, but I've been running into a similar problem as described above with a new dataset on which I've been working: I run the formatted dataset and get no significantly differential features between classes even though running a basic Mann-Whitney shows that certain groups (e.g., the genus Treponema) have significantly different abundance between classes (p<0.0001). I've attached one of the datasets with which I've had this problem if you wouldn't mind taking a look. Thanks so much!

-Allison

[Nicola Segata]

Hi Allison,

 looking at your data, it seems that row one is "Bacteria" which is almost always 0. Probably you meant ("Bacteria|OTU_XXX")? Otherwise when LEfSe sums up the internal values there will be problem. I think you have two possibilities: avoid the per sample normalization (-o in the command line) and lower the LDA threshold, or attach the OTU level to all taxa.

I hope this helps,

thanks

Nicola