基于单倍型频率数据做AMOVA遇参数错误，求解决方案

问题背景 #

基于线粒体DNA单倍型开展AMOVA分析以探究种群结构，最初尝试用poppr包的amova函数（基于genind类）失败，推测因该方法更适配二倍体数据。转用pegas包的amova函数，基于单倍型频率计算距离矩阵时出现错误：

Error in amova(dist_matrix ~ populations, nperm = 1000) : unused argument (nperm = 1000)

虽能打印AMOVA结果，但不确定其可信度。

错误原因分析 #

1. 参数误用：pegas::amova函数不支持nperm参数，该参数属于poppr包的amova函数。pegas中需通过permutest()单独做置换检验获取显著性。
2. 数据结构错误：原始代码中haplo_freq的Population列仅4个种群，但每个单倍型列有6个数值，行数不匹配，会导致数据逻辑错误。
3. 距离矩阵选择不当：用欧氏距离基于单倍型频率计算距离，无法反映线粒体单倍型的遗传差异，AMOVA需基于遗传距离（如K2P等序列进化模型距离）。

修正后的代码示例 #

# 1. 修正单倍型频率数据：确保种群数与样本数行数匹配
haplo_freq <- data.frame(
  Population = c("A", "B", "C", "D", "E", "F"),
  Hap_1 = c(1, 0, 0, 0, 0, 0),
  Hap2 = c(1, 0, 0, 0, 11, 5),
  Hap3 = c(0, 5, 2, 0, 5, 12),
  Hap4 = c(5, 0, 0, 3, 19, 14),
  Hap5 = c(0, 0, 0, 0, 5, 0),
  Hap6 = c(0, 0, 0, 0, 4, 0),
  Hap7 = c(0, 0, 0, 0, 9, 1),
  Hap8 = c(0, 0, 0, 0, 1, 0)
)
rownames(haplo_freq) <- haplo_freq$Population
haplo_freq <- haplo_freq[, -1]

# 2. 将频率数据转换为个体水平的单倍型分配（适配AMOVA要求）
ind_hap_list <- list()
for (pop in rownames(haplo_freq)) {
  # 按频率重复单倍型名称，生成个体级数据
  hap_reps <- rep(colnames(haplo_freq), times = haplo_freq[pop, ])
  ind_hap_list[[pop]] <- data.frame(
    Population = rep(pop, length(hap_reps)),
    Haplotype = hap_reps
  )
}
ind_hap_df <- do.call(rbind, ind_hap_list)

# 3. 基于单倍型序列计算遗传距离（假设hap_seqs是DNAbin格式的单倍型序列）
# 若已有单倍型序列，先加载并命名为Hap_1到Hap8
# hap_seqs <- read.dna("haplotypes.fasta", format = "fasta")
dist_ind <- dist.dna(hap_seqs[ind_hap_df$Haplotype], model = "K80")

# 4. 运行AMOVA并做置换检验
amova_result <- amova(dist_ind ~ Population, data = ind_hap_df)
perm_result <- permutest(amova_result, nperm = 1000)

# 查看结果
print(amova_result)
print(perm_result)

关键说明 #

• 原始打印的AMOVA结果不可信：R会忽略不识别的nperm参数，但数据本身存在结构错误，且未做置换检验无法判断分化显著性。
• 必须使用遗传距离：基于序列的遗传距离才能准确反映线粒体单倍型的进化差异，是AMOVA分析的合理输入。

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