Статистика и R

Занятие 5

Статистика

Применение t-теста
Проверка на нормальность
Непараметрические тесты

R

Добавление текста на графики
Использование пакета ggpubr
Длинный и широкий формат
Объединение датафреймов

Широкий формат

Несколько переменных, похожих по своей природе, в разных столбцах

df_wide

# A tibble: 5 × 4
  sample  gene1 gene2 gene3
  <chr>   <dbl> <dbl> <dbl>
1 sample1    10     5     1
2 sample2     5     3     0
3 sample3     7     4    10
4 sample4     8     6     9
5 sample5    11     4     7

Длинный формат

Столбец с названиями переменных
Столбец со значениями переменных

df_long

# A tibble: 15 × 3
   sample  gene  expression
   <chr>   <chr>      <dbl>
 1 sample1 gene1         10
 2 sample1 gene2          5
 3 sample1 gene3          1
 4 sample2 gene1          5
 5 sample2 gene2          3
 6 sample2 gene3          0
 7 sample3 gene1          7
 8 sample3 gene2          4
 9 sample3 gene3         10
10 sample4 gene1          8
11 sample4 gene2          6
12 sample4 gene3          9
13 sample5 gene1         11
14 sample5 gene2          4
15 sample5 gene3          7

df_long %>% 
  ggplot(aes(
    x = gene, 
    y = expression, 
    fill = gene)) +
  geom_boxplot()

Длинный и широкий формат

Из широкого в длинный - pivot_longer()
Из длинного в широкий - pivot_wider()

`pivot_longer()`

Из широкого в длинный формат

df_wide %>% 
  pivot_longer(
    cols = gene1:gene3,        # в каких столбцах нужные переменные
    names_to = "gene",         # название столбца с именами переменных
    values_to = "expression")  # название столбца со значениями переменных

# A tibble: 15 × 3
   sample  gene  expression
   <chr>   <chr>      <dbl>
 1 sample1 gene1         10
 2 sample1 gene2          5
 3 sample1 gene3          1
 4 sample2 gene1          5
 5 sample2 gene2          3
 6 sample2 gene3          0
 7 sample3 gene1          7
 8 sample3 gene2          4
 9 sample3 gene3         10
10 sample4 gene1          8
11 sample4 gene2          6
12 sample4 gene3          9
13 sample5 gene1         11
14 sample5 gene2          4
15 sample5 gene3          7

`pivot_wider()`

df_long %>%
  pivot_wider(
    names_from = gene,        # в каком столбце имена переменных
    values_from = expression) # в каком столбце значения переменных

# A tibble: 5 × 4
  sample  gene1 gene2 gene3
  <chr>   <dbl> <dbl> <dbl>
1 sample1    10     5     1
2 sample2     5     3     0
3 sample3     7     4    10
4 sample4     8     6     9
5 sample5    11     4     7

Объединение датафреймов

Объедиение проискходит по одному или нескольким столбцам с “ключами”.

full_join() - сохранить все строчки из обеих таблиц
inner_join() - сохранить только строчки с общими “ключами”
left_join() - сохранить все строчки из “первой” таблицы
right_join() - сохранить все строчки из “второй” таблицы

Нужно объединить таблицу с результатом эксперимента и таблицу с метаданными:

df_meta

# A tibble: 6 × 2
  sample  source   
  <chr>   <chr>    
1 sample1 tail     
2 sample2 left_ear 
3 sample3 tail     
4 sample4 right_ear
5 sample5 left_ear 
6 sample6 right_ear

df_joined <- full_join(df_wide, df_meta, by = "sample")

df_joined

# A tibble: 6 × 5
  sample  gene1 gene2 gene3 source   
  <chr>   <dbl> <dbl> <dbl> <chr>    
1 sample1    10     5     1 tail     
2 sample2     5     3     0 left_ear 
3 sample3     7     4    10 tail     
4 sample4     8     6     9 right_ear
5 sample5    11     4     7 left_ear 
6 sample6    NA    NA    NA right_ear

df_joined %>% 
  pivot_longer(cols = gene1:gene3, names_to = "gene", values_to = "expression") %>% 
  ggplot(aes(
    x = gene, 
    y = expression, 
    fill = source)) +
  geom_boxplot()

Пример - данные про пингвинов

library(tidyverse)

penguins <- read_csv("https://rstats-at-bio-msu.netlify.app/data/penguins.csv")
penguins

