<p>I am trying to understand metagenomeSeq, I will be very much thankful if you could help me to understand this:</p>

<p>I do have three treatment condition (T1, T2, and T3) and want to see how these treatments affect the microbial diversity in rice roots. I am interested in the pairwise comparison between them.</p>

<p>I used below code but I am not sure whether it is right or wrong, I do also have few queries regarding below code. I initially imported whole data into phyloseq and then filter the root data and convert into metagenomeSeq object.</p>

<pre><code>root_data &lt;- phyloseq_to_metagenomeSeq(root_M_filter)
</code></pre>

<p>DO I need to used normalize count for further processing? what is setting?</p>

<pre><code>normalizedMatrix &lt;- MRcounts(root_data, norm = TRUE)
settings = zigControl(maxit = 1, verbose = FALSE)
</code></pre>

<p>Is it correct way to make model matrix?</p>

<pre><code>mod = model.matrix(~root_data<span class="math-container">$Treatment)
colnames(mod) = levels(root_data$</span>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
</code></pre>

<p>What is design here?</p>

<pre><code>finalMod = slot(res, "fit")$design
</code></pre>

<p>It is correct option for contract.make, or should I also involved some
more options?</p>

<pre><code>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
</code></pre>

CC BY-SA 4.0
2019-05-22T10:13:12.560Z
You get warnings, not errors. The problem is you keep trying to load ph_mgs as a dataset from a package, while it is an object you create. However, this warnings do not prevent you to analyze the ph_mgs as you wish
llrs
CC BY-SA 4.0
2019-05-23T07:26:07.410Z
@DrYogesh If the problem is unrelated then I would recommend deleting this post. But you could also [edit] this one to include the relevant section (if your problem is highly related to this code).
llrs
CC BY-SA 4.0
2019-05-24T07:43:17.657Z
Please use code formatting not quotation formatting for the question.
llrs
<bioconductor>
Queries regarding MetagenomeSeq code
# Queries regarding MetagenomeSeq code

(Asked by: username_0 on 2019-05-18)

I am trying to understand metagenomeSeq, I will be very much thankful if you could help me to understand this:
I do have three treatment condition (T1, T2, and T3) and want to see how these treatments affect the microbial diversity in rice roots. I am interested in the pairwise comparison between them.
I used below code but I am not sure whether it is right or wrong, I do also have few queries regarding below code. I initially imported whole data into phyloseq and then filter the root data and convert into metagenomeSeq object.
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)

Is it correct way to make model matrix?
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

What is design here?
finalMod = slot(res, "fit")$design

It is correct option for contract.make, or should I also involved some
more options?
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

**Tags:** bioconductor
> **Comments:**
> **username_1**: You get warnings, not errors. The problem is you keep trying to load ph_mgs as a dataset from a package, while it is an object you create. However, this warnings do not prevent you to analyze the ph_mgs as you wish  
> **username_1**: @DrYogesh If the problem is unrelated then I would recommend deleting this post. But you could also [edit] this one to include the relevant section (if your problem is highly related to this code).  
> **username_1**: Please use code formatting not quotation formatting for the question.  

---

## 0 answers

