MetagenomeSeq for amplicon differential abundance analysis
0
0
Entering edit mode
4.9 years ago
Bioinfonext ▴ 460

Hi

I have imported phyloseq object into metagenomeSeq, I do have three treatment group:T1, T2 and T3. I want to see how these treatment group affect the microbial diversity in plant roots.

I want to do pairwise comparison, T1_VS_T2 , T2_VS_T3 and T1_VS_T3. I have tried below code to extract T1_vS_T2, IS IT CORRECT?

here I use the code but need some correction:

> root_data <- phyloseq_to_metagenomeSeq(root_M_filter)

**DO I need to used normalize count for further processing? what is setting? 

> normalizedMatrix <- MRcounts(root_data, norm = TRUE)


> settings = zigControl(maxit = 1, verbose = FALSE)



> mod = model.matrix(~root_data$Treatment)

> colnames(mod) = levels(root_data$Treatment)

> res = fitZig(obj = root_data, mod = mod, control = settings)
> res
An object of class "fitZigResults"
Slot "fit":
An object of class "MArrayLM"
$coefficients
                                         T1           T2          T3
71f84e7910006f22684121564206e8ca -2.8129112 -0.321366092  0.06113294
54e204fb99c80764e964456dadd6a0e5  0.5643061  0.019018800  1.01351768
55cd1fc570879d645bbf7a3642e9b0a8 -3.1159264  0.072348076  0.02205983
65f5c31e12c277aec319e2096463f9d2 -0.1355134 -0.004567206  0.53119153
7bed62f0fef250fd831dcf13bf43f4fc  1.9469670 -0.596315540 -0.35124338
                                 scalingFactor
71f84e7910006f22684121564206e8ca    1.62101725
54e204fb99c80764e964456dadd6a0e5   -0.24292482
55cd1fc570879d645bbf7a3642e9b0a8    1.50462644
65f5c31e12c277aec319e2096463f9d2    0.05833637
7bed62f0fef250fd831dcf13bf43f4fc   -0.56247683

> finalMod = slot(res, "fit")$design
> contrast.matrix = makeContrasts(T1 - T2, levels = finalMod)
> fit2 = contrasts.fit(zigFit, contrast.matrix)
> fit2 = eBayes(fit2)
> topTable(fit2)
                                     logFC   AveExpr         t      P.Value
3147790f0d5a78316fb9dd64f53b9473 11.262064 14.091887 16.331949 3.482815e-15
d2ec9f3b77975c0f457e4b7413b217ff  9.466087 12.212476 10.606122 6.159311e-11
7b908379d1fd2e75d84fd9aaecfd3f77  6.907182  9.784971  9.163780 1.267984e-09
d5fdc1df6dfd31acce9a34b9d8a76590  5.404459  8.373729  7.035926 1.802316e-07
19ddf8e59673c884734cf6093c672da0  5.381863  8.449413  7.011910 1.912363e-07

Kind Regards

R bioconductor • 1.3k views
ADD COMMENT

Login before adding your answer.

Traffic: 3041 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