Небольшой пример обработки микробных данных на примере анализа измкенчивости микробиома литторин (данные изменены)

Осторожно, встречаются неоптимизированные действия, неединообразный код и конфликты пакетов.

1 Основные действия до фильтрации

1.1 Загрузка данных

## phyloseq-class experiment-level object
## otu_table()   OTU Table:         [ 7920 taxa and 80 samples ]
## sample_data() Sample Data:       [ 80 samples by 12 sample variables ]
## tax_table()   Taxonomy Table:    [ 7920 taxa by 8 taxonomic ranks ]
##  [1] "L1_L_fab_gut2"    "L1_L_nod_BF2"     "L1_L_obt_gut1"    "L1_L_obt_gut2"   
##  [5] "L1_L_obt_gut3"    "L1_L_obt_gut5"    "L1_L_obt_gut6"    "L1_L_obt_gut7"   
##  [9] "L1_L_obt_gut8"    "L1_L_ser_BF1"     "L1_L_ser_BF2"     "L1_L_stones_BF1" 
## [13] "L1_L_ves_BF1"     "L1_L_ves_BF2"     "L1_L_ves_BF3"     "L1_U_obt_gut1"   
## [17] "L1_U_obt_gut2"    "L1_U_obt_gut3"    "L1_U_obt_gut7"    "L11_L_fab_gut2"  
## [21] "L11_L_fab_gut3"   "L11_L_fab_gut4"   "L11_L_fab_gut5"   "L11_L_fab_gut6"  
## [25] "L11_L_fab_gut7"   "L11_L_fab_gut8"   "L11_L_fab_tent1"  "L11_L_fab_tent2" 
## [29] "L11_L_fill_BF1"   "L11_L_fill_BF2"   "L11_L_fill_BF3"   "L11_L_nod_BF2"   
## [33] "L11_L_nod_BF3"    "L11_L_obt_gut3"   "L11_L_obt_gut4"   "L11_L_obt_gut5"  
## [37] "L11_L_obt_gut6"   "L11_L_obt_gut7"   "L11_L_obt_gut8"   "L11_L_obt_tent1" 
## [41] "L11_L_obt_tent2"  "L11_L_stones_BF1" "L11_L_stones_BF2" "L11_L_stones_BF3"
## [45] "L11_U_fill_BF1"   "L11_U_nod_BF2"    "L11_U_obt_gut1"   "L11_U_obt_gut2"  
## [49] "L11_U_obt_gut3"   "L11_U_obt_gut4"   "L11_U_obt_gut5"   "L11_U_obt_gut8"  
## [53] "L11_U_obt_tent2"  "L11_U_stones_BF1" "L11_U_stones_BF3" "L5_L_fab_gut1"   
## [57] "L5_L_fab_gut3"    "L5_L_fab_gut5"    "L5_L_fab_gut6"    "L5_L_fab_gut7"   
## [61] "L5_L_fab_gut8"    "L5_L_nod_BF2"     "L5_L_obt_gut1"    "L5_L_obt_gut3"   
## [65] "L5_L_obt_gut4"    "L5_L_obt_gut7"    "L5_L_obt_gut8"    "L5_L_stones_BF2" 
## [69] "L5_L_ves_BF3"     "L5_U_fill_BF1"    "L5_U_fill_BF2"    "L5_U_nod_BF1"    
## [73] "L5_U_nod_BF3"     "L5_U_obt_gut1"    "L5_U_obt_gut2"    "L5_U_obt_gut4"   
## [77] "L5_U_obt_gut5"    "L5_U_obt_gut6"    "L5_U_obt_gut8"    "L5_U_obt_tent2"
## [1] "Domain"  "Phylum"  "Class"   "Order"   "Family"  "Genus"   "Species"
## [8] "X"
##  [1] "site"              "species"           "sp_type"          
##  [4] "subtype"           "type"              "sufficient"       
##  [7] "level"             "sea"               "season"           
## [10] "season_time"       "season_site"       "season_site_level"

2 Фильтрация и нормализация образцов

2.1 Общая фильтрация

Фильтрация мягкая - удаляем ASV, встреченные 1 раз во всем наборе данных (singletones), ASV, присутствующие менее чем в 5% образцов. Порог обилия моджет быть иным, но даже если удалить фичи, встреченные менее, чем в 10%, это может помешать найти различия между образцами по редким микробам.

##  L11_L_obt_gut5  L11_L_obt_gut7 L1_L_stones_BF1    L5_L_ves_BF3  L11_L_fill_BF2 
##               0               0            2715            4108            4621 
##   L1_L_obt_gut7  L11_L_obt_gut6   L5_U_obt_gut2   L5_L_obt_gut7  L11_L_obt_gut3 
##            5210            5366            5494            6242            6543 
##   L5_L_obt_gut8  L11_U_obt_gut1 L11_L_fab_tent2   L5_U_obt_gut5   L5_L_fab_gut6 
##            7368            8125            8395            8413            8469
## L11_L_stones_BF2    L1_L_obt_gut5    L1_L_fab_gut2    L1_L_obt_gut2 
##            34988            34010            29782            28611 
##    L1_U_obt_gut3    L1_U_obt_gut1    L1_L_obt_gut8    L1_U_obt_gut7 
##            26706            25122            23868            23568 
## L11_L_stones_BF3     L5_U_nod_BF3    L1_L_obt_gut1    L11_U_nod_BF2 
##            22115            21893            21624            21035 
##    L1_L_obt_gut3 L11_U_stones_BF3 L11_U_stones_BF1 
##            20357            19611            19070

