如何基于R2相关系数值重新分组SNP数据？

核心逻辑梳理 #

要实现需求，关键要精准区分两种SNP：

• 需分配新组的SNP：在原分组内，与所有其他SNP的R2均<0.4（且必须存在至少一条与组内SNP的LD记录，无LD记录的情况归为保留原组）
• 保留原组的SNP：要么和组内至少一个SNP的R2≥0.4，要么没有任何与组内SNP的LD记录

假设数据结构 #

先定义示例数据方便测试验证：

library(data.table)

# 原遗传变异分组数据
snps <- data.table(
  snp = c("rs1", "rs2", "rs3", "rs4", "rs5", "rs6"),
  group = c(1, 1, 1, 2, 2, 3)  # rs6是组3唯一SNP，无LD记录
)

# 变异间R2相关系数数据（单向存储）
ld <- data.table(
  snp1 = c("rs1", "rs1", "rs2", "rs3", "rs4"),
  snp2 = c("rs2", "rs3", "rs3", "rs1", "rs5"),
  r2 = c(0.5, 0.3, 0.6, 0.3, 0.2)
)

修正后的代码实现 #

# 1. 处理LD数据：将单向SNP对转为双向，确保每个SNP的组内关联都能被匹配到
ld_both <- rbind(ld, ld[, .(snp1 = snp2, snp2 = snp1, r2)])

# 2. 关联原分组数据，匹配每个SNP对应的组内其他SNP的R2值
snps_ld <- snps[ld_both, on = .(snp = snp1), 
                .(snp = snp1, group, partner_snp = snp2, r2)]

# 3. 按SNP+原分组聚合，判断是否符合新组条件并生成结果
snps_result <- snps[, {
  # 获取当前SNP在原分组内的所有R2值
  snp_r2 <- snps_ld[snp == .BY$snp & group == .BY$group, r2]
  
  # 核心条件：存在LD记录 且 所有R2均小于0.4
  need_new_group <- length(snp_r2) > 0 & all(snp_r2 < 0.4)
  
  # 生成新分组：符合条件则用原组号加后缀标识，否则保留原组
  new_group <- if (need_new_group) paste0(group, "_new") else as.character(group)
  
  .(new_group = new_group)
}, by = .(snp, group)]

结果说明 #

运行代码后，snps_result的输出完全符合需求：

• rs1：与rs2的R2=0.5≥0.4 → 保留组1
• rs2：与rs1、rs3的R2均≥0.4 → 保留组1
• rs3：与rs2的R2=0.6≥0.4 → 保留组1
• rs4：与rs5的R2=0.2<0.4（组内唯一LD记录）→ 分配新组2_new
• rs5：与rs4的R2=0.2<0.4（组内唯一LD记录）→ 分配新组2_new
• rs6：无任何LD记录 → 保留组3

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