-
-
Save claczny/3415270a6c6919969bff79d2c246a527 to your computer and use it in GitHub Desktop.
Call: | |
adonis(formula = GPfr_phylum_bray ~ SampleType, data = sampledf) | |
Permutation: free | |
Number of permutations: 999 | |
Terms added sequentially (first to last) | |
Df SumsOfSqs MeanSqs F.Model R2 Pr(>F) | |
SampleType 8 8.2618 1.03272 5.81 0.7322 0.001 *** | |
Residuals 17 3.0218 0.17775 0.2678 | |
Total 25 11.2836 1.0000 | |
--- | |
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 | |
Call: | |
adonis(formula = GPfr_phylum_bray ~ Description, data = sampledf) | |
Permutation: free | |
Number of permutations: 999 | |
Terms added sequentially (first to last) | |
Df SumsOfSqs MeanSqs F.Model R2 Pr(>F) | |
Description 24 11.1949 0.46646 5.2634 0.99215 0.027 * | |
Residuals 1 0.0886 0.08862 0.00785 | |
Total 25 11.2836 1.00000 | |
--- | |
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 | |
Call: | |
adonis(formula = GPfr_phylum_bray ~ SampleType + Description, data = sampledf) | |
Permutation: free | |
Number of permutations: 999 | |
Terms added sequentially (first to last) | |
Df SumsOfSqs MeanSqs F.Model R2 Pr(>F) | |
SampleType 8 8.2618 1.03272 11.6530 0.73220 0.022 * | |
Description 16 2.9331 0.18332 2.0686 0.25995 0.060 . | |
Residuals 1 0.0886 0.08862 0.00785 | |
Total 25 11.2836 1.00000 | |
--- | |
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 |
# Adjusted from http://deneflab.github.io/MicrobeMiseq/demos/mothur_2_phyloseq.html#permanova | |
library(phyloseq) | |
library(dplyr) | |
library(ggplot2) | |
library(vegan) | |
data(GlobalPatterns) | |
GPr = transform_sample_counts(GlobalPatterns, function(x) x / sum(x) ) | |
GPfr = filter_taxa(GPr, function(x) mean(x) > 1e-5, TRUE) | |
GPfr_phylum <- GPfr %>% | |
tax_glom(taxrank = "Phylum") %>% # agglomerate at phylum level | |
transform_sample_counts(function(x) {x/sum(x)} ) %>% # Transform to rel. abundance | |
psmelt() %>% # Melt to long format | |
filter(Abundance > 0.02) %>% # Filter out low abundance taxa | |
arrange(Phylum) | |
### | |
# Plot ##### | |
### | |
p.phyla <- ggplot(GPfr_phylum, aes(x = Sample, y = Abundance, fill = Phylum)) + | |
geom_bar(stat = "identity") + | |
# Remove x axis title | |
theme(axis.title.x = element_blank()) + | |
# | |
guides(fill = guide_legend(reverse = TRUE, keywidth = 1, keyheight = 1)) + | |
ylab("Relative Abundance (Phyla > 2%) \n") + | |
ggtitle("Phylum Composition of GlobalPatterns Communities by Sample") | |
print(p.phyla) | |
### | |
# PERMANOVA ##### | |
### | |
set.seed(1) | |
# Calculate bray curtis distance matrix | |
GPfr_phylum_bray <- phyloseq::distance(GPfr, method = "bray") | |
# make a data frame from the sample_data | |
sampledf <- data.frame(sample_data(GPfr)) | |
# Adonis test | |
adonis(GPfr_phylum_bray ~ SampleType, data = sampledf) | |
adonis(GPfr_phylum_bray ~ Description, data = sampledf) | |
adonis(GPfr_phylum_bray ~ SampleType + Description, data = sampledf) |
Note that the order of the covariates (a.k.a., independent variables) MATTERS in this case, i.e., different results are obtained using SampleType + Description
or Description + SampleType
(results not shown). For more info, s. https://stats.stackexchange.com/questions/188519/adonis-in-vegan-order-of-variables-or-use-of-strata
[UPDATE] This is for ANOVA http://www.statmethods.net/stats/anova.html says:
In a nonorthogonal design with more than one term on the right hand side of the equation order will matter (i.e., A+B and B+A will produce different results)! We will need use the
drop1( )
function to produce the familiar Type III results.
[UPDATE2] For the adonis
function, and maybe related to aov
, the testing is done sequentially, s. vegandevs/vegan#229. Hence, the effect of the ordering. Unsure what the "marginal" means, instead of "sequential".
[Update3] Note the comment in vegandevs/vegan#229:
I do not recommend using all possible orders. Either you do marginal tests and accept the statistics it can deliver, or you have an a priori meaningful order.
-> See the
a priori meaningful order
Not sure how to get this though.
Note, as described in http://deneflab.github.io/MicrobeMiseq/demos/mothur_2_phyloseq.html#permanova, the data should be scaled to have similar sequencing depth. This is omitted in this example, but should be performed when using your own data.