微生物群落关键种识别:一种不依赖于网络的自上而下的方法

本文主要是介绍微生物群落关键种识别:一种不依赖于网络的自上而下的方法,希望对大家解决编程问题提供一定的参考价值,需要的开发者们随着小编来一起学习吧!

  微生物群落在促进养分循环、协助植物生长、维持人体健康等方面发挥着重要的作用。群落关键种对维持微生物群落稳定性具有重要影响,识别关键种一直是微生物生态学中的热点话题。识别关键种主要有两种框架:数据驱动的方法(data driven method)去除实验(perturbation experiment)。其中数据驱动的方法主要有三种:

  • 基于共现网络的方法
  • top-down方法
  • 基于深度学习的方法

注意:数据驱动的方法确定的关键种只是可能的关键种,还需要通过去除实验进一步地验证。

  • 基于共现网络的方法主要包括:构建共现网络→划分模块→计算模块间连通度和模块内连通度→确定关键种,该方法已在之前的博客中有所介绍:计算网络节点模块内连通度(within modular degree)和模块间连通度(between modular degree)。
  • 基于深度学习的方法:这里先做个预告,代码和数据都整理好了,预计下周上线,具体可参考论文Identifying keystone species in microbial communities using deep learning
  • 本文主要介绍top-down方法,该方法源于论文:Top-down identification of keystone taxa in the microbiome。该方法通过计算Empirical Presence-abundance Interrelation (EPI)来衡量物种的重要性。

EPI指标计算的流程是:

  1. 根据物种i的有-无划分为两组:
  2. 将该物种去除,并将剩余物种的相对丰度标准化,使其和为1;
  3. 然后计算组和组的距离,即该物种的重要性,EPI;
  4. 物种EPI高于平均值+两个标准差的物种可以确定为关键种。

这里的某物种 i i i 的EPI有三种衡量方法:
在这里插入图片描述
D 1 i {D}_{1}^{i} D1i 的计算:

  1. 根据物种 i i i 的有-无划分为两组:
  2. 将该物种去除,并将剩余物种的相对丰度标准化,使其和为1;
  3. 计算组和组样品的两两间的Bray-Crutis距离。假设有5个样品A、B、C、D、E,其中组:A、B、C, 组: D、E。组和组样品的两两间的距离矩阵为:
IDABC
Dxxxxxxxxx
Exxxxxxxxx
  1. 然后取该矩阵的平均值,即为 D 1 i {D}_{1}^{i} D1i

计算 D 1 i {D}_{1}^{i} D1i R代码如下:

EPI_D1 <- function(S) {library(vegan)# InitializationN <- nrow(S)M <- ncol(S)S_01 <- ifelse(S>0,1,0)D1 <- rep(NA, N)for (i in 1:N) {# If the species is always present/absent, D1 is undefinedif (sum(S_01[i, ], na.rm = TRUE) != 0 & sum(S_01[i, ], na.rm = TRUE) != M) {print(i)ind_pres <- S_01[i, ] != 0S2 <- S[-i, , drop = FALSE]S2 <- S2 / colSums(S2)bc <- as.matrix(vegdist(t(S2)))bc2 <- bc[ind_pres,!ind_pres]D1[i] <- sum(bc2) / (sum(ind_pres) * sum(!ind_pres))}}return(D1)
}

D 2 i {D}_{2}^{i} D2i 的计算:

  1. 根据物种 i i i 的有-无划分为两组:
  2. 将该物种去除,并将剩余物种的相对丰度标准化,使其和为1;
  3. 分别计算组和组样品的平均物种组成,获得 P ‾ \overline P P (P: Presence)和 A ‾ \overline A A (A: Absence),然后计算两者的平均值。假设有5个样品A、B、C、D、E,其中组:A、B、C, 组: D、E。组和组样品平均值如下:
IDABC P ‾ \overline P P
taxa1x1x2x3average(x1,x2,x3)
taxa2y1y2y3average(y1,y2,y3)
taxa3z1z2z3average(z1,z2,z3)
IDCD A ‾ \overline A A
taxa1x1x2average(x1,x2)
taxa2y1y2average(y1,y2)
taxa3z1z2average(z1,z2)
  1. 然后计算 P ‾ \overline P P A ‾ \overline A A的Bray-Crutis距离,即为 D 2 i {D}_{2}^{i} D2i

计算 D 2 i {D}_{2}^{i} D2i R代码如下:

