Need help in applying LRT to my count abundance data
0
0
Entering edit mode
14 months ago
rishav513 ▴ 30

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 • 446 views
ADD COMMENT
0
Entering edit mode

Someone kindly provide valuable suggestions

ADD REPLY

Login before adding your answer.

Traffic: 2594 users visited in the last hour
Help About
FAQ
Access RSS
API
Stats

Use of this site constitutes acceptance of our User Agreement and Privacy Policy.

Powered by the version 2.3.6