---
layout: post
title: "Testing for significance between groups using adonis from the vegan R package"
blurb: "Are the differences between groups significant?"
author: "PD Schloss"
date: 2021-03-18 11:45
comments: false
youtube: 1ETBgbXl-BM
---

<iframe style="margin: 0 auto;display:block;" width="560" height="315" src="https://www.youtube.com/embed/{{ page.youtube }}" frameborder="0" allow="accelerometer; autoplay; encrypted-media; gyroscope; picture-in-picture" allowfullscreen></iframe>

## Code

This is where we started before the episode

```R
library(tidyverse)
library(readxl)
library(ggtext)

nmds <- read_tsv(file="raw_data/schubert.braycurtis.nmds_2d.axes")
metadata <- read_excel(path="raw_data/schubert.metadata.xlsx")
metadata_nmds <- inner_join(metadata, nmds, by=c('sample_id'='group')) %>%
  mutate(disease_stat = factor(disease_stat,
                               levels=c("DiarrhealControl",
                                        "Case",
                                        "NonDiarrhealControl")
  )
)

my_legend <- tibble(x = c(-0.75, 0.55, -0.85),
                    y = c(-0.55, 0.8, 0.7),
                    label = c("<strong>Healthy</strong>", "<strong>Diarrhea</strong>", "<strong>*C. difficile*<br>positive<strong>"),
                    color=c("NonDiarrhealControl", "DiarrhealControl", "Case")
)

ggplot(metadata_nmds,
       aes(x=axis1, y=axis2, color=disease_stat, fill=disease_stat)) +
  stat_ellipse(geom="polygon",type="norm", level=0.75, alpha=0.2, show.legend=F) +
  geom_point(show.legend=FALSE) +
  geom_richtext(data=my_legend,
                aes(x=x, y=y, label=label, color=color),
                hjust=0, lineheight=0.75,
                inherit.aes = FALSE, show.legend=FALSE,
                fill = NA, label.color = NA
  ) +
  coord_cartesian(xlim=c(-0.8, 0.8), ylim=c(-0.8, 0.8)) +
  labs(title="<span style='color:#999999'>**Healthy individuals**</span> have a different microbiota<br>from those with <span style='color:#0000FF'>**diarrhea**</span> and those with diarrhea<br>who are <span style='color: #ff0000'><strong>positive for *C. difficile*<strong></span>",
       x="NMDS Axis 1",
       y="NMDS Axis 2") +
  scale_color_manual(name=NULL,
                     breaks=c("NonDiarrhealControl", "DiarrhealControl", "Case"),
                     values=c("gray", "blue", "red"),
                     labels=c("Healthy", "Diarrhea",
                              "*C. difficile* positive")
  )+
  scale_fill_manual(name=NULL,
                    breaks=c("NonDiarrhealControl", "DiarrhealControl", "Case"),
                    values=c("lightgray", "dodgerblue", "pink"),
                    labels=c("Healthy", "Diarrhea",
                             "*C. difficile* positive")
  )+
  theme_classic() +
  theme(
    axis.title.y = element_text(hjust=1),
    axis.title.x = element_text(hjust=0),
    legend.key.size = unit(0.25, "cm"),
    legend.position = c(1.05, 0.95),
    legend.background = element_rect(fill="white", color="black"),
    legend.margin = margin(t=-2, r=3, b=3, l=3),
    legend.text = element_markdown(),
    plot.title.position = "plot",
    plot.title = element_markdown(margin=margin(b=1, unit="lines")),
    plot.margin = margin(l=0.5, r=0.5, t=0.5, b=0.5, unit="lines")
  )

ggsave("schubert_nmds.tiff", width=4.5, height=4.5)
```

## Installations

* [See video](https://www.youtube.com/watch?v=D6CunpqF04E) for instructions on getting set up
* [R](https://r-project.org)
* [RStudio](https://rstudio.com)
* [Raw data](https://github.com/riffomonas/raw_data/releases/latest)
