如何用betadisper计算样本与异组质心的多元距离？

单个站点的处理方法 #

假设你已经整理好单个站点数据：包含物种数据框comm（行=样本，列=物种）、分组变量year_period（取值"before"/"after"）、年份year、站点标识site。

1. 计算群落距离并生成PCoA坐标

library(vegan)
# 计算Bray-Curtis距离（可按需替换为其他距离类型，如"euclidean"）
dist_mat <- vegdist(comm, method = "bray")
# 执行PCoA分析，提取样本坐标（取前2轴，如需更多轴可调整索引）
pcoa_result <- cmdscale(dist_mat, eig = TRUE, add = TRUE)
sample_coords <- as.data.frame(pcoa_result$points[, 1:2])

2. 计算before组的质心坐标

# 筛选before组样本的PCoA坐标，取各轴均值作为质心
before_centroid <- colMeans(sample_coords[year_period == "before", ])

3. 计算所有样本到before质心的欧氏距离

# 逐行计算每个样本坐标与质心的欧氏距离
sample_coords$dist_to_before_centroid <- apply(sample_coords, 1, function(x) {
  sqrt(sum((x - before_centroid)^2))
})
# 合并站点、年份、分组信息
single_site_result <- cbind(site = unique(site), year, year_period, sample_coords)

多站点批量处理 #

如果全数据集包含多个站点，可通过dplyr+purrr按站点分组批量计算：

假设全数据集full_data结构：site（站点）、year（年份）、year_period（分组）、后续列是物种数据。

library(dplyr)
library(purrr)

full_site_result <- full_data %>%
  group_by(site) %>%
  nest() %>%
  mutate(
    analysis = map(data, function(df) {
      # 提取当前站点的物种数据（假设物种列以"sp_"开头，按需修改匹配规则）
      comm <- df %>% select(starts_with("sp_"))
      # 重复单个站点的分析流程
      dist_mat <- vegdist(comm, method = "bray")
      pcoa_result <- cmdscale(dist_mat, eig = TRUE, add = TRUE)
      sample_coords <- as.data.frame(pcoa_result$points[, 1:2])
      
      # 处理无before样本的站点（避免报错）
      before_idx <- df$year_period == "before"
      if (sum(before_idx) == 0) {
        sample_coords$dist_to_before_centroid <- NA
      } else {
        before_centroid <- colMeans(sample_coords[before_idx, ])
        sample_coords$dist_to_before_centroid <- apply(sample_coords, 1, function(x) {
          sqrt(sum((x - before_centroid)^2))
        })
      }
      
      # 合并结果
      cbind(df %>% select(site, year, year_period), sample_coords)
    })
  ) %>%
  unnest(analysis) %>%
  ungroup()

用betadisper直接提取坐标的替代方法 #

如果你一定要基于betadisper的结果计算，可直接从对象中提取样本和质心坐标：

bd <- betadisper(dist_mat, group = year_period)
# 提取样本PCoA坐标
sample_coords_bd <- as.data.frame(bd$points)
# 提取before组质心（注意核对质心行名与分组名一致）
before_centroid_bd <- bd$centroids[rownames(bd$centroids) == "before", ]
# 计算距离
sample_coords_bd$dist_to_before <- apply(sample_coords_bd, 1, function(x) {
  sqrt(sum((x - before_centroid_bd)^2))
})

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