Same result regardless of the value of rho

[wviechtb]

In a meta-regression model with correlated effects (i.e., modelweights = "CORR"), I get the same results no matter what value for rho I specify. Shouldn't the results change depending on the value of rho? Here is an example:

library(metafor)
library(robumeta)

dat <- dat.riley2003
dat$vi <- dat$sei^2

rhos <- seq(0, 1, by=.1)

res <- lapply(rhos, function(rho) 
   robu(yi ~ outcome - 1, data=dat, modelweights="CORR", studynum=study, 
        var.eff.size=vi, rho=rho, small=TRUE))

cbind(rho  = sapply(res, function(x) c(x$mod_info$rho)),
      tau2 = sapply(res, function(x) c(x$mod_info$tau.sq)))

This yields:

      rho      tau2
 [1,] 0.0 0.4283492
 [2,] 0.1 0.4283492
 [3,] 0.2 0.4283492
 [4,] 0.3 0.4283492
 [5,] 0.4 0.4283492
 [6,] 0.5 0.4283492
 [7,] 0.6 0.4283492
 [8,] 0.7 0.4283492
 [9,] 0.8 0.4283492
[10,] 0.9 0.4283492
[11,] 1.0 0.4283492

Is this a bug or am I misunderstanding how the method works?

[zackfisher]

Hi Wolfgang,

I don't think you are misunderstanding how the method works. For example, see the following:

library(metafor)
library(robumeta)

rhos <- seq(0, 1, by=.1)

res <- lapply(rhos, function(rho) 
   robu(effectsize ~ 1, data=corrdat, modelweights="CORR", studynum=studyid, 
        var.eff.size=var, rho=rho, small=TRUE))

cbind(rho  = sapply(res, function(x) c(x$mod_info$rho)),
      tau2 = sapply(res, function(x) c(x$mod_info$tau.sq)))

which results in

      rho      tau2
 [1,] 0.0 0.1539968
 [2,] 0.1 0.1540746
 [3,] 0.2 0.1541524
 [4,] 0.3 0.1542302
 [5,] 0.4 0.1543080
 [6,] 0.5 0.1543858
 [7,] 0.6 0.1544637
 [8,] 0.7 0.1545415
 [9,] 0.8 0.1546193
[10,] 0.9 0.1546971
[11,] 1.0 0.1547749

That said I'm not yet entirely sure what is going on in your data example. I'll look into it and get back to you with more specific information.

[jepusto]

Wolfgang, the invariance you've illustrated is a property of the tau-squared estimator used for the correlated effects model in robumeta. It occurs whenever you're estimating a strictly multivariate model, i.e., separate intercepts for each dimension of the model, with no more than one effect size estimate per dimension in a given study.

For such models, the tau-squared estimator is calculated by taking a moment estimator for tau-squared (similar to DL estimator) within each dimension and then pooling them across dimensions. Because the moment estimator for each dimension is unbiased, so is the weighted average thereof. The calculation ends up being very similar (though not quite identical) to this:

DL_info <- function(x) {
  res <- rma(yi = yi, vi = v_avg, data = dat, method = "DL", subset = outcome == x)
  wi <- 1 / res$vi
  c_const <- sum(wi) - sum(wi^2) / sum(wi)

  data.frame(tau2 = res$tau2, c_const = c_const, QE = res$QE)
}

outcome_levels <- unique(dat$outcome)
names(outcome_levels) <- outcome_levels
tau2_components <- do.call(rbind, lapply(outcome_levels, DL_info))

with(tau2_components, sum(tau2 * c_const) / sum(c_const))

[wviechtb]

Thanks for the replies. But help(robu) says: "For the correlated effects model the method-of-moments estimar depends on the user-specified value of rho." So is this not correct? Also, not just the estimate of tau^2 is the same -- all of the results are the same (e.g., compare res[[1]] and res[[9]]). So as far as I can tell, rho has no influence whatsoever on the results. That seems a bit strange to me, but again maybe I am not understanding correctly what this parameter is supposed to do here.

[marianolake]

Just piping in to say that rho only affects the results through the estimate of tau^2. So if tau^2 isn’t affected by rho, the regression estimates, etc won’t be affected either.

Elizabeth Tipton Associate Professor of Statistics Faculty Fellow, Institute for Policy Research Northwestern University www.bethtipton.com

On Jul 20, 2020, 10:48 AM -0500, Wolfgang Viechtbauer notifications@github.com, wrote:

Thanks for the replies. But help(robu) says: "For the correlated effects model the method-of-moments estimar depends on the user-specified value of rho." So is this not correct? Also, not just the estimate of tau^2 is the same -- all of the results are the same (e.g., compare res[[1]] and res[[9]]). So as far as I can tell, rho` has no influence whatsoever on the results. That seems a bit strange to me, but again maybe I am not understanding correctly what this parameter is supposed to do here. — You are receiving this because you are subscribed to this thread. Reply to this email directly, view it on GitHub, or unsubscribe.

[zackfisher]

I will update the documentation (and fix that typo) to make this more clear.

[jepusto]

Thanks for the replies. But help(robu) says: "For the correlated effects model the method-of-moments estimar depends on the user-specified value of rho." So is this not correct? Also, not just the estimate of tau^2 is the same -- all of the results are the same (e.g., compare res[[1]] and res[[9]]). So as far as I can tell,rho` has no influence whatsoever on the results. That seems a bit strange to me, but again maybe I am not understanding correctly what this parameter is supposed to do here.

Right. Generally, the moment estimator will depend to some extent on the user-specified rho. It is only in this particular class of models (strictly multivariate) that it does not depend on rho at all.

As an aside, this is also one of the cases where I think the "correlated effects" weights used in robumeta are most troublesome (because they don't use off-diagonal terms at all) and where the results will tend to be most discrepant with rma.mv() %>% robust() or rma.mv() %>% clubSandwich::vcovCR().

[wviechtb]

Ok, thanks for the clarification. Given that this is intended behavior, feel free to close the issue.