廿TT

譬如水怙牛過窓櫺 頭角四蹄都過了 因甚麼尾巴過不得

腸内細菌叢のデータで遊ぼう(Kostic et al, 2015)

今日の川柳

Kostic et al, 2015 (https://www.cell.com/cell-host-microbe/fulltext/S1931-3128(15)00021-9) では、乳幼児の腸内細菌叢のコホート研究をおこなっています。

その結果として、腸内細菌の種の多様性は年齢とともに増加する傾向があるが、糖尿病を発症する群は2歳くらいで種の多様性の増加がサチるという現象を発見しました。

f:id:abrahamcow:20160623203028j:plain
https://www.cell.com/cell-host-microbe/fulltext/S1931-3128(15)00021-9 より)

データセット
DIABIMMUNE
に公開されています。

本当は統計モデリングでこの説を検証したかったけど、いまのところノーアイデアなので、今回は可視化だけ。

Kostic らは種の多様性の指標として Chao 1 (https://www.ism.ac.jp/editsec/toukei/pdf/60-2-263.pdf)というのを使っているけど、Chao 1がなんでこの式になるのかよくわからなかったので、ぼくはエントロピーを使うことにします。

腸内細菌のでデータはDNAを便のサンプルからまとめて回収してきて、データベースに照会することで菌種を特定します。

このときDNAのひとつの断片が菌種 k に該当する確率を p_k と表します。

エントロピー

 - \sum_{k}p_k \log p_k

で定義されます。

このエントロピーをサンプルごとに計算してプロットしたのがこちらです。

f:id:abrahamcow:20190809212941p:plain

T1Dが糖尿病です。たしかにちょっとサチってる気がする。

サンプルごとの繰り返し測定を表現したくて、こんな図も書いてみた。

f:id:abrahamcow:20190809213421p:plain

R のコードを貼ります。

library(tidyverse)
library(readxl)

tab0 <- read_tsv("~/data/Kostic/diabimmune_t1d_16s_otu_table.txt",skip=1)
meta <- read_excel("~/data/Kostic/diabimmune_t1d_wgs_metadata.xlsx")
abundance <- tab0 %>% 
  dplyr::select(starts_with("G")) %>% 
  mutate_all(as.numeric) %>% 
  t()


S_obs <- ncol(abundance)

S0 <- rowSums(abundance>0)
F1 <- rowSums(abundance==1)
F2 <- rowSums(abundance==2)

P <- (abundance+1)/rowSums(abundance+1)

dfalpha <- tibble(Gid_16S = row.names(abundance),diversity=S0+(F1^2)/(2*F2),
                  entropy=-rowSums(P*log(P)))
datout <- left_join(dfalpha,meta) %>% 
  drop_na()

theme_set(theme_bw(16))

ggplot(datout,aes(x=Age_at_Collection,y=entropy,colour=T1D_Status))+
  geom_point()+
  stat_smooth()

ggplot(datout,aes(x=Age_at_Collection, y = Subject_ID, colour=entropy))+
  geom_point(size=10, pch=15,alpha=0.8)+
  scale_color_gradient2(mid="gray", midpoint = mean(datout$entropy))+
  facet_grid(T1D_Status~.,scales = "free_y", space = "free")

だれか種の多様性の変化のモデリングやりませんか?