This is a test version of Biostars. For the public version, visit https://www.biostars.org.
Need help in applying LRT to my count abundance data

Hello all,

I am little confused in using LRT of DESeq2

I have ASV count table of microbiome data on which i want to do differential abundance analysis among different groups.

Here is the overview of my sample metadata

enter image description here

So, first of all, I load the required libraries and phyloseq object containing count information

library(phyloseq)
library(DESeq2)
ps10

Then, i try to compare the differences of effect of soil origin on plant compartments vs effect of plant genotype on plant compartments

diagdds = phyloseq_to_deseq2(ps10, ~ Plant_compartments + Soil_origin + Plant_genotype + Plant_compartments:Soil_origin + Plant_compartments:Plant_genotype)
dds = DESeq(diagdds, test="LRT", reduced=~Soil_origin + Plant_compartments:Soil_origin, fitType="local")
resultsNames(dds)
 [1] "Intercept"                            
 [2] "Plant_compartments_N_vs_L"            
 [3] "Plant_compartments_R_vs_L"            
 [4] "Plant_compartments_S_vs_L"            
 [5] "Soil_origin_F_vs_C"                   
 [6] "Plant_genotype_G27_vs_G1"             
 [7] "Plant_genotype_G96_vs_G1"             
 [8] "Plant_compartmentsN.Soil_originF"     
 [9] "Plant_compartmentsR.Soil_originF"     
[10] "Plant_compartmentsS.Soil_originF"     
[11] "Plant_compartmentsN.Plant_genotypeG27"
[12] "Plant_compartmentsR.Plant_genotypeG27"
[13] "Plant_compartmentsS.Plant_genotypeG27"
[14] "Plant_compartmentsN.Plant_genotypeG96"
[15] "Plant_compartmentsR.Plant_genotypeG96"
[16] "Plant_compartmentsS.Plant_genotypeG96"

Afterwards when i tried to compare between different groups such as

 res = results(dds, name = "Plant_compartmentsN.Plant_genotypeG27")
 sigtab = res[which(res$padj < 0.1), ]
 sigtab = cbind(as(sigtab, "data.frame"), as(tax_table(ps10)[rownames(sigtab), ], "matrix"))
 dim(sigtab)


    res1 = results(dds, name = "Plant_compartmentsR.Plant_genotypeG27")
    sigtab1 = res1[which(res$padj < 0.1), ]
    sigtab1 = cbind(as(sigtab1, "data.frame"), as(tax_table(ps10)[rownames(sigtab1), ], "matrix"))
    dim(sigtab1)

I am getting the same abundance of microbes that is approx 1500 microbes over every comparison . I want to know, what I am doing wrong

Thank you in advance

lrt r deseq2

Someone kindly provide valuable suggestions

0 answers

No answers yet.

Log in to answer this question.