[R-meta] simulating multivariate random effect meta-analysis

[Filippo Gambarota]

Hi,
I'm trying to simulate a very simple multivariate meta-analysis model.
I'm following section 5.3 of Cheung (2015, meta-analysis, a structural
equation modeling approach) and the blog post by James
(https://www.jepusto.com/simulating-correlated-smds/). I would like to
have some suggestions about the general approach. In particular, I'm
not sure if the way I use the d + random effect for calculating the
sampling variance and then sampling the error term is correct.
This is my code:

# Packages

library(metafor)
library(tidyr)
library(MASS)

# Parameters

set.seed(2023)

k <- 500 # number of studies
p <- 3 # number of outcomes
d <- 0.5 # outcomes effect size
minn <- 50 # minimum sample size
meann <- 100 # average sample size
rs <- 0.6 # sampling error correlation
r <- 0.6 # pre-post correlation (for the vi)
rho <- 0.7 # the between outcomes correlation
tau <- 0.4 # the outcomes tau

```
# Setup

sim <- tidyr::expand_grid(
  paper = 1:k,
  outcome = paste0("y", 1:p),
  d = d
)

n <- minn + rpois(k, meann - minn) # generate sample size from a poisson
sim$n <- rep(n, each = p) # append to the dataset

# between-outcomes variance covariance matrix
Cmat <- rho + diag(1 - rho, nrow = p)
Tmat <- diag(tau, nrow = p) %*% Cmat %*% diag(tau, nrow = p)

# study-specific adjustment to the average effect
sim$u <- c(t(MASS::mvrnorm(k, rep(0, p), Tmat)))

# known sampling variance computed using d + u
sim$vi <- with(sim, 2*(1-rpp) * (1/n + 1/n) + (d + u)^2 / (2*(n + n)))

# block variance covariance matrix of sampling errors
V <- metafor::vcalc(vi = vi, cluster = paper, obs = outcome, rho = rs,
data = sim)

# errors
sim$e <- MASS::mvrnorm(1, rep(0, nrow(V)), Sigma = V)

# observed effect size yi
sim$yi <- with(sim, d + u + e)

# fitting the correct model according to the parameters
fit <- rma.mv(yi, V,
              mods = ~ 0 + outcome,
              random = ~outcome|paper,
              struct = "CS",
              data = sim,
              sparse = TRUE)

summary(fit)
```
What do you think? I'm aware that the simulation is very simple.
Thanks
-- 
Filippo Gambarota
PhD Student - University of Padova
Department of Developmental and Social Psychology
Website: filippogambarota.xyz
Research Groups: Colab   Psicostat