GMSB文章七：微生物整合分析

最新推荐文章于 2024-08-15 22:29:28 发布
最新推荐文章于 2024-08-15 22:29:28 发布
阅读量647 收藏 14
点赞数 16
分类专栏： 复现文章 文章标签： 数据分析 r语言
版权声明：本文为博主原创文章，遵循 CC 4.0 BY-SA 版权协议，转载请附上原文出处链接和本声明。
本文链接：https://blog.csdn.net/H20230717/article/details/140058296
版权

欢迎大家关注全网生信学习者系列：

  • WX公zhong号：生信学习者
  • Xiao hong书：生信学习者
  • 知hu：生信学习者
  • CDSN：生信学习者2

介绍

本文通过多元方差分析和典型相关分析研究微生物（species）、细胞因子（cytokine）和短链脂肪酸（SCFA）之间的相关关系。以下是两种分析的定义：

**多元方差分析（Multivariate Analysis of Variance，简称MANOVA）是一种统计方法，用于同时分析多个因变量（dependent variables）**对一个或多个自变量（independent variables）的影响。它是一种扩展了单变量方差分析（ANOVA）的技术，允许研究者检验多个响应变量是否受到一个或多个分类自变量的影响。

  • 多维数据：MANOVA处理的是多维数据集，即每个观测值都有多个响应变量的测量值。

  • 线性模型：它基于线性模型，其中每个因变量可以表示为自变量的线性组合加上误差项。

  • 假设检验：MANOVA检验的核心是假设检验，主要检验自变量对因变量的总体影响是否显著。这包括检验自变量的主效应、交互效应以及它们对因变量的联合效应。

  • 协方差矩阵：MANOVA考虑了因变量之间的相关性，通过分析协方差矩阵来评估这种相关性。

  • Wilks’ Lambda, Pillai’s Trace, Hotelling’s Trace, Roy’s Largest Root：这些都是MANOVA中常用的统计量，用于检验自变量对因变量的影响。

加载R包

library(readr)
library(openxlsx)
library(compositions)
library(tidyverse) 
library(mia)
library(ggpubr)

导入数据

大家通过以下链接下载数据：

  • 百度网盘链接：https://pan.baidu.com/s/1fz5tWy4mpJ7nd6C260abtg
  • 提取码: 请关注WX公zhong号_生信学习者_后台发送 复现gmsb 获取提取码
df_v1 <- read_csv("./data/GMSB-data/df_v1.csv", show_col_types = FALSE)
bias_corr_species <- read_csv("./data/GMSB-data/results/outputs/bias_corr_species.csv")
raw_species <- readRDS("./data/GMSB-data/rds/primary_ancombc2_species.rds")
sig_species <- read.xlsx("./data/GMSB-data/results/outputs/res_ancombc2.xlsx", sheet = 1)

数据预处理

  • 提取差异物种丰度表

  • 合并分组变量和差异物种丰度表

df_v1 <- df_v1 %>%
  dplyr::filter(group1 != "missing")

# Microbiome data
bias_corr_species <- bias_corr_species %>%
  dplyr::rowwise() %>%
  dplyr::filter(grepl("Species:", species)|grepl("Genus:", species)) %>%
  dplyr::mutate(species = ifelse(grepl("Genus:", species), 
                        paste(strsplit(species, ":")[[1]][2], "spp."),
                        strsplit(species, ":")[[1]][2])) %>%
  dplyr::ungroup() 

# Significant taxa from ANCOM-BC2
sig_species <- sig_species %>%
  dplyr::filter(p_val < 0.05) %>%
  .$taxon

# Subset significant taxa
df_da_species <- bias_corr_species %>%
  dplyr::filter(species %in% sig_species)
df_da_species <- t(df_da_species)
colnames(df_da_species) <- df_da_species[1, ]
df_da_species <- data.frame(df_da_species[-1, , drop = FALSE], check.names = FALSE) %>%
  rownames_to_column("sampleid") %>%
  dplyr::mutate(across(-1, as.numeric))
df_da_species[is.na(df_da_species)] <- 0

# 合并分组标签和物种丰度表
df_add_species <- df_v1 %>%
  dplyr::left_join(df_da_species, by = "sampleid")
df_add_species <- data.frame(df_add_species)

sig_species <- make.names(sig_species)

head(sig_species)

