分组因子成对比较所得p值绘制热图的实现方法咨询

首先加载所需依赖包：

library(tidyverse)
library(rstatix)
library(ggplot2)

步骤1：计算成对检验结果并格式化 #

你原有代码得到的检验结果是长格式，需要先转为对称的热图适配格式，同时完成显著性区间划分：

# 先保存成对检验结果
p_result <- data %>%
  group_by(site) %>%
  t_test(variable ~ organism) %>%
  select(site, group1, group2, p.adj) # 选用校正后p值做展示

# 转为热图所需对称格式+划分显著性区间
heatmap_data <- p_result %>%
  # 补全反向比较结果（A和B比较的p值与B和A一致）
  bind_rows(
    p_result %>% rename(group1 = group2, group2 = group1)
  ) %>%
  # 补全对角线项（自身比较的p值统一设为1，也可按需设为NA）
  bind_rows(
    expand_grid(site = unique(p_result$site),
                group1 = unique(p_result$group1),
                group2 = unique(p_result$group1)) %>%
      filter(group1 == group2) %>%
      mutate(p.adj = 1)
  ) %>%
  # 按要求划分显著性区间
  mutate(sig_level = case_when(
    p.adj < 0.001 ~ "p < 0.001",
    p.adj < 0.01 ~ "p < 0.01",
    p.adj < 0.05 ~ "p < 0.05",
    TRUE ~ "p ≥ 0.05"
  )) %>%
  # 固定行列顺序，保证热图展示逻辑一致
  mutate(
    group1 = factor(group1, levels = c("Insects", "Mammals", "Reptiles")),
    group2 = factor(group2, levels = c("Insects", "Mammals", "Reptiles")),
    sig_level = factor(sig_level, levels = c("p ≥ 0.05", "p < 0.05", "p < 0.01", "p < 0.001"))
  )

步骤2：绘制显著性热图 #

可以按site分面同时展示两个站点的结果，也可以筛选单个站点单独绘图：

ggplot(heatmap_data, aes(x = group1, y = group2, fill = sig_level)) +
  geom_tile(color = "white", size = 1) +
  # 叠加p值文本，按需保留
  geom_text(aes(label = round(p.adj, 4)), size = 3.5) +
  # 按站点分面，不需要分面的话删掉此行+filter筛选单个site即可
  facet_wrap(~site) +
  # 自定义显著性颜色，越显著颜色越深
  scale_fill_brewer(palette = "Reds", direction = 1) +
  labs(x = "物种分类", y = "物种分类", fill = "显著性水平") +
  theme_bw(base_size = 12) +
  theme(axis.text.x = element_text(angle = 45, hjust = 1),
        panel.grid = element_blank())

可选：导出纯矩阵格式 #

如果需要输出类似pwpm函数的矩阵结果，可通过以下代码转换：

# 以SITE1为例输出p值矩阵
site1_p_matrix <- heatmap_data %>%
  filter(site == "SITE1") %>%
  pivot_wider(id_cols = group1, names_from = group2, values_from = p.adj) %>%
  column_to_rownames("group1") %>%
  as.matrix()

内容的提问来源于stack exchange，提问作者Valen78