大数跨境
0
0

ggpicrust2优雅的探索PICRUSt2分析数据

ggpicrust2优雅的探索PICRUSt2分析数据 R语言数据分析指南
2023-08-19
1

欢迎关注R语言数据分析指南

本节来介绍一款R包「ggpicrust2」,主要用于对「PICRUSt2」输出的结果进行深度操作,ggpicrust2集成了ko丰度到kegg通路丰度转换、通路注释、差异丰度(DA)分析及结果图的可视化等。下面小编来简单介绍下,更多详细的案例内容请参考作者官方文档。「此款R包安装依赖较多,请耐心安装」

官方文档

https://github.com/cafferychen777/ggpicrust2

关注下方公众号下回更新不迷路

安装R包

if (!requireNamespace("BiocManager", quietly = TRUE))
  install.packages("BiocManager")

pkgs <- c("phyloseq""ALDEx2""SummarizedExperiment""Biobase""devtools"
          "ComplexHeatmap""BiocGenerics""BiocManager""metagenomeSeq"
          "Maaslin2""edgeR""lefser""limma""KEGGREST""DESeq2")

for (pkg in pkgs) {
  if (!requireNamespace(pkg, quietly = TRUE))
    BiocManager::install(pkg)
}

devtools::install_github("cafferychen777/ggpicrust2")

加载R包

library(ggpicrust2)
library(tidyverse)
library(GGally)
library(ggprism)
library(patchwork)
library(ggh4x)

差异丰度分析

data("metacyc_abundance")
data("metadata")
metacyc_daa_results_df <- pathway_daa(abundance = metacyc_abundance %>% 
column_to_rownames("pathway"), metadata = metadata, 
group = "Environment", daa_method = "LinDA",
select = NULL, p.adjust = "BH", reference = NULL)

多模型结果比较

methods <- c("ALDEx2""DESeq2""edgeR")
daa_results_list <- lapply(methods, function(method) {
  pathway_daa(abundance = metacyc_abundance %>% 
  column_to_rownames("pathway"), metadata = metadata, group = "Environment", daa_method = method)
})

method_names <- c("ALDEx2_Welch's t test","ALDEx2_Wilcoxon rank test","DESeq2""edgeR")
# Compare results across different methods
comparison_results <- compare_daa_results(daa_results_list = daa_results_list, method_names = method_names)
> comparison_results
                     method num_features num_common_features num_diff_features
1     ALDEx2_Welch's t test          124                  10               124
2 ALDEx2_Wilcoxon rank test           41                  10                41
3                    DESeq2           41                  10                41

kegg通路注释

data("kegg_abundance")
data("metadata")

daa_results_df <- pathway_daa(abundance = kegg_abundance,
metadata = metadata, group = "Environment", daa_method = "LinDA")

daa_annotated_results_df <- pathway_annotation(pathway = "KO",
daa_results_df = daa_results_df, ko_to_kegg = TRUE)

Ko号注释

data("ko_abundance")
data("metadata")

daa_results_df <- pathway_daa(abundance = ko_abundance %>%
column_to_rownames("#NAME"), metadata = metadata, 
group = "Environment", daa_method = "LinDA")

daa_annotated_results_df <- pathway_annotation(pathway = "KO",
daa_results_df = daa_results_df, ko_to_kegg = FALSE)

pathway_errorbar

data("ko_abundance")
data("metadata")
kegg_abundance <- ko2kegg_abundance(data = ko_abundance) 

daa_results_df <- pathway_daa(kegg_abundance, metadata = metadata,
group = "Environment", daa_method = "LinDA")
daa_annotated_results_df <- pathway_annotation(pathway = "KO",
daa_results_df = daa_results_df, ko_to_kegg = TRUE)

p <- pathway_errorbar(abundance = kegg_abundance,
                      daa_results_df = daa_annotated_results_df,
                      Group = metadata$Environment,
                      ko_to_kegg = TRUE,
                      p_values_threshold = 0.05,
                      order = "pathway_class",
                      select = NULL,
                      p_value_bar = TRUE,
                      colors = NULL,
                      x_lab = "pathway_name")
data("metacyc_abundance")
data("metadata")
metacyc_daa_results_df <- pathway_daa(abundance = metacyc_abundance %>% 
column_to_rownames("pathway"), metadata = metadata,
group = "Environment", daa_method = "LinDA")

metacyc_daa_annotated_results_df <- pathway_annotation(pathway = "MetaCyc",
daa_results_df = metacyc_daa_results_df, ko_to_kegg = FALSE)

p <- pathway_errorbar(abundance = metacyc_abundance %>% 
column_to_rownames("pathway"),
                      daa_results_df = metacyc_daa_annotated_results_df,
                      Group = metadata$Environment,
                      ko_to_kegg = FALSE,
                      p_values_threshold = 0.05,
                      order = "group",
                      select = NULL,
                      p_value_bar = TRUE,
                      colors = NULL,
                      x_lab = "description")

通路热图

data("metacyc_abundance")
data("metadata")

# Perform differential abundance analysis
metacyc_daa_results_df <- pathway_daa(
  abundance = metacyc_abundance %>% column_to_rownames("pathway"),
  metadata = metadata,
  group = "Environment",
  daa_method = "LinDA"
)

# Annotate the results

annotated_metacyc_daa_results_df <- pathway_annotation(
  pathway = "MetaCyc",
  daa_results_df = metacyc_daa_results_df,
  ko_to_kegg = FALSE
)

feature_with_p_0.05 <- metacyc_daa_results_df %>% 
  filter(p_adjust < 0.05)

pathway_heatmap(
  abundance = metacyc_abundance %>% 
    right_join(
      annotated_metacyc_daa_results_df %>% select(all_of(c("feature","description"))),
      by = c("pathway" = "feature")
    ) %>% 
    filter(pathway %in% feature_with_p_0.05$feature) %>% 
    select(-"pathway") %>% 
    column_to_rownames("description"),
  metadata = metadata, 
  group = "Environment"
)

Pathway_PCA

pathway_pca(abundance = metacyc_abundance %>% column_to_rownames("pathway"),
            metadata = metadata, group = "Environment")

本节介绍到此结束,可以看到其内置函数还是很实用的。有需要学习个性化数据可视化的朋友,欢迎到小编的「淘宝店铺」 「R语言数据分析指南」购买「2023年度会员文档」同步更新中「售价149元」,内容主要包括各种「高分论文的图表分析复现以及一些个性化图表的绘制」均包含数据+代码;按照往年数据小编年产出约在150+以上

购买后微信发小编订单截图即邀请进新的会员交流群,小编的文档为按年售卖,只包含当年度的「除系列课程外」的文档,有需要往年文档的朋友也可下单购买,需要了解更多信息的朋友欢迎交流咨询。

淘宝购买方式

2023会员群案例展示


【声明】内容源于网络
0
0
R语言数据分析指南
R语言重症爱好者,喜欢绘制各种精美的图表,喜欢的小伙伴可以关注我,跟我一起学习
内容 1180
粉丝 0
R语言数据分析指南 R语言重症爱好者,喜欢绘制各种精美的图表,喜欢的小伙伴可以关注我,跟我一起学习
总阅读551
粉丝0
内容1.2k