SIBER包绘制小样本校正SEAc椭圆的方法求助

核心问题排查 #

你遇到的addEllipse只生成单个椭圆、plotSiberObject仅输出SEA椭圆的问题，本质是未正确构建多分组的SIBER对象，以及未针对SEAc指定正确的计算参数。以下是分步解决方法：

1. 重新构建包含group-community组合的SIBER对象 #

首先需要将group和community合并为一个独立分组列，让SIBER识别出5×2=10个独立组：

library(SIBER)
library(dplyr)

# 假设你的数据集为df，包含列：group, community, d13C, d15N
df <- df %>% mutate(combined_group = paste(group, community, sep = "_"))

# 构建SIBER对象，指定合并后的分组列
siber_obj <- createSiberObject(
  df,
  group = "combined_group",
  isoNames = c("d13C", "d15N") # 替换为你的同位素列名
)

2. 计算SEAc椭圆参数 #

SEAc是小样本校正后的椭圆，需要调用calcEllipse时指定p = NULL（SEA对应p=0.95这类置信水平）：

# 批量计算所有10个组的SEAc椭圆坐标
seac_ellipses <- lapply(1:siber_obj$n.groups, function(i) {
  calcEllipse(
    siber_obj$ML.mu[[i]],    # 第i组的均值
    siber_obj$ML.cov[[i]],   # 第i组的协方差矩阵
    p = NULL,                # 指定计算SEAc
    n = 100                  # 椭圆的点数，越大越平滑
  )
})

3. 手动绘制SEAc椭圆 #

先初始化画布，再循环绘制所有椭圆，用颜色/线型区分group和community：

# 初始化画布，可选择绘制散点
plotSiberObject(
  siber_obj,
  plot.type = 1,
  bty = "L",
  xlab = expression(paste(delta^{13}, "C (\u2030)")),
  ylab = expression(paste(delta^{15}, "N (\u2030)")),
  col = "darkgrey" # 散点颜色
)

# 定义样式：5个group用不同颜色，2个community用不同线型
group_cols <- rep(c("#E69F00", "#56B4E9", "#009E73", "#F0E442", "#0072B2"), each = 2)
community_ltys <- rep(c(1, 2), 5)

# 循环绘制所有SEAc椭圆
for (i in 1:length(seac_ellipses)) {
  lines(
    seac_ellipses[[i]],
    col = group_cols[i],
    lty = community_ltys[i],
    lwd = 2
  )
}

# 可选：添加图例
legend("topright",
       legend = unique(df$combined_group),
       col = group_cols,
       lty = community_ltys,
       bty = "n",
       cex = 0.8)

4. 修复addEllipse的用法 #

如果坚持用addEllipse，需要循环为每个组传入对应的均值和协方差：

plotSiberObject(siber_obj, plot.type = 1, bty = "L")

for (i in 1:siber_obj$n.groups) {
  addEllipse(
    siber_obj$ML.mu[[i]],
    siber_obj$ML.cov[[i]],
    p = NULL,
    col = group_cols[i],
    lty = community_ltys[i],
    lwd = 2
  )
}

关键注意事项 #

• 确保createSiberObject的group参数是合并后的combined_group，否则SIBER无法识别10个独立组。
• SEAc的小样本校正逻辑已内置在calcEllipse(p = NULL)中，会自动根据每组样本量调整椭圆大小。
• 若需要填充椭圆，将lines替换为polygon，并通过adjustcolor设置透明度：polygon(seac_ellipses[[i]], col = adjustcolor(group_cols[i], alpha.f = 0.2), border = group_cols[i])

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