EPI_D2 <- function(S) {N <- nrow(S)M <- ncol(S)S_01 <- ifelse(S>0,1,0)D2 <- rep(NA, N)for (i in 1:N) {# If the species is always present/absent, D2 is undefinedif (sum(S_01[i, ], na.rm = TRUE) != 0 & sum(S_01[i, ], na.rm = TRUE) != M) {print(i)# Dividing into the two groupsind_pres <- S_01[i, ] != 0S_pres <- as.matrix(S[, ind_pres])S_abs <- as.matrix(S[, !ind_pres])# Removing the i speciesS_pres <- S_pres[-i, , drop = FALSE]S_abs <- S_abs[-i, , drop = FALSE]# NormalizingS_pres <- S_pres / colSums(S_pres)S_abs <- S_abs / colSums(S_abs)# Calculating D2D2[i] <- vegdist(rbind(rowMeans(S_pres), rowMeans(S_abs)))[1]}}return(D2)
}

Q i {Q}^{i} Qi 的计算:

  1. 根据物种 i i i 的有-无划分为两组:
  2. 将该物种去除,并将剩余物种的相对丰度标准化,使其和为1;
  3. 计算样品间的Bray-Crutis距离;
  4. 设定一定的阈值,构建样品-样品的网络,这里网络中的节点代表样品;
  5. 对网络中的节点(代表样品)赋予模块,例如:模块1代表模块2代表
  6. 计算该网络的模块度(modularity),即为 Q i {Q}^{i} Qi

计算 Q i {Q}^{i} Qi 的R代码如下:

EPI_Q <- function(S, threshold_net) {N <- nrow(S)M <- ncol(S)S_01 <- ifelse(S > 0,1,0)Q <- rep(NA, N)modularity <- function(B, s) {library(igraph)B_graph <- graph.adjacency(B, mode = "undirected")d <- degree(B_graph) # Degree of each sampleq <- sum(B) / 2Qmod <- (t(s) %*% (B - (d %*% t(d)) / (2 * q)) %*% s) / (4 * q)return(Qmod)}for (i in 1:N) {# If the species is always present/absent, Q is undefinedif (sum(S_01[i, ], na.rm = TRUE) != 0 & sum(S_01[i, ], na.rm = TRUE) != M) {print(i)# Removing the i speciesS_i <- S[-i,]# NormalizingS_i <- S_i / colSums(S_i)# Building the networkdistances_i <- as.matrix(vegdist(t(S_i)))dist_threshold <- quantile(distances_i, threshold_net)B_i <- as.matrix(distances_i <= dist_threshold)diag(B_i) <- 0s_i <- as.numeric(S_01[i, ])s_i[s_i == 0] <- -1# CalculatingQ[i] <- modularity(B_i, s_i)}}return(Q)
}

更多测试数据及R代码可参考如下连接:https://mbd.pub/o/bread/ZZ2bm5hx

这篇关于微生物群落关键种识别:一种不依赖于网络的自上而下的方法的文章就介绍到这儿,希望我们推荐的文章对编程师们有所帮助!



http://www.chinasem.cn/article/874689

相关文章

Python中反转字符串的常见方法小结

《Python中反转字符串的常见方法小结》在Python中,字符串对象没有内置的反转方法,然而,在实际开发中,我们经常会遇到需要反转字符串的场景,比如处理回文字符串、文本加密等,因此,掌握如何在Pyt... 目录python中反转字符串的方法技术背景实现步骤1. 使用切片2. 使用 reversed() 函

Python中将嵌套列表扁平化的多种实现方法

《Python中将嵌套列表扁平化的多种实现方法》在Python编程中,我们常常会遇到需要将嵌套列表(即列表中包含列表)转换为一个一维的扁平列表的需求,本文将给大家介绍了多种实现这一目标的方法,需要的朋... 目录python中将嵌套列表扁平化的方法技术背景实现步骤1. 使用嵌套列表推导式2. 使用itert

Python使用pip工具实现包自动更新的多种方法

《Python使用pip工具实现包自动更新的多种方法》本文深入探讨了使用Python的pip工具实现包自动更新的各种方法和技术,我们将从基础概念开始,逐步介绍手动更新方法、自动化脚本编写、结合CI/C... 目录1. 背景介绍1.1 目的和范围1.2 预期读者1.3 文档结构概述1.4 术语表1.4.1 核

在Linux中改变echo输出颜色的实现方法

《在Linux中改变echo输出颜色的实现方法》在Linux系统的命令行环境下,为了使输出信息更加清晰、突出,便于用户快速识别和区分不同类型的信息,常常需要改变echo命令的输出颜色,所以本文给大家介... 目python录在linux中改变echo输出颜色的方法技术背景实现步骤使用ANSI转义码使用tpu

Conda与Python venv虚拟环境的区别与使用方法详解

《Conda与Pythonvenv虚拟环境的区别与使用方法详解》随着Python社区的成长,虚拟环境的概念和技术也在不断发展,:本文主要介绍Conda与Pythonvenv虚拟环境的区别与使用... 目录前言一、Conda 与 python venv 的核心区别1. Conda 的特点2. Python v

Spring Boot中WebSocket常用使用方法详解

《SpringBoot中WebSocket常用使用方法详解》本文从WebSocket的基础概念出发,详细介绍了SpringBoot集成WebSocket的步骤,并重点讲解了常用的使用方法,包括简单消... 目录一、WebSocket基础概念1.1 什么是WebSocket1.2 WebSocket与HTTP

SQL Server配置管理器无法打开的四种解决方法

《SQLServer配置管理器无法打开的四种解决方法》本文总结了SQLServer配置管理器无法打开的四种解决方法,文中通过图文示例介绍的非常详细,对大家的学习或者工作具有一定的参考学习价值,需要的... 目录方法一:桌面图标进入方法二:运行窗口进入检查版本号对照表php方法三:查找文件路径方法四:检查 S

MyBatis-Plus 中 nested() 与 and() 方法详解(最佳实践场景)

《MyBatis-Plus中nested()与and()方法详解(最佳实践场景)》在MyBatis-Plus的条件构造器中,nested()和and()都是用于构建复杂查询条件的关键方法,但... 目录MyBATis-Plus 中nested()与and()方法详解一、核心区别对比二、方法详解1.and()

golang中reflect包的常用方法

《golang中reflect包的常用方法》Go反射reflect包提供类型和值方法,用于获取类型信息、访问字段、调用方法等,文中通过示例代码介绍的非常详细,对大家的学习或者工作具有一定的参考学习价值... 目录reflect包方法总结类型 (Type) 方法值 (Value) 方法reflect包方法总结

C# 比较两个list 之间元素差异的常用方法

《C#比较两个list之间元素差异的常用方法》:本文主要介绍C#比较两个list之间元素差异,本文通过实例代码给大家介绍的非常详细,对大家的学习或工作具有一定的参考借鉴价值,需要的朋友参考下吧... 目录1. 使用Except方法2. 使用Except的逆操作3. 使用LINQ的Join,GroupJoin