[1] "A.muciniphila"     "A.onderdonkii"     "Anaerovibrio.spp." "B.adolescentis"    "B.caccae"         
[6] "B.fragilis"

函数

  • lm_eqn：提取线性模型结果

  • plot_scatter：两个变量的散点图，关联关系

lm_eqn <- function(m){
  
  a <- unname(coef(m))[1]
  b <- unname(coef(m))[2]
  p_val <- summary(m)$coef[2, "Pr(>|t|)"]
  
  b <- ifelse(sign(b) >= 0, 
             paste0(" + ", format(b, digits = 2)), 
             paste0(" - ", format(-b, digits = 2)))
  
  eq <- substitute(paste(italic(y) == a, b, italic(x), ", ", italic(p) == p_val),
                  list(a = format(a, digits = 2), b = b,
                       p_val = format(p_val, digits = 2)))
  
  return(as.character(as.expression(eq)))              
}

plot_scatter <- function(df, var_name, tax_name, y_lab) {
  
  df_fig <- df %>%
    dplyr::select(all_of(c(var_name, tax_name))) %>%
    dplyr::rename(y = var_name) %>%
    tidyr::pivot_longer(-y, names_to = "tax", values_to = "x")
  
  df_lm <- df_fig %>% 
    dplyr::group_by(tax) %>%
    dplyr::do(fit = lm(formula = y ~ x, data = .)) %>%
    dplyr::summarise(eq = lm_eqn(m = fit))
  
  df_ann <- df_fig %>% 
    dplyr::group_by(tax) %>%
    dplyr::summarise(y = ifelse(mean(y, na.rm = TRUE) > 0, 
                         0.5 * max(y, na.rm = TRUE),
                         0.2 * abs(mean(y, na.rm = TRUE))),
              x = median(x, na.rm = TRUE)) %>%
      dplyr::mutate(eq = df_lm$eq,
             y_max = 1.05 * y)
  
  fig <- df_fig %>% 
    ggplot(aes(x = x, y = y)) + 
    geom_point(alpha = 0.8, color = "#BEAED4") +
    geom_smooth(method = "lm", se = TRUE, color = "skyblue", 
                formula = y ~ x) +
    facet_wrap(.~tax, scales = "free") +
    labs(x = "Micorbial abundances", 
         y = y_lab, 
         title = NULL) + 
    geom_blank(data = df_ann, aes(y = y_max)) +
    geom_text(data = df_ann, aes(x = x, y = y, label = eq), size = 3, 
              parse = TRUE, inherit.aes = FALSE) + 
    theme_bw() +
    theme(strip.background = element_rect(fill = "white"),
          legend.position = "bottom",
          plot.title = element_text(hjust = 0.5))
  
  return(fig)
}

Microbiome vs. cytokines

差异物种和细胞因子的关联分析，采用**多重协方差分析（MANCOVA, Multivariate Analysis of Covariance）**方法来评估细胞因子和微生物物种之间的多变量关系

  • 因变量：细胞因子

  • 自变量：差异菌

t_formula <- as.formula(paste0("cbind(crp, cd14, cd163) ~ ",
                              paste0(sig_species, collapse = " + ")))
fit <- manova(t_formula, data = df_add_species)

summ <- summary(fit, test = "Pillai")

df_summ <- data.frame(summ$stats) %>%
  rownames_to_column("Taxon") %>%
  dplyr::transmute(Taxon,
                   approx.F = signif(approx.F, 2),
                   num.Df, den.Df, 
                   P = round(`Pr..F.`, 2)) %>%
  dplyr::filter(Taxon != "Residuals") %>%
  dplyr::arrange(P)

head(df_summ)

cowplot::plot_grid(
  plot_scatter(
    df = df_add_species,
    var_name = "crp",
    tax_name = "Lachnospira.spp.",
    y_lab = "crp"),
  plot_scatter(
    df = df_add_species,
    var_name = "cd14",
    tax_name = "Lachnospira.spp.",
    y_lab = "cd14"),
  plot_scatter(
    df = df_add_species,
    var_name = "cd163",
    tax_name = "Lachnospira.spp.",
    y_lab = "cd163"),
  ncol = 2, align = "hv")
Taxonapprox.Fnum.Dfden.DfP
1Lachnospira.spp.3.232120.02
2RFN20.spp.2.832120.04
3Succinivibrio.spp.1.932120.13
4B.uniformis1.432120.25
5Bifidobacterium.spp.1.432120.25
6B.fragilis1.332120.28

