GUniFrac系统发育距离结合多重距离统群落差异
写在前面
系统发育距离介绍
在多变量微生物组组成分析中，加权和未加权的 UniFrac distance是测量两种微生物组群体的的最常用方法。UniFrac距离，也被称为未加权UniFrac距离，由Lozupone等人在2005年提出(Lozupone and Knight 2005)。它是根据系统发育距离估计微生物组样品之间差异的一种方法。UniFrac距离度量的目标是使来自不同条件的微生物组样本之间的客观比较成为可能。
安装R包
# install.packages("GUniFrac")
library(GUniFrac)
数据获取
library(phyloseq)
metadata = read.delim("https://raw.githubusercontent.com/taowenmicro/R-_function/main/metadata.tsv",row.names = 1)
otutab = read.delim("https://raw.githubusercontent.com/taowenmicro/R-_function/main/otutab.txt", row.names=1)
taxonomy = read.table("https://raw.githubusercontent.com/taowenmicro/R-_function/main/taxonomy.txt", row.names=1)
tree = read_tree("https://raw.githubusercontent.com/taowenmicro/R-_function/main/otus.tree")
library(Biostrings)
rep = readDNAStringSet("https://raw.githubusercontent.com/taowenmicro/R-_function/main/otus.fa")
ps = phyloseq(sample_data(metadata),
otu_table(as.matrix(otutab), taxa_are_rows=TRUE),
tax_table(as.matrix(taxonomy)), phy_tree(tree),refseq(rep)
)
library(tidyverse)
ps = filter_taxa(ps, function(x) sum(x ) > 500, TRUE);ps #筛选序列数量大于1的
提取otu表格
为什么我要多此一举，先构建phyloseq对象，然后在进行提取呢？因为phyloseq对象可以对不同类型的一个微生物组数据进行重排和整理，保证每个部分都可以匹配上。
otu_tab <- t(otu_table(ps))
library(vegan)
otu_tab_rarefy <- rrarefy(otu_tab, min(apply(otu_tab,1,sum)))
提取进化树
otu_tree <- phy_tree(ps)
otu_tree
# Phylogenetic tree with 187 tips and 186 internal nodes.
#
# Tip labels:
# ASV_3, ASV_4, ASV_100, ASV_6, ASV_96, ASV_28, ...
#
# Rooted; includes branch lengths.
计算系统发育距离
这里运算的很慢，可能和这个函数不能调用多线程有关吧，所以建议大家进行去除低丰度OTU后减少OTU数量来做。可以加快速度。
unifracs <- GUniFrac((otu_tab_rarefy), otu_tree, alpha=c(0, 0.5, 1))
提取对应的不同系统发育距离矩阵
# class(unifracs)
# str(unifracs)
dw <- unifracs$unifracs[,,"d_1"]# Weighted UniFrac
du <- unifracs$unifracs[,,"d_UW"] # Unweighted UniFrac
dv <- unifracs$unifracs[,, "d_VAW"]# Variance adjusted weighted UniFrac
d0 <- unifracs$unifracs[,, "d_0"] # GUniFrac with alpha 0
d5 <- unifracs$unifracs[,, "d_0.5"]# GUniFrac with alpha 0.5
使用UniFrac距离检测差异
set.seed(123)
adonis(as.dist(d5) ~ metadata$Group)
# Call:
# adonis(formula = as.dist(d5) ~ metadata$Group)
#
# Permutation: free
# Number of permutations: 999
#
# Terms added sequentially (first to last)
#
# Df SumsOfSqs MeanSqs F.Model R2 Pr(>F)
# metadata$Group 2 0.10458 0.052288 2.4643 0.24731 0.001 ***
# Residuals 15 0.31828 0.021218 0.75269
# Total 17 0.42285 1.00000
# ---
# Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
使用多个距离矩阵进行差异检测
结合多个距离矩阵可以增加假设检验的能力。下面的R代码使用d(0)，d(0.5),d(1)和函数PermanovaG()进行排列多元方差分析:
PermanovaG(unifracs$unifracs[, , c("d_0", "d_0.5", "d_1")] ~ metadata$Group)
# Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
#
# $aov.tab.list[[3]]
# Permutation: free
# Number of permutations: 999
#
# Terms added sequentially (first to last)
#
# Df SumsOfSqs MeanSqs F.Model R2 Pr(>F)
# metadata$Group 2 0.08836 0.044178 2.7813 0.27052 0.001 ***
# Residuals 15 0.23826 0.015884 0.72948
# Total 17 0.32661 1.00000
# ---
# Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1
根际互作生物学研究室 简介
根际互作生物学研究室是沈其荣教授土壤微生物与有机肥团队下的一个关注于根际互作的研究小组。本小组由袁军副教授带领，主要关注：1.植物和微生物互作在抗病过程中的作用；2 环境微生物大数据整合研究；3 环境代谢组及其与微生物过程研究体系开发和应用。团队在过去三年中在 isme J， Microbiome， PCE，SBB，Horticulture Research等期刊上发表了多篇文章。欢迎关注 微生信生物 公众号对本研究小组进行了解。