...

For each beta diversity method, both overall significance and pairwise significance can be calculated using a Permutational Multivariate Analysis of Variance (PERMANOVA), a non-parametric multivariate statistical test. The R-squared value represents the percentage of variance explained by the examined groups. E.g. if R = 0.23 then 23% of the total diversity is explained by groupwise differences. PERMANOVA is based on groupwise differences, thus cannot be applied to continuous data.

...

5a. Choosing a variable to analyse

In the Alpha Diversity section, you already imported your nfcore/ampliseq results and converted them to an ampvis2 object. To have another quick look at your ampvis object:

...

group <- "Nose_size"

Ordering your variable

The plotting done in ampvis2 is done by the ggplot2 package. ggplot factorises variables and automatically orders them on the plot by alphabetical order.

This can cause your groups to be ordered incorrectly on the plot axes (e.g. a time series may not be plotted sequentially).

You can manually set the order of your variable here.

First have a look at how ggplot will order your variable.

...

And once again you’ll need to order the groups in your variable. See the ‘Ordering your variable’ section in the Alpha diversity section for more details.

Choose how you want to order your groups:

lev <- c("Small", "Medium", "Big")

Then run the following to apply the levels to your data:

ampvisdata$metadata[[group]] <- factor(ampvisdata$metadata[[group]], levels = lev))

If these are in the order you want to see them on your plot axes, nothing needs to be done. If they are in the wrong order you need to order them manually by setting the levels.

Choose how you want to order your groups here:

lev <- c("Small", "Medium", "Big")

To order your variable you need to put all the variable levels into the lev = c(..). Make sure each level is in double quotes and separated by a comma.

Then run the following to apply the levels to your data:

ampvisdata$metadata[[group]] <- factor(ampvisdata$metadata[[group]], levels = lev)

5b. PCoA plots and statistics

The overview section described (with links and references) the ordination methods that can be used to estimate and plot beta diversity.

Briefly, these are: Bray–Curtis dissimilarity, Cao index, and Jaccard similarity index. Each of these has strengths and weaknesses. It's up to you, the researcher, to explore the literature and decide which is the best index to use for your data.

Choose the ordination method you want to use to estimate and plot beta diversity.

Bray–Curtis dissimilarity is used by default ("bray") Change this to "cao" for Cao index, or "jaccard" for Jaccard similarity index.

index <- "bray"

Now create the PCoA plot.

p <- amp_ordinate(ampvisdata, type = "pcoa", transform = "none", distmeasure = index, sample_color_by = group, sample_point_size = 3, sample_colorframe = TRUE) +
scale_color_manual(values=c("Red", "Green", "Blue")) +
scale_fill_manual(values=c("Red", "Green", "Blue")) +
theme(text = element_text(size = 16))
p

Save your plot as a 300dpi (i.e. publication quality) tiff or pdf file:

tiff_exp <- paste0("PCoA_beta_div_", group, "", index, "_samples.tif")
ggsave(file = tiff_exp, dpi = 300, compression = "lzw", device = "tiff", plot = p, width = 20, height = 20, units = "cm")

pdf_exp <- paste0("PCoA_beta_div_", group, "_", index,, "_samples.pdf")
ggsave(file = pdf_exp, device = "pdf", plot = p, width = 20, height = 20, units = "cm")

Statistical analysis

To compare the overall differences between groups within your chosen variable, a PERMANOVA test can be performed and similarly a pairwise PERMANOVA test can be performed to compare differences between each group.

PERMANOVA:

# Need to remove rows (from ASV abundance table) with all 0 counts first
asvmatrix <- ampvisdata$abund
asvmatrix <- asvmatrix[rowSums(asvmatrix) > 0, ]
# Also need to transpose (samples need to be as rows, asv's as columns)
asvmatrix <- t(asvmatrix)
# Then generate pairwise distance matrix
sampdist <- vegdist(asvmatrix, method="bray")
# Use adonis function (vegan package: "Permutational Multivariate Analysis of Variance Using Distance Matrices") to run PERMANOVA on distances
pathotype.adonis <- adonis2(sampdist ~ get(group), data = samples_table)
# Output the r squared and p values as variables
r2 <- pathotype.adonis$R2[1]
pval <- pathotype.adonis$`Pr(>F)`[1]