# A tibble: 344 × 8
   species island    bill_length_mm bill_depth_mm flipper_…¹ body_…² sex    year
   <chr>   <chr>              <dbl>         <dbl>      <dbl>   <dbl> <chr> <dbl>
 1 Adelie  Torgersen           39.1          18.7        181    3750 male   2007
 2 Adelie  Torgersen           39.5          17.4        186    3800 fema…  2007
 3 Adelie  Torgersen           40.3          18          195    3250 fema…  2007
 4 Adelie  Torgersen           NA            NA           NA      NA <NA>   2007
 5 Adelie  Torgersen           36.7          19.3        193    3450 fema…  2007
 6 Adelie  Torgersen           39.3          20.6        190    3650 male   2007
 7 Adelie  Torgersen           38.9          17.8        181    3625 fema…  2007
 8 Adelie  Torgersen           39.2          19.6        195    4675 male   2007
 9 Adelie  Torgersen           34.1          18.1        193    3475 <NA>   2007
10 Adelie  Torgersen           42            20.2        190    4250 <NA>   2007
# … with 334 more rows, and abbreviated variable names ¹flipper_length_mm,
#   ²body_mass_g

Действительно ли длины крыльев у пингвинов Adelie и Chinstrap различаются?

# Вектор с длинами крыльев Adelie
adelie <- penguins %>% 
  filter(species == "Adelie") %>% 
  pull(flipper_length_mm)

# Вектор с длинами крыльев Chinstrap
chinstrap <- penguins %>% 
  filter(species == "Chinstrap") %>% 
  pull(flipper_length_mm)

t.test(adelie, chinstrap)


    Welch Two Sample t-test

data:  adelie and chinstrap
t = -5.7804, df = 119.68, p-value = 6.049e-08
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -7.880530 -3.859244
sample estimates:
mean of x mean of y 
 189.9536  195.8235

Добавить текст на график

Любой текстовый элемент можно добавить на график с помощью функции annotate()

# Получили p-значение
pval <- t.test(adelie, chinstrap)$p.value
pval <- format(pval, digits = 3, scientific = T)

penguins %>% filter(species %in% c("Adelie", "Chinstrap")) %>% 
  ggplot(aes(x = species, y = flipper_length_mm, fill = species)) +
  geom_boxplot() + theme(legend.position = "none") +
  annotate("text", x = 1, y = 220, label = paste("T-test, p-val =", pval), size = 4)# Получили p-значение
pval <- t.test(adelie, chinstrap)$p.value
pval <- format(pval, digits = 3, scientific = T)

penguins %>% filter(species %in% c("Adelie", "Chinstrap")) %>% 
  ggplot(aes(x = species, y = flipper_length_mm, fill = species)) +
  geom_boxplot() + theme(legend.position = "none") +
  annotate("text", x = 1, y = 220, label = paste("T-test, p-val =", pval), size = 4)

Добавить результат теста на график

Функция stat_compare_means() из пакета ggpubr проводит тест и выводит реузльтат на график

library(ggpubr)

penguins %>% filter(species %in% c("Adelie", "Chinstrap")) %>% 
  ggplot(aes(x = species, y = flipper_length_mm, fill = species)) +
  geom_boxplot() + theme(legend.position = "none") +
  stat_compare_means(method = "t.test")library(ggpubr)

penguins %>% filter(species %in% c("Adelie", "Chinstrap")) %>% 
  ggplot(aes(x = species, y = flipper_length_mm, fill = species)) +
  geom_boxplot() + theme(legend.position = "none") +
  stat_compare_means(method = "t.test")

Добавить текст с помощью `geom_text()`

penguin_species <- penguins %>% 
  count(species)

penguin_species

# A tibble: 3 × 2
  species       n
  <chr>     <int>
1 Adelie      152
2 Chinstrap    68
3 Gentoo      124

ggplot(penguin_species, 
  aes(
    x = species, 
    y = n,
    fill = species)) +
  geom_col() +
  geom_text(aes(label = n, y = n + 5)) +
  theme(legend.position = "none")ggplot(penguin_species, 
  aes(
    x = species, 
    y = n,
    fill = species)) +
  geom_col() +
  geom_text(aes(label = n, y = n + 5)) +
  theme(legend.position = "none")

Добавить текст близко к координате

Возьмем небольшую подвыборку:

set.seed(42)

ten_penguins <- penguins %>% 
  slice_sample(n = 10) # случайные 10 строк из датафрейма

ggplot(ten_penguins, aes(x = body_mass_g, y = flipper_length_mm)) +
  geom_point() +
  geom_text(aes(label = species))

Добавить текст близко к координате

С помощью функции geom_text_repel() из пакета ggrepel можно добавить надпись к точке

library(ggrepel)

ggplot(ten_penguins, aes(x = body_mass_g, y = flipper_length_mm)) +
  geom_point() +
  geom_text_repel(aes(label = species))library(ggrepel)

ggplot(ten_penguins, aes(x = body_mass_g, y = flipper_length_mm)) +
  geom_point() +
  geom_text_repel(aes(label = species))