---
title: "Figure 6A PCR settings"
author: "Mads Albertsen"
date: "`r format(Sys.time(), '%d-%m-%Y')`"
output: html_document
---

## Load packages
```{r Load_packages, message=FALSE, warning=FALSE, results='hide'}
library("ampvis")
```

## Load data
```{r load_data}
data(DNAext_1.0)
```

## Subset to the relevant dataset
```{r subset, message=FALSE}
pcr <- subset_samples(V13, Exp.PCR == "YES") %>%
  rarefy_even_depth(sample.size = 25000, rngseed = 712) %>%
  filter_taxa(function(x) max(x) >= 10, TRUE)
```

## Overall differences between samples using PCA
PCA with square root transformed OTU abundances.
```{r PCA}
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.
```{r pca_plot, fig.align='center', fig.height=7, fig.width=8}
pca$plot + theme(legend.position = "none")
```

The model reports a p-value of `r pca$eff.model$factors$pvals`, hence there is an effect of different PCR settings.
```{r pca_model}
pca$eff.model
```

## Overall differences between samples using beta diversity
The Bray-Curtis dissimilarity index is used as an alternative method to test for significant groupings in the dataset. 
```{r beta}
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 `r beta$adonis$aov.tab$"Pr(>F)"[1]`.
```{r beta_adonis}
beta$adonis
```

## Figure 6A: Influence of PCR settings
```{r Fig6A, fig.align='center', fig.width=2, fig.height=2}
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")
        )
```

## Save the plot
```{r save, eval=FALSE}
ggsave("plots/Fig6A.eps", width = 50, height = 50, units = "mm")
```