[R-meta] multilevel and multivariate model

Tue Jun 27 15:53:26 CEST 2023

Yefeng's observations are right on target.

Regarding aggregating over the lowest level, if the model does not include
predictors that vary within the sample, then aggregating is equivalent to
estimating a model with the disaggregated effect sizes but dropping the
lowest-level random effects. More detailed explanation here:
https://www.jepusto.com/sometimes-aggregating-effect-sizes-is-fine/

If the model does include predictors that vary within sample, then
aggregating is not a good idea because you're removing relevant variation
in the outcome.

James

On Tue, Jun 27, 2023 at 5:24 AM Filippo Gambarota via R-sig-meta-analysis <
r-sig-meta-analysis using r-project.org> wrote:

> thank you Yefeng. You are right, also robust variance is a good
> option. My only concern is about omitting the 4th level as presented
> by Wolfgang here
> (https://www.metafor-project.org/doku.php/analyses:konstantopoulos2011)
> about including or not the nested effect. In fact, my data structure
> in the lowest level has only 2 (sometimes 1) effect within each
> experiment thus maybe that level can be problematic to estimate. The
> only way to reduce model complexity should be to aggregate the lowest
> level reducing the model to a 3-level.
>
> On Tue, 27 Jun 2023 at 10:11, Yefeng Yang <yefeng.yang1 using unsw.edu.au>
> wrote:
> >
> > Hi Filippo,
> >
> > Besides what James and Wolfgang said, you can try different random
> structures and choose the so-called "best" (3 levels vs. 4 levels; nested
> vs. cross-classified) from the perspective of information criteria like
> (c)AIC. But it should be noted that some researchers think the structure
> should be informed by your understanding of the data generation mechanism
> rather than relying on some sort of measures. The point estimate is usually
> not biased. The important thing that needs to take care of is the sampling
> variance of the point estimate. In this regard, either a 3- or 4-level
> model can work well generally. Robust variance estimation also works well.
> >
> > Best,
> > Yefeng
> > ________________________________
> > From: R-sig-meta-analysis <r-sig-meta-analysis-bounces using r-project.org>
> on behalf of Filippo Gambarota via R-sig-meta-analysis <
> r-sig-meta-analysis using r-project.org>
> > Sent: Tuesday, 27 June 2023 17:46
> > To: R Special Interest Group for Meta-Analysis <
> r-sig-meta-analysis using r-project.org>
> > Cc: Filippo Gambarota <filippo.gambarota using gmail.com>
> > Subject: Re: [R-meta] multilevel and multivariate model
> >
> > Thank you James, Yes my parameters seem to be correctly recovered. My
> > only doubt was about the 4-level model that looks a little bit
> > complicated but in fact if I have effects nested within experiments
> > nested within papers seems to be the only solution. Beyond
> > aggregating, I'm not aware of any other solution for this situation.
> > What do you think?
> > Thank you!
> > Filippo
> >
> > On Tue, 27 Jun 2023 at 08:52, Viechtbauer, Wolfgang (NP) via
> > R-sig-meta-analysis <r-sig-meta-analysis using r-project.org> wrote:
> > >
> > > Hi James,
> > >
> > > No particular reason -- added to my to-do list (probably will make it
> an option, like 'sparse=TRUE' as in rma.mv()).
> > >
> > > Best,
> > > Wolfgang
> > >
> > > >-----Original Message-----
> > > >From: R-sig-meta-analysis [mailto:
> r-sig-meta-analysis-bounces using r-project.org] On
> > > >Behalf Of James Pustejovsky via R-sig-meta-analysis
> > > >Sent: Monday, 26 June, 2023 19:25
> > > >To: R Special Interest Group for Meta-Analysis
> > > >Cc: James Pustejovsky
> > > >Subject: Re: [R-meta] multilevel and multivariate model
> > > >
> > > >Filippo, Your code looks right to me. To double check it for
> yourself, try
> > > >increasing K to something much larger and verify that rma.mv()
> returns
> > > >estimates that are very close to the specified parameter values. (This
> > > >might take some time to compute because of the vcalc() call, which
> returns
> > > >a square matrix with the same number of rows as dat. If nrow(dat) is
> large,
> > > >then V can get very big.)
> > > >
> > > >Wolfgang, is there a reason that V returns a dense k x k matrix? As
> far as
> > > >I can see, V should always have block diagonal form since the cluster
> > > >argument is required. Given that, would it be feasible for vcalc() to
> > > >return a Matrix::bdiag (or the metafor equivalent)?
> > > >
> > > >James
> > > >
> > > >On Mon, Jun 26, 2023 at 10:22 AM Filippo Gambarota via
> R-sig-meta-analysis <
> > > >r-sig-meta-analysis using r-project.org> wrote:
> > > >
> > > >> Hi,
> > > >> I'm working on a meta-analysis with a multilevel (effects within the
> > > >> same paper with independent groups) and multivariate (effects within
> > > >> the same paper/experiment measured on the same participants i.e.
> > > >> different outcomes). I'm not sure if my workflow is at least
> capturing
> > > >> the effect size dependency appropriately.
> > > >>
> > > >> I have some simulated data with the same structure:
> > > >>
> > > >> ```
> > > >> library(metafor)
> > > >>
> > > >> seqw <- function(x){
> > > >>   unlist(sapply(x, function(i) 1:i))
> > > >> }
> > > >>
> > > >> set.seed(2023)
> > > >>
> > > >> K <- 50 # number of papers
> > > >> J <- sample(1:3, K, replace = TRUE) # number of studies, within
> each paper
> > > >> Z <- sample(1:2, sum(J), replace = TRUE) # number of outcomes per
> > > >> study/paper
> > > >>
> > > >> dat <- data.frame(paper = rep(rep(1:K, J), Z),
> > > >>                   exp = rep(seqw(J), Z),
> > > >>                   effect = seqw(Z))
> > > >> head(dat)
> > > >> ```
> > > >> the `paper` variable is the paper, the `exp` is the experiment
> > > >> (different experiments have different subjects) and `effect` is the
> > > >> outcome within each experiment (1 and/or 2).
> > > >>
> > > >> Then I simulate a 4-level model:
> > > >>
> > > >> ```
> > > >> set.seed(2023)
> > > >> # residual variance components
> > > >> tau2 <- 0.3
> > > >> omega2 <- 0.1
> > > >> zeta2 <- 0.1
> > > >>
> > > >> # random effects
> > > >> b0_i <- rnorm(K, 0, sqrt(tau2))
> > > >> b0_ij <- rnorm(sum(J), 0, sqrt(omega2))
> > > >> b0_ijz <- rnorm(nrow(dat), 0, sqrt(zeta2))
> > > >>
> > > >> # add to dataframe
> > > >> dat$b0_i <- b0_i[dat$paper]
> > > >> dat$b0_ij <- rep(b0_ij, Z)
> > > >> dat$b0_ijz <- b0_ijz
> > > >> dat$vi <- runif(nrow(dat), 0.05, 0.1)
> > > >> ```
> > > >> Now I create the block variance-covariance matrix where sampling
> > > >> errors are correlated within each experiment and independent across
> > > >> experiments and papers:
> > > >>
> > > >> ```
> > > >> set.seed(2023)
> > > >> # create block-matrix
> > > >> V <- vcalc(vi, cluster = paper, subgroup = exp, obs = effect, rho =
> > > >> 0.7, data = dat)
> > > >> # sampling errors
> > > >> e_ij <- MASS::mvrnorm(1, mu = rep(0, nrow(V)), Sigma = V)
> > > >> ```
> > > >> Finally I add a dummy variable for the outcome 1 or 2 and simulate
> the
> > > >> observed effects:
> > > >>
> > > >> ```
> > > >> b0 <- 0.1
> > > >> b1 <- 0.1
> > > >> # moderator
> > > >> dat$x <- ifelse(dat$effect == 1, 1, 0)
> > > >>
> > > >> # simulate effect
> > > >> dat$yi <- with(dat, (b0 + b0_i + b0_ij + b0_ijz) + b1*x + e_ij)
> > > >> dat$x <- factor(dat$x)
> > > >> dat$exp <- factor(dat$exp)
> > > >> ```
> > > >> Finally my model should be written as:
> > > >>
> > > >> ```
> > > >> fit <- rma.mv(yi, V, mods = ~0 + x, random = ~1|paper/exp/effect,
> data
> > > >> = dat, sparse = TRUE)
> > > >> ```
> > > >> My question regards if the simulated data structure is correctly
> > > >> captured by the proposed model.
> > > >> Thank you!
> > > >>
> > > >> Filippo
> > > >>
> > > >> --
> > > >> Filippo Gambarota, PhD
> > > >> Postdoctoral Researcher - University of Padova
> > > >> Department of Developmental and Social Psychology
> > > >> Website: filippogambarota.xyz
> > > >> Research Groups: Colab   Psicostat
> > > _______________________________________________
> > > R-sig-meta-analysis mailing list @ R-sig-meta-analysis using r-project.org
> > > To manage your subscription to this mailing list, go to:
> > > https://stat.ethz.ch/mailman/listinfo/r-sig-meta-analysis
> >
> >
> >
> > --
> > Filippo Gambarota, PhD
> > Postdoctoral Researcher - University of Padova
> > Department of Developmental and Social Psychology
> > Website: filippogambarota.xyz
> > Research Groups: Colab   Psicostat
> >
> > _______________________________________________
> > R-sig-meta-analysis mailing list @ R-sig-meta-analysis using r-project.org
> > To manage your subscription to this mailing list, go to:
> > https://stat.ethz.ch/mailman/listinfo/r-sig-meta-analysis
>
>
>
> --
> Filippo Gambarota, PhD
> Postdoctoral Researcher - University of Padova
> Department of Developmental and Social Psychology
> Website: filippogambarota.xyz
> Research Groups: Colab   Psicostat
>
> _______________________________________________
> R-sig-meta-analysis mailing list @ R-sig-meta-analysis using r-project.org
> To manage your subscription to this mailing list, go to:
> https://stat.ethz.ch/mailman/listinfo/r-sig-meta-analysis
>

	[[alternative HTML version deleted]]