在这里插入图片描述

结果：自变量species对因变量细胞因子的检验结果

  • 自变量Lachnospira.spp.p值小于0.05，这表示它对至少一个因变量（crp, cd14, cd163）产生了影响，可以通过散点图查看结果；

  • 自变量Lachnospira.spp.和因变量cd14存在显著负相关。

Microbiome vs. SCFAs

差异物种和短链脂肪酸的关联分析，采用**多重协方差分析（MANCOVA, Multivariate Analysis of Covariance）**方法来评估短链脂肪酸和微生物物种之间的多变量关系

  • 因变量：短链脂肪酸

  • 自变量：差异菌

t_formula <- as.formula(paste0("cbind(acetate, valerate) ~ ",
                              paste0(sig_species, collapse = " + ")))
fit <- manova(t_formula, data = df_add_species)

summ <- summary(fit, test = "Pillai")

df_summ <- data.frame(summ$stats) %>%
  rownames_to_column("Taxon") %>%
  dplyr::transmute(Taxon,
                   approx.F = signif(approx.F, 2),
                   num.Df, den.Df, 
                   P = round(`Pr..F.`, 2)) %>%
  dplyr::filter(Taxon != "Residuals") %>%
  dplyr::arrange(P)

head(df_summ)

cowplot::plot_grid(
  plot_scatter(
    df = df_add_species,
    var_name = "acetate",
    tax_name = "B.uniformis",
    y_lab = "acetate"),
  plot_scatter(
    df = df_add_species,
    var_name = "valerate",
    tax_name = "B.uniformis",
    y_lab = "valerate"),
  ncol = 2, align = "hv")
Taxonapprox.Fnum.Dfden.DfP
1B.uniformis6.321960.00
2Megasphaera.spp.3.621960.03
3Succinivibrio.spp.3.421960.04
4Butyricimonas.spp.2.521960.09
5Paraprevotella.spp.2.421960.09
6B.fragilis2.021960.13

在这里插入图片描述

结果：自变量species对因变量短链脂肪酸的检验结果

  • 自变量B.uniformisp值小于0.05，这表示它对至少一个因变量（acetate, valerate）产生了影响，可以通过散点图查看结果；

  • 自变量B.uniformis和两个因变量acetate, valerate分别存在显著正负相关。

Cytokines vs. SCFAs

细胞因子和短链脂肪酸的关联分析，采用**多重协方差分析（MANCOVA, Multivariate Analysis of Covariance）**方法来评估细胞因子和短链脂肪酸之间的多变量关系

  • 因变量：细胞因子

  • 自变量：短链脂肪酸

fit <- manova(cbind(crp, cd14, cd163) ~ 
               acetate + valerate, 
             data = df_v1)

summ <- summary(fit, test = "Pillai")

df_summ <- data.frame(summ$stats) %>%
  rownames_to_column("Taxon") %>%
  dplyr::transmute(Taxon,
                   approx.F = signif(approx.F, 2),
                   num.Df, den.Df, 
                   P = round(`Pr..F.`, 2)) %>%
  dplyr::filter(Taxon != "Residuals") %>%
  dplyr::arrange(P)

head(df_summ)

cowplot::plot_grid(
  plot_scatter(
    df = df_add_species,
    var_name = "crp",
    tax_name = "acetate",
    y_lab = "crp"),
  plot_scatter(
    df = df_add_species,
    var_name = "cd14",
    tax_name = "acetate",
    y_lab = "cd14"),
  ncol = 2, align = "hv")
Taxonapprox.Fnum.Dfden.DfP
1acetate2.532160.06
2valerate1.232160.30

结果：自变量短链脂肪酸对因变量细胞因子的检验结果

  • 自变量acetatep = 0.06，这表示它对至少一个因变量（crp, cd14, cd163）产生了轻微影响，可以通过散点图查看结果；

  • 自变量acetate和因变量crp存在显著正相关。

生信学习者2
关注
  • 16
    点赞
  • 踩
  • 14
    收藏
    觉得还不错? 一键收藏
  • 打赏
    打赏
  • 0
    评论
