基于RStudio生成样本Variant Caller的AAchange.refGene列韦恩图求助

步骤1：环境准备与文件读取 #

先确保所需R包已安装，再读取目标目录下的所有xlsx文件并解析文件名：

# 首次运行安装所需包
install.packages(c("readxl", "VennDiagram", "dplyr", "purrr", "writexl"))

# 加载依赖包
library(readxl)
library(VennDiagram)
library(dplyr)
library(purrr)
library(writexl)

# 获取/Samples/目录下所有xlsx文件路径
file_paths <- list.files(path = "/Samples/", pattern = "\\.xlsx$", full.names = TRUE)

# 解析文件名，提取样本名与Variant Caller类型
file_metadata <- tibble(file_path = file_paths) %>%
  mutate(file_name = basename(file_path)) %>%
  # 按点分割文件名，提取样本名、caller，剩余后缀忽略
  separate(file_name, into = c("sample_id", "caller", "suffix"), sep = "\\.", extra = "merge") %>%
  # 筛选目标caller
  filter(caller %in% c("bcftools", "freebayes", "mutect2.filtered", "strelka.variants"))

步骤2：按样本分组并提取AAchange数据 #

将每个样本对应的4个caller数据分组，提取AAchange.refGene列的唯一非空值：

# 按样本分组，读取数据并整理AAchange集合
sample_aa_data <- file_metadata %>%
  group_by(sample_id) %>%
  nest() %>%
  mutate(
    caller_aa_sets = map(data, function(group_df) {
      # 读取每个caller的文件，提取AAchange列并去重、去NA
      group_df %>%
        mutate(aa_unique = map(file_path, ~ read_xlsx(.x) %>% 
                                 pull(AAchange.refGene) %>% 
                                 na.omit() %>% 
                                 unique())) %>%
        select(caller, aa_unique) %>%
        deframe() # 转换为以caller为名称的列表
    })
  )

步骤3：生成每个样本的韦恩图 #

遍历每个样本，绘制4组数据的韦恩图并保存到指定目录：

# 批量生成韦恩图
walk2(sample_aa_data$sample_id, sample_aa_data$caller_aa_sets, function(sample, aa_sets) {
  venn.diagram(
    x = aa_sets,
    filename = paste0("/Samples/Venn_", sample, ".png"),
    fill = c("#E41A1C", "#377EB8", "#4DAF4A", "#984EA3"), # 自定义配色
    alpha = 0.6,
    main = paste("AAchange Intersection - Sample", sample),
    main.cex = 1.2,
    cat.cex = 0.9,
    cex = 1.1,
    margin = 0.1
  )
})

步骤4：计算并导出交集数据 #

提取所有可能的交集组合（4组全交、两两交、三者交），并保存为汇总文件：

# 计算各样本的交集统计
sample_intersection_results <- sample_aa_data %>%
  mutate(
    # 4个caller的共同交集
    four_way_intersect = map(caller_aa_sets, ~ Reduce(intersect, .x)),
    # 所有交集组合的详细统计
    all_intersections = map(caller_aa_sets, function(sets) {
      caller_names <- names(sets)
      # 生成所有非空子集组合
      all_subsets <- lapply(1:length(caller_names), function(k) {
        combn(caller_names, k, simplify = FALSE)
      }) %>% unlist(recursive = FALSE)
      
      # 计算每个子集的交集大小与具体值
      map_dfr(all_subsets, function(subset) {
        intersect_values <- Reduce(intersect, sets[subset])
        tibble(
          combination = paste(subset, collapse = " & "),
          variant_count = length(intersect_values),
          aa_changes = list(intersect_values)
        )
      })
    })
  )

# 导出结果到xlsx文件
sample_intersection_results %>%
  select(sample_id, four_way_intersect, all_intersections) %>%
  write_xlsx(path = "/Samples/Sample_AAchange_Intersections.xlsx")

关键说明 #

• 代码自动过滤AAchange.refGene列的NA值与重复值，确保每个变异仅被统计一次
• 韦恩图文件以Venn_<样本名>.png格式保存到/Samples/目录
• 交集汇总文件包含每个样本的4组全交数据，以及所有子集组合的交集数量和具体AAchange值

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