我在几个班级上运行了以下 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 .
我的下一步工作是对这些类进行元回归并添加协变量。为此，我需要获得每个类的平均同质系数和 std.errors 并将它们添加到元回归中：
rma(homophily_estimates, homophily_errors, mods = ~ minor_perc+class_size+mean_income, data=dfmeta2)
但是，我不确定如何继续，因为这是我第一次进行元回归。你有什么建议吗？
我认为您的目标是进行类似于以下的分析：吕贝尔斯和斯奈德斯在https://doi.org/10.1016/j.socnet.2007.03.002？我不确定你所说的“如何进行？”到底是什么意思，但这里是并尝试猜测它:)如果猜测不完全正确，我将编辑我的答案。步骤基本上是：
例如类似的事情：
library(ergm)
#> Loading required package: network
#>
#> 'network' 1.18.2 (2023-12-04), part of the Statnet Project
#> * 'news(package="network")' for changes since last version
#> * 'citation("network")' for citation information
#> * 'https://statnet.org' for help, support, and other information
#>
#> 'ergm' 4.6.0 (2023-12-17), part of the Statnet Project
#> * 'news(package="ergm")' for changes since last version
#> * 'citation("ergm")' for citation information
#> * 'https://statnet.org' for help, support, and other information
#> 'ergm' 4 is a major update that introduces some backwards-incompatible
#> changes. Please type 'news(package="ergm")' for a list of major
#> changes.
library(purrr)
library(dplyr)
#>
#> Attaching package: 'dplyr'
#> The following objects are masked from 'package:stats':
#>
#> filter, lag
#> The following objects are masked from 'package:base':
#>
#> intersect, setdiff, setequal, union
# Let's pretend the three Sampson's networks are three classes
data(samplk)
networks <- list(samplk1, samplk2, samplk3)
# Fit the models to every network
l <- map(networks, ~ ergm(.x ~ edges + nodematch("group")))
#> Starting maximum pseudolikelihood estimation (MPLE):
#> Obtaining the responsible dyads.
#> Evaluating the predictor and response matrix.
#> Maximizing the pseudolikelihood.
#> Finished MPLE.
#> Evaluating log-likelihood at the estimate.
#> Starting maximum pseudolikelihood estimation (MPLE):
#> Obtaining the responsible dyads.
#> Evaluating the predictor and response matrix.
#> Maximizing the pseudolikelihood.
#> Finished MPLE.
#> Evaluating log-likelihood at the estimate.
#> Starting maximum pseudolikelihood estimation (MPLE):
#> Obtaining the responsible dyads.
#> Evaluating the predictor and response matrix.
#> Maximizing the pseudolikelihood.
#> Finished MPLE.
#> Evaluating log-likelihood at the estimate.
names(l) <- paste0("cl", seq(along = l))
# Data frame of coefficients
d <- map(l, broom::tidy) |>
bind_rows(.id = "cl")
d
#> # A tibble: 6 × 7
#> cl term estimate std.error mcmc.error statistic p.value
#> <chr> <chr> <dbl> <dbl> <dbl> <dbl> <dbl>
#> 1 cl1 edges -2.11 0.212 0 -9.98 1.79e-23
#> 2 cl1 nodematch.group 1.73 0.318 0 5.45 5.06e- 8
#> 3 cl2 edges -2.42 0.239 0 -10.1 5.86e-24
#> 4 cl2 nodematch.group 2.47 0.334 0 7.40 1.34e-13
#> 5 cl3 edges -2.48 0.245 0 -10.1 6.28e-24
#> 6 cl3 nodematch.group 2.53 0.338 0 7.48 7.34e-14
# Artificial class-level data to be added
classd <- data.frame(
cl = c("cl1", "cl2", "cl3"),
income = rnorm(3),
class_size = c(15, 21, 18)
)
# Join model coefs with class-level data
dd <- left_join(d, classd, by = "cl")
# Fit the meta-regression
metafit <- metafor::rma(
estimate ~ income,
std.error^2,
data = dd,
subset = term == "nodematch.group"
)
summary(metafit)
#>
#> Mixed-Effects Model (k = 3; tau^2 estimator: REML)
#>
#> logLik deviance AIC BIC AICc
#> -0.5423 1.0846 7.0846 1.0846 31.0846
#>
#> tau^2 (estimated amount of residual heterogeneity): 0.0625 (SE = 0.2450)
#> tau (square root of estimated tau^2 value): 0.2500
#> I^2 (residual heterogeneity / unaccounted variability): 36.08%
#> H^2 (unaccounted variability / sampling variability): 1.56
#> R^2 (amount of heterogeneity accounted for): 35.50%
#>
#> Test for Residual Heterogeneity:
#> QE(df = 1) = 1.5646, p-val = 0.2110
#>
#> Test of Moderators (coefficient 2):
#> QM(df = 1) = 1.3664, p-val = 0.2424
#>
#> Model Results:
#>
#> estimate se zval pval ci.lb ci.ub
#> intrcpt 2.0185 0.3004 6.7190 <.0001 1.4297 2.6073 ***
#> income 0.3102 0.2654 1.1689 0.2424 -0.2099 0.8304
#>
#> ---
#> Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1