如何从ERGM提取平均同质性系数以开展元回归分析？

我已在多个班级样本上运行以下ERGM模型：

mFriend<-ergm(networks$`network_id`~edges+mutual+ 
                  nodeofactor('sex') + nodeofactor('minor') +
                  nodemix ('sex', base = c(2,3))+nodemix ('minor', base = c(2,3)) + absdiff('gpa'),
                control = control.ergm(seed=1, MCMC.burnin=50000, MCMC.interval=5000), verbose=TRUE)

某一班级的模型结果摘要如下：

Monte Carlo Maximum Likelihood Results:

                    Estimate Std. Error MCMC % z value Pr(>|z|)    
edges                -3.3709     0.7549      0  -4.466  < 1e-04 ***
mutual                2.7828     0.4667      0   5.963  < 1e-04 ***
nodeofactor.sex.1    -0.8513     1.0916      0  -0.780 0.435492    
nodeofactor.minor.1   0.2849     0.9734      0   0.293 0.769728    
mix.sex.0.0           2.2518     0.6215      0   3.623 0.000291 ***
mix.sex.1.1           3.2358     0.8129      0   3.980  < 1e-04 ***
mix.minor.0.0         0.1571     0.5563      0   0.282 0.777615    
mix.minor.1.1           -Inf     0.0000      0    -Inf  < 1e-04 ***
absdiff.gpa          -0.5670     0.3446      0  -1.645 0.099939 .  

下一步我计划结合班级协变量开展元回归分析，需获取每个班级的平均同质性系数及标准误，用于执行以下元回归代码：

rma(homophily_estimates, homophily_errors, mods = ~ minor_perc+class_size+mean_income, data=dfmeta2)

由于是首次开展元回归，不清楚具体操作流程，恳请提供相关建议。

1. 提取每个班级的同质性系数及标准误 #

你的ERGM模型中，nodemix项直接对应同质性相关的对数优势比系数：

• 性别同质性：mix.sex.0.0（女性-女性连接）、mix.sex.1.1（男性-男性连接）
• 少数群体同质性：mix.minor.0.0（非少数群体-非少数群体连接）、mix.minor.1.1（少数群体-少数群体连接）

计算平均同质性系数 #

你可以选择两种方式计算平均：

• 简单平均：直接取两个同质性系数的均值
• 加权平均：按班级内对应群体的节点占比加权（更贴合实际网络结构）

代码示例（批量提取） #

假设networks是存储所有班级网络的列表，covariates是存储班级协变量的数据框（按class_id匹配）：

# 初始化结果数据框
dfmeta2 <- data.frame(
  class_id = character(),
  avg_sex_homophily = numeric(),
  avg_sex_se = numeric(),
  avg_minor_homophily = numeric(),
  avg_minor_se = numeric(),
  minor_perc = numeric(),
  class_size = numeric(),
  mean_income = numeric(),
  stringsAsFactors = FALSE
)

# 遍历每个班级
for (id in names(networks)) {
  # 运行ERGM模型
  mFriend <- ergm(networks[[id]] ~ edges + mutual + 
                    nodeofactor('sex') + nodeofactor('minor') +
                    nodemix('sex', base = c(2,3)) + nodemix('minor', base = c(2,3)) + absdiff('gpa'),
                  control = control.ergm(seed=1, MCMC.burnin=50000, MCMC.interval=5000), verbose=FALSE)
  
  # 提取系数和标准误
  coefs <- coef(mFriend)
  ses <- sqrt(diag(vcov(mFriend)))
  
  # 处理性别同质性
  sex_coefs <- coefs[c("mix.sex.0.0", "mix.sex.1.1")]
  sex_ses <- ses[c("mix.sex.0.0", "mix.sex.1.1")]
  avg_sex_homophily <- mean(sex_coefs, na.rm=TRUE)
  avg_sex_se <- mean(sex_ses, na.rm=TRUE) # 简单平均标准误，加权需用方差加权后开方
  
  # 处理少数群体同质性（解决-Inf问题）
  minor_coefs <- coefs[c("mix.minor.0.0", "mix.minor.1.1")]
  minor_ses <- ses[c("mix.minor.0.0", "mix.minor.1.1")]
  # 替换-Inf为极小值（如-10，需匹配模型系数范围），0标准误替换为极小值
  minor_coefs[is.infinite(minor_coefs)] <- -10
  minor_ses[minor_ses == 0] <- 0.0001
  avg_minor_homophily <- mean(minor_coefs, na.rm=TRUE)
  avg_minor_se <- mean(minor_ses, na.rm=TRUE)
  
  # 匹配班级协变量
  covar <- covariates[covariates$class_id == id, ]
  
  # 写入结果
  dfmeta2 <- rbind(dfmeta2, data.frame(
    class_id = id,
    avg_sex_homophily = avg_sex_homophily,
    avg_sex_se = avg_sex_se,
    avg_minor_homophily = avg_minor_homophily,
    avg_minor_se = avg_minor_se,
    minor_perc = covar$minor_perc,
    class_size = covar$class_size,
    mean_income = covar$mean_income,
    stringsAsFactors = FALSE
  ))
}

2. 元回归执行与分析 #

选择分析指标 #

确定要分析的同质性维度（性别/少数群体），分别执行元回归：

library(metafor)

# 性别同质性元回归
rma_sex <- rma(avg_sex_homophily, avg_sex_se, mods = ~ minor_perc+class_size+mean_income, data=dfmeta2)
summary(rma_sex)

# 少数群体同质性元回归
rma_minor <- rma(avg_minor_homophily, avg_minor_se, mods = ~ minor_perc+class_size+mean_income, data=dfmeta2)
summary(rma_minor)

关键检查点 #

• ERGM收敛性：提取系数前必须用mcmc.diagnostics(mFriend)检查模型收敛情况，不收敛的结果直接丢弃
• 异质性检验：查看结果中的I²统计量，若异质性高，可考虑调整模型（如加入更多协变量或使用混合效应模型）
• 敏感性分析：针对少数群体同质性的异常值处理方式，尝试不同替换值（如-8、-12）验证结果稳健性

3. 注意事项 #

• ERGM中的nodemix系数是对数优势比，元回归结果解释为协变量对同组连接相对跨组连接优势的影响
• 若班级数量较少，元回归的统计效力有限，需谨慎解读结果
• 对于mix.minor.1.1为-Inf的班级，也可选择直接排除，需在论文中说明理由

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