专栏目录
微生物网络分析：网络鲁棒性评价
01-24
对微生物网络进行鲁棒性评价有利于我们探索微生物群落的稳定性。该资源包含exe文件和python源代码及测试数据。exe文件可在不安装任何其他资源情况下，直接在windows环境中对网络进行鲁棒性评价；有兴趣的同学还可...
miaViz:微生物组分析绘图和可视化
03-19
微生物组分析绘图和可视化 该软件包的范围是微生物组数据的绘制和可视化。接口的主要类是TreeSummarizedExperiment类。 贡献 随意通过分叉和打开拉取请求来做出贡献。请确保将所需的数据整理设计为尽可能可重用，并...
NICU_Microbiome_Study:微生物组分析和相关数据
03-15
微生物组分析和相关数据。 此回购包含手稿“表征非常早产儿的细菌肠道微生物组”的补充数据。 这包括： 完整的工作流程，从原始读取到分析。 元数据（已标识）。 元数据摘要。 对于上下文，以下是本研究中使用...
RibosomeProfiling:表征微生物的分析管道
05-05
在“RibosomeProfiling:表征微生物的分析管道”项目中，重点在于利用Python编程语言构建和优化数据分析流程。Python因其强大的库支持和易用性，在生物信息学领域广泛应用。以下是一些关键的知识点： 1. **样本准备*...
微生物组：微生物组R包
02-05
微生物组分析是生物学领域的一个重要研究方向，主要关注在特定环境或生物体内共生的各种微生物的种类、数量及其相互作用。R语言作为一种强大的统计分析和图形生成工具，为微生物组学研究提供了丰富的软件包。"微生物...
Python 绘图入门
qq_39889473的博客
08-08 740
ggplot是基于matplotlib并旨在以简单方式提高matplotlib可视化感染力的库，它采用叠加图层的形式绘制图形，例如先绘制坐标轴所在的图层，再绘制点所在的图层，最后绘制线所在的图层，但其并不适用于个性化定制图形。seaborn是基于matplotlib进行高级封装的可视化库，它支持交互式界面，使得绘制图表的功能变得愈加容易，且图表的色彩更具吸引力，可以画出丰富多样的统计图表。的过程，它可以简单地理解为将一个不易描述的事物形成一个可感知画面的过程，也就是从数据空间到图形空间的。
Python 数据可视化，怎么选出合适数据的图表
qq_39889473的博客
08-12 875
掌握数据可视化的最佳实践有助于创建既美观又有效的图表。通过选择合适的图表类型，遵循简洁明了的设计原则，并避免常见错误，你可以显著提高数据可视化的质量。希望本文提供的建议和示例能帮助你在数据可视化的过程中做出更好的决策，提升图表的表现力和清晰度。如果你对数据可视化有更多的兴趣，欢迎持续关注我们的 Python 绘图专栏。1. 使用Matplotlib绘制柱状图2. 使用 Matplotlib 绘制折线图3. Python 绘图入门4. Python绘图入门：使用Matplotlib绘制饼状图。
14.2 Pandas数据处理
在路上的专栏
08-13 641
在数据科学的世界里，Pandas 就像是你厨房里的万能工具，无论是切菜、搅拌还是调味，它都能轻松应对。本文将带你走进 Pandas 的世界，揭示它如何帮助我们在数据处理的海洋中游刃有余。从基本的数据操作，到更高级的数据分析，我们将一步步揭开 Pandas 的神秘面纱。你会发现，这不仅仅是一门技术，更是一场数据世界里的冒险之旅。
科研绘图系列：Python语言时间趋势图
H20230717的博客
08-12 110
科研绘图系列：Python语言时间趋势图
近十年首款！肝病“first-in-class”PPAR激动剂获FDA批准
北京义翘神州官方账号
08-12 944
2024年6月11日，FDA加速批准Ipsen和Genfit的Iqirvo（Elafibranor）上市，用于治疗成人原发性胆汁性胆管炎（PBC）。Elafibranor是近十年首款用于治疗这类肝病的“first-in-class”小分子药物，为患者提供了新的治疗选择。PBC是一种自身免疫性胆汁淤积肝病，特征是小胆管炎症和纤维化。由于胆管的慢性损伤，导致肝脏清除体内毒素的能力下降，可能引发肝硬化和肝功能衰竭。PPAR激动剂可通过抗炎和抗纤维化作用来减缓疾病的进展。
R语言中的命名规则
2301_79425796的博客
08-12 569
在R语言中，为变量、元素、函数或对象命名时，遵循一定的规则是至关重要的。这些规则确保代码的可读性、可维护性，并避免与R语言的内部函数和保留字产生冲突。下面我们将详细探讨R中命名的具体规则和注意事项。
4章9节：用R做数据重塑，增加变量和赋值修改，和mutate()函数的复杂用法
2301_79425796的博客
08-12 491
如何使用mutate()函数创建新变量和修改现有变量，如何通过中括号和赋值符号直接修改变量，以及如何处理因子类型的变量。
如何利用R语言学习机器学习？
人邮异步社区
08-13 1012
机器学习是近年来非常热门的一个学科，R 语言经过一段时间的发展也逐渐成为主流的编程语言之一。机器学习是人工智能的一个分支，旨在通过构建和训练计算机算法和模型，使计算机能够 从数据中学习并自主进行决策和预测。传统的计算机程序是由开发者编写特定的规则和指令以 实现特定的任务。然而，机器学习不同，它依赖数据驱动的学习过程。机器学习算法会分析大 量的输入数据，并根据这些数据的模式、趋势和统计特征来提取信息并进行预测。机器学习在 许多领域都有广泛的应用，例如，、图像和语音识别、推荐系统、金融预测、医 疗诊断等。
R语言的下载和安装
2301_79425796的博客
08-12 811
内容来自：专栏《R 语言与数据科学的终极指南：2章1节：R和RStudio的下载和安装（Windows 和 Mac）》的部分
R语言的算数运算
2301_79425796的博客
08-15 558
自然指数运算是指以自然常数ee为底的指数运算。在R语言中，我们使用exp()函数来进行自然指数运算。自然常数ee的近似值为 2.71828。这个常数在数学和科学中有广泛的应用，如在复利计算、统计学中的分布函数等。
R语言的矩阵运算
最新发布
2301_79425796的博客
08-15 432
矩阵运算是数据处理和分析的重要组成部分。矩阵是一种二维数据结构，由行和列组成。R语言提供了多种函数来操作和处理矩阵，使得数据分析更加便捷和高效。矩阵可以通过多种方式创建，最常见的方法是使用matrix()函数。我们可以通过指定数据、行数和列数来创建矩阵。
如何选择适合自己的编程语言，为什么R是非计算机专业数据分析的最佳选择，五大点告诉你
2301_79425796的博客
08-12 833
对于非计算机专业背景的学者和学生来说，选择适合自己的编程语言可能看似困难。本文将探讨为什么对于那些需要进行本科生论文、研究生论文、或者发表学术成果的用户来说，R语言是最好的选择。
4章8节：用R做数据重塑，行列命名和数据类型转换
2301_79425796的博客
08-12 342
在R语言中，行列命名和数据类型转换是数据处理中的两个基础性操作。它们不仅对数据的可读性和组织性至关重要，而且在执行数据分析、模型构建和结果解释时也扮演着重要的角色。
R语言统计分析——独立性检验
maizeman126的博客
08-12 306
Fisher精确检验的原假设是：边界固定的列联表中行和列是相互独立的。vcd包也提供了kappa()函数，可以计算混淆矩阵的Cohen's kappa值以及加权的kappa值（举例说明，混淆矩阵可以表示两位判断者对于一系列对象进行分类所得结果的一致程度）。结果解读：卡方值为13.055，对应的p值为0.001463，小于0.05，说明治疗方案和改善水平之间存在显著的相关关系。值得注意的是：这里的fisher.test()函数可以在任意行列数大于等于2的联列表上使用，但不能用于2×2的列联表。
非均匀Chemostat模型：微生物竞争与共存状态分析
本文探讨了一种特殊的微生物竞争模型，该模型基于非均匀的Chemostat环境，涉及携带质粒的微生物与无质粒微生物之间的竞争。Chemostat是一种在微生物生态学中常用的实验装置，用于培养单一或多种微生物种群。模型通过...
评论
添加红包

请填写红包祝福语或标题

个

红包个数最小为10个

元

红包金额最低5元

当前余额3.43元 前往充值 >
需支付：10.00元
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

生信学习者2

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付：¥1
获取中
扫码支付

您的余额不足，请更换扫码支付或充值

打赏作者

实付元
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明：

1.余额是钱包充值的虚拟货币，按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载，可以购买VIP、付费专栏及课程。

余额充值