2.2 Удаление низкой глубины прочтения

первые два образца (с наиболее низкой глубиной прочтения, она НОЛЬ) далее в анализе не используются

кроме того, имеет смысл удалить лишние летние образцы (Лето 2021, неполные L.fabalis)

было принято волевое решение не удалять единственную пробу F.vesiculosus

##    L5_L_ves_BF3  L11_L_fill_BF2  L11_L_obt_gut6   L5_U_obt_gut2   L5_L_obt_gut7 
##            4108            4621            5366            5494            6242 
##  L11_L_obt_gut3   L5_L_obt_gut8  L11_U_obt_gut1 L11_L_fab_tent2   L5_U_obt_gut5 
##            6543            7368            8125            8395            8413 
##   L5_L_fab_gut6   L5_L_obt_gut1   L11_L_nod_BF2   L5_L_obt_gut3  L11_L_obt_gut4 
##            8469            8531            8536            8813            8861

Значения размеров библиотек для самых больших образцов если что такие:

##    L5_U_obt_gut1 L11_L_stones_BF1   L11_U_obt_gut3   L11_U_obt_gut5 
##            15521            15837            17379            17694 
##   L11_U_obt_gut8   L11_L_fill_BF3   L11_U_obt_gut4     L5_L_nod_BF2 
##            17959            18104            18279            18738 
##   L11_L_fill_BF1 L11_U_stones_BF1 L11_U_stones_BF3    L11_U_nod_BF2 
##            18896            19070            19611            21035 
##     L5_U_nod_BF3 L11_L_stones_BF3 L11_L_stones_BF2 
##            21893            22115            34988

2.3 Нормализация

Мы в своем анализе делили на медианный размер библиотеки (здесь - 11539) для всех образцов, но можно хоть логарифмировать

3 Интегральные характеристики

Давайте сначала опишем простые вещи

3.1 Таблица

альфа-ранообразие, полный датасет, объединено по типу образца

3.2 Выровненность и Шеннон общие

3.3 Ранее построенный Шеннон

3.3.0.1 Шеннон в кишках, по сезонам на Белом

3.3.0.2 Выровненность в кишках, по сезонам на Белом

4 Ординации и PERMANOVA

Наконец-то многомерный подход

Базово я привыкла делать nMDS по матрице различий Брея-Куртиса, но с некоторыми (напр, логарифмическими нормализациями он не дружит) Можно делать и Жаккара в качестве расстояний, и хоть главные компоненты в качестве ординации Можно делать хитмап или кластеризацию, главное понять, на какие расстояния относит композиция микробиома ваши образцы друг ото друга

## Square root transformation
## Wisconsin double standardization
## Run 0 stress 0.200135 
## Run 1 stress 0.2243491 
## Run 2 stress 0.2147016 
## Run 3 stress 0.2089865 
## Run 4 stress 0.2105712 
## Run 5 stress 0.205666 
## Run 6 stress 0.2051715 
## Run 7 stress 0.223335 
## Run 8 stress 0.2075067 
## Run 9 stress 0.2032072 
## Run 10 stress 0.2146564 
## Run 11 stress 0.2043048 
## Run 12 stress 0.2136925 
## Run 13 stress 0.2033042 
## Run 14 stress 0.2151472 
## Run 15 stress 0.2036829 
## Run 16 stress 0.1987382 
## ... New best solution
## ... Procrustes: rmse 0.02500797  max resid 0.165155 
## Run 17 stress 0.2072265 
## Run 18 stress 0.204302 
## Run 19 stress 0.1969736 
## ... New best solution
## ... Procrustes: rmse 0.03306505  max resid 0.2160054 
## Run 20 stress 0.2051777 
## *** Best solution was not repeated -- monoMDS stopping criteria:
##     20: stress ratio > sratmax
## [1] 0.1969736

4.1 Картинки

4.1.1 Все точки на одной панели автоматом

4.1.2 Все точки на одной панели руками

4.1.3 Фасетки по сезонам

4.1.4 Кишечная ординация

4.1.5 Для одного вида

4.2 Циферки

будем тестировать с помощью PERMANOVA, есть ли влияние вида для кишечного микробиома и есть ли влияние сезона

4.2.1 кишки, 2 сезона, 2 вида, 2 фактора

## [1] FALSE
## 
## May-June_2022 November_2021 
##            17            17
## 
##  L_fabalis L_obtusata 
##         13         21

4.2.2 кишки летнего сезона, 2 вида, однофакторно

4.2.3 кишки летнего сезона, 2 вида, однофакторно, иная перманова

## 
## Analysis of Variance, using Residual Randomization
## Permutation procedure: Randomization of null model residuals 
## Number of permutations: 1000 
## Estimation method: Ordinary Least Squares 
## Sums of Squares and Cross-products: Type I 
## Effect sizes (Z) based on F distributions
## 
##           Df     SS      MS     Rsq      F      Z Pr(>F)  
## species    1 0.3341 0.33409 0.09513 1.5769 1.3603  0.089 .
## Residuals 15 3.1779 0.21186 0.90487                       
## Total     16 3.5120                                       
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Call: lm.rrpp(f1 = guts_white_bray_snails ~ species, data = guts_white_df)

5 Барплоты обилия

5.1 для всего датасета

5.1.1 встроенное в филосек

5.1.2 специальный пакет для вложенных столбиков

5.2 Барплоты для отдельных видов

5.2.1 Pseudoalteromonas

5.2.2 Streptococcus

6 Бонусный пример хитмапа