library("ampvis")
data(DNAext_1.0)
pcr <- subset_samples(V13, Exp.PCR == "YES") %>%
rarefy_even_depth(sample.size = 25000, rngseed = 712) %>%
filter_taxa(function(x) max(x) >= 10, TRUE)
PCA with square root transformed OTU abundances.
pca <- amp_ordinate(data = pcr,
plot.color = "Add.Label",
plot.point.size = 5,
plot.group = "chull",
plot.group.label = "Add.Label",
output = "complete",
envfit.factor = "Add.Label",
envfit.show = F,
plot.theme = "clean")
Plot the PCA. It looks like there might be some significant grouping.
pca$plot + theme(legend.position = "none")
The model reports a p-value of 0.001, hence there is an effect of different PCR settings.
pca$eff.model
##
## ***FACTORS:
##
## Centroids:
## PC1 PC2
## Add.Label1 ng -2.2568 1.6732
## Add.Label25 cycles -2.1403 0.6196
## Add.Label35 cycles 0.0822 0.0188
## Add.Label50 ng 2.8735 1.1218
## Add.Label52 C 0.0454 -4.4422
## Add.Label58 C 2.3664 1.2207
## Add.LabelStandard -0.9431 -0.2056
##
## Goodness of fit:
## r2 Pr(>r)
## Add.Label 0.8995 0.001 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## Permutation: free
## Number of permutations: 999
The Bray-Curtis dissimilarity index is used as an alternative method to test for significant groupings in the dataset.
beta <- amp_test_cluster(data = pcr, group = "Add.Label", method = "bray", plot.color = "Add.Label", plot.label = "Add.Label")
Using adonis we also find a significant effect as the p-value is 0.001.
beta$adonis
##
## Call:
## adonis(formula = test_formula, data = sample)
##
## Permutation: free
## Number of permutations: 999
##
## Terms added sequentially (first to last)
##
## Df SumsOfSqs MeanSqs F.Model R2 Pr(>F)
## Add.Label 6 0.068761 0.0114601 2.217 0.50574 0.001 ***
## Residuals 13 0.067200 0.0051692 0.49426
## Total 19 0.135960 1.00000
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
amp_ordinate(data = pcr,
plot.color = "Add.Label",
plot.point.size = 2,
plot.group = "chull",
plot.theme = "clean"
) +
xlim(-4,4) +
annotate("text", x = -2, y = 3, label = "1 ng", size = 2) +
annotate("text", x = -1, y = -0.9, label = "Standard", size = 2) +
annotate("text", x = -1, y = -1.4, label = "56*' '*degree*C*', 5 ng, 30 cyc'", size = 2, parse = T) +
annotate("text", x = -3.4, y = 0.5, label = "25 cyc", size = 2) +
annotate("text", x = 0.8, y = -4, label = "52*' '*degree*C", size = 2, parse = T) +
annotate("text", x = 0.9, y = 0.3, label = "35 cyc", size = 2) +
annotate("text", x = 2.7, y = 0.3, label = "50 ng", size = 2) +
annotate("text", x = 2, y = 1.8, label = "58*' '*degree*C", size = 2, parse = T) +
theme(legend.position = "none",
axis.text = element_text(size = 6, color = "black"),
text = element_text(size = 8, color = "black")
)
ggsave("plots/Fig6A.eps", width = 50, height = 50, units = "mm")