waldronlab / lefser

R implementation of the LEfSe method
https://waldronlab.io/lefser/
38 stars 6 forks source link

appear to be constant with groups #32

Closed celinechen1885 closed 4 months ago

celinechen1885 commented 9 months ago

Hi!

I'm working with a shotgun metagenomics data, and we aim to analyze the functional profile. we have 224 samples in total, with 356 different pathways. below are the codes that I used metagenomic_lefse <- read.csv se <- SummarizedExperiment(assays = list(counts = as.matrix(metagenomic_lefse))) colData(se) <- DataFrame(metadata) se$group <- as.factor(se$group) se$DI <- as.factor(se$DI) se_relab <- relativeAb(se)

res <- lefser(se_relab, groupCol = "group")

I got an error message: Error in lda.default(x, grouping, ...) : variables 85 145 appear to be constant within groups

Could anyone help me with it? Thanks!

lwaldron commented 5 months ago

This seems to happen with normalization to column sums adding to 1 instead of 1 million as in the original implementation. The following code after normalization accomplishes this, and we will make this the default of the relativeAb function.

assay(rse) = assay(rse)*1e6 # will become the default of relativeAb
lwaldron commented 5 months ago

@Peacesandy append * 1e6 to the line https://github.com/waldronlab/lefser/blob/db8f36ad48b05d5749518650ae490234cdfa7e1d/R/utils.R#L28 .

The following code currently reproduces the problem, so use it to create a testthat unit test of the fix:

library(lefser)
data(zeller14)
zeller14sub <- zeller14[, zeller14$study_condition != "adenoma"]
zeller14ra <- relativeAb(zeller14sub)
lefser(zeller14ra, groupCol = "study_condition")
Peacesandy commented 5 months ago

Okay, I'll make a pull request as soon as I'm done

Peacesandy commented 5 months ago

I just made a pull request @lwaldron

https://github.com/waldronlab/lefser/pull/34

lwaldron commented 4 months ago

Fixed by 8864568581444d1dfbcb12f8ef3348e19968ea08 (v1.13.2) and this error no longer appears. Thank you for the report, @celinechen1885 . @asyakhl will you push to bioc-devel ASAP? @celinechen1885, the easiest way to access this fix for now is to install from GitHub (BiocManager::install("waldronlab/lefser")), or just multiply your relative abundance values by 1e6.

lwaldron commented 4 months ago

@asyakhl please push to bioc-devel and cherry-pick to bioc-release. This is an important bugfix.

LiNk-NY commented 4 months ago

@lwaldron I've updated both. Thanks!

lwaldron commented 4 months ago

Thanks Marcel!