Exercise 4.2: Exploring the armpit dataset - solution

In this exercise, we will explore the “smelly armpit” dataset.

1 Smelly armpit dataset

Smelly armpits are not caused by sweat, itself. The smell is caused by specific micro-organisms belonging to the group of Corynebacterium spp. that metabolise sweat. Another group of abundant bacteria are the Staphylococcus spp., these bacteria do not metabolise sweat in smelly compounds.

The CMET-group at Ghent University does research to on transplanting the armpit microbiome to save people with smelly armpits.

• Proposed Therapy:

  1. Remove armpit-microbiome with antibiotics
  2. Influence armpit microbiome with microbial transplant, see this 2 minute talk on youtube

• Experiment:

  • 20 students with smelly armpits are attributed to one of two treatment groups
  • placebo (only antibiotics)
  • transplant (antibiotica followed by microbial transplant).
  • The microbiome is sampled 6 weeks upon the treatment
  • The relative abundance of Staphylococcus spp. on Corynebacterium spp. + Staphylococcus spp. in the microbiome is measured via DGGE (Denaturing Gradient Gel Electrophoresis).

Load the libraries

library(tidyverse)

2 Import the dataset

ap <- read_csv("https://raw.githubusercontent.com/statOmics/PSLS21/data/armpit.csv")

## Rows: 20 Columns: 2

## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): trt
## dbl (1): rel

## 
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.

glimpse(ap)

## Rows: 20
## Columns: 2
## $ trt <chr> "placebo", "placebo", "placebo", "placebo", "placebo", "placebo", …
## $ rel <dbl> 54.99208, 31.84466, 41.09948, 59.52064, 63.57341, 41.48649, 30.440…

3 Goal

The overarching goal of this study was to assess if the relative abundance Staphylococcus spp. in the microbiome of the armpit is affected by transplanting the microbiome. To this end, the researchers randomized patients to two treatments: A treatment with antibiotics only and a treatment with antibiotics and a microbial transplant.

In the tutorial on hypotheses testing we will use a formal statistical test to generalize the results from the sample to that of the population. For this test to be valid, we have to assess following assumptions:

1. The data in each treatment group is normally distributed.

2. The data from the two treatment groups has the same variance.

A statistical analysis always starts with a exploratory data analysis to get insight in the nature and distribution of the data, and to assess the assumptions of the downstream data analysis. Mastering this data exploration step is the purpose of the current exercise

4 Data visualization

A crucial first step in a data analysis is to visualize and to explore the raw data.

First, try to make a histogram of the data. Fill in the missing parts in the chunk of code below to get a good-looking visualization:

ap %>%
  ggplot(aes(x=rel,fill=trt)) + # fill in the correct values for x and fill 
  geom_histogram() +
  facet_grid(rows = vars(trt)) + # fill in to put the histograms for both treatment conditions in a separate row
  theme_bw() +
  xlab("relative abundance (%)")

## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.

Based on this plot, it seems that the relative abundance is higher for subjects who had the transplants. However, given the small sample size the histograms are not optimally informative. A better option for these data would be to show the data in a boxplot:

ap %>%  ggplot(aes(x=trt,y=rel,fill=trt)) +
  geom_boxplot(outlier.shape=NA) + 
  geom_point(position="jitter") +
  ylab("relative abundance (%)") +
  xlab("treatment group") + 
  stat_summary(fun.y=mean, geom="point", shape=5, size=3, color="black", fill="black")

## Warning: `fun.y` is deprecated. Use `fun` instead.

What do you observe?

5 Descriptive statistics

Here, we will generate some informative descriptive statistics for the dataset.

We first summarize the data and calculate the mean, standard deviation, number of observations and standard error and store the result in an object apRelSum via ’apRelSum<-`

1. We pipe the ap dataframe to the group_by function to group the data by treatment trt group_by(trt)
2. We pipe the result to the summarize_at function to summarize the “rel” variable and calculate the mean, standard deviation and the number of observations
3. We pipe the result to the mutate function to make a new variable in the data frame that is named se for which we calculate the standard error

apRelSum<-ap%>%
  group_by(trt)%>%
  summarize_at("rel",
               list(mean=~mean(.,na.rm=TRUE),
                    sd=~sd(.,na.rm=TRUE),
                    n=function(x) x%>%is.na%>%`!`%>%sum)) %>%
  mutate(se=sd/sqrt(n))

apRelSum

This concludes the data exploration. Tomorrow, we will learn how to formally test if the observed difference is statistically significant.

---

This work is licensed under the CC BY-NC-SA 4.0 licence.