Example using PERMANOVA with single and multiple covariates (a.k.a., independent variables).

adonis_output.txt

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

permanova_demo.R

# 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, 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.

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.