論文

Environmental factors shaping the gut microbiome in a Dutch population

https://www.nature.com/articles/s41586-022-04567-7

s41586-022-04567-7.pdf

數(shù)據(jù)和代碼下載鏈接

https://github.com/GRONINGEN-MICROBIOME-CENTRE/DMP

論文中提供的是模擬數(shù)據(jù)集

這個(gè)分析的具體原理暫時(shí)還看不明白，當(dāng)前只能試著把代碼跑通

輸入數(shù)據(jù)集部分截圖

image.png

讀取數(shù)據(jù)集

inDFmeta <- read.table('Mock_data/taxa.txt')
inDF <- inDFmeta

對(duì)數(shù)據(jù)集進(jìn)行過(guò)濾

他這里自定義了一個(gè)函數(shù)，很長(zhǎng)很長(zhǎng)，這里把他自定義的函數(shù)準(zhǔn)備到一個(gè)文件里，然后加載

source("filterMetaGenomeDF.R")

對(duì)數(shù)據(jù)集過(guò)濾

inDFmm <- filterMetaGenomeDF(inDF,
                             presPerc = -1, 
                             minMedRelAb = -1, 
                             minMRelAb = -1, 
                             keepDomains = "All",
                             keepLevels = c("S","G","F","O","C","P","K"))
dag3S <- filterMetaGenomeDF(inDFmm,keepLevels = "S",presPerc = -1,minMRelAb = 0.0,minMedRelAb = -1)

這個(gè)是物種水平的操作

對(duì)數(shù)據(jù)集進(jìn)行操作

dag3S.t <- t.data.frame(dag3S)
dag3S.t.pa <- dag3S.t
dag3S.t.pa[dag3S.t.pa > 0] <- 1
dag3S.t.pa.rs <- rowSums(dag3S.t.pa)

使用iNEXT包進(jìn)行計(jì)算

iNEXT包的幫助文檔 https://cran.r-project.org/web/packages/iNEXT/vignettes/Introduction.html

#install.packages("iNEXT")
library(iNEXT)
D_abund <- iNEXT(dag3S.t.pa.rs, 
                 datatype = 'abundance',
                 knots = 250,
                 endpoint = sum(dag3S.t.pa.rs)*1.25)
D_abund$DataInfo$n <- 2000
D_abund$iNextEst$m <- D_abund$iNextEst$m/sum(dag3S.t.pa.rs)*2000

作圖代碼

library(ggplot2)
gg.s <- ggiNEXT(D_abund, 
              type=1, 
              se=TRUE, 
              facet.var="none", 
              color.var="site", 
              grey=FALSE) + 
  theme_classic() + 
  ylab("Number of Species") + 
  xlab("Sample size") + 
  theme(text = element_text(size = 18))
print(gg.s)

image.png

屬水平的操作

dag3G <- filterMetaGenomeDF(inDFmm,keepLevels = "G",presPerc = -1,minMRelAb = 0.0000,minMedRelAb = -1)
dag3G.t <- t.data.frame(dag3G)
dag3G.t.pa <- dag3G.t
dag3G.t.pa[dag3G.t.pa > 0] <- 1
dag3G.t.pa.rs <- rowSums(dag3G.t.pa)
D_abundG <- iNEXT (dag3G.t.pa.rs, datatype = 'abundance',knots = 250,endpoint = sum(dag3G.t.pa.rs)*1.25)
D_abundG$DataInfo$n <- 2000
D_abundG$iNextEst$m <- D_abundG$iNextEst$m/sum(dag3G.t.pa.rs)*2000

gg.g <- ggiNEXT(D_abundG, 
              type=1, 
              se=TRUE, 
              facet.var="none", 
              color.var="site", 
              grey=FALSE) + 
  theme_classic() + ylab("Number of Genera") + 
  xlab("Sample size") + theme(text = element_text(size = 18))
print(gg.g)

image.png

把兩個(gè)圖拼接到一起

gg.s + 
  theme(legend.position = c(0.8,0.2))+
  scale_color_manual(values = "red")+
  guides(color="none",shape="none",fill="none") -> p1

gg.g + 
  theme(legend.position = c(0.8,0.2))+
  scale_color_manual(values = "darkgreen")+
  guides(color="none",shape="none",fill="none") -> p2

library(patchwork)

p1+p2 +
  plot_annotation(tag_levels = "a",tag_suffix = ".")

image.png

示例數(shù)據(jù)和代碼可以在公眾號(hào)后臺(tái)留言20220614獲取

歡迎大家關(guān)注我的公眾號(hào)

小明的數(shù)據(jù)分析筆記本

小明的數(shù)據(jù)分析筆記本 公眾號(hào) 主要分享：1、R語(yǔ)言和python做數(shù)據(jù)分析和數(shù)據(jù)可視化的簡(jiǎn)單小例子；2、園藝植物相關(guān)轉(zhuǎn)錄組學(xué)、基因組學(xué)、群體遺傳學(xué)文獻(xiàn)閱讀筆記；3、生物信息學(xué)入門(mén)學(xué)習(xí)資料及自己的學(xué)習(xí)筆記！