[R-meta] A potential addition to metafor random-effect structures
Interesting question, Reza. I've also wondered about using factor-analytic vcov structures like this. I think they could be potentially quite useful. As Reza noted, one application could be for multivariate meta-analysis (multivariate in the strict sense <https://www.jepusto.com/what-does-multivariate-mean/>), where each study could in principle measure effect sizes on a set of p outcomes, but in practice not every study reports all outcomes. With complete reporting for a large number of studies, using unstructured random effects variances works, but with missingness and/or a limited number of studies, struct = "UN" can be hard to fit. In my experience, the solutions end up returning correlations at the boundaries of the parameter space (e.g., r = 0.999 or r = -0.999 for a bivariate random effects model, which is equivalent to a one-factor model). For a p-dimensional structure, a d-dimensional factor model has sum(p + 1 - 1:d) parameters. So these structures might be useful just as an atheoretical model-building tool, which bridges between the low-dimensional structures like CS (2 parameters) or HCS (p + 1 parameters) and the totally unconstrained UN structure (p x (p + 1) / 2 parameters). I could also see applications where such models have a meaningful theoretical interpretation. For example, perhaps there are p outcomes, which vary in their degree of sensitivity to intervention. Studies might vary along a single latent factor of intervention potency, so strong interventions have relatively large effect sizes for all outcomes, weak interventions have relatively small effects for all outcomes. The random effect for outcome j in study i might then be described by u_ij = L_j X f_i, where f_i is the latent factor of intervention potency and L_j is the sensitivity to intervention of outcome j. I could also imagine extending this further to two or more factors---maybe intervention potency and population risk level, with u_ij = L_1j X f_1j + L_2j x f_2j? James On Sun, Feb 5, 2023 at 2:31 PM Reza Norouzian via R-sig-meta-analysis <
r-sig-meta-analysis at r-project.org> wrote:
Hi Wolfgang, Thank you for your interest. Yes, potentially we can lower G's rank but it may no longer be invertible. I haven't looked at the guts of glmmTMB but obviously they use TMB in the back end for higher speed for larger models. The other thing about rr() in glmmTMB is that my quick search didn't return any simulation studies testing how approximate this approximation can be, especially given that in practice *d* is pretty much determined by consulting the information-criteria-type model fit indices. But overall, there is some potential for this modification to help users test multivariate-multilevel models currently difficult or nearly impossible to fit. I've not been lucky enough to come across a large number of such datasets, but in the few cases where this was the case, I had to drop a few of the assumptions I had in mind which eventually led me to finding about the rank-reduced structure recently added to the glmmTMB package. I may also be looking to see if I can have such models actually fit using glmmTMB, if it allows flexibility in its `dispformula=` and `control=` arguments. Kind regards, Reza On Sun, Feb 5, 2023, 7:29 AM Viechtbauer, Wolfgang (NP) via R-sig-meta-analysis <r-sig-meta-analysis at r-project.org> wrote:
I have been doing a bit more thinking about this (can't help myself). One might consider using one of the various decompositions (e.g., SVD) to accomplish this. In fact: https://en.wikipedia.org/wiki/Low-rank_approximation Something even simpler might be to use the Cholesky decomposition, that is, if G is a p*p symmetric positive-definite var-cov matrix, then t(chol(G)) %*% chol(G) == G. So, we could use t(chol(G[1:r,])) %*% chol(G[1:r,]) as a lower rank approximation to G, with r < p. In fact,
for
struct="UN", rma.mv() uses the Cholesky decomposition anyway for
ensuring
that G is positive-definite. So it might be possible to implement this without too much difficulty. Problems might creep in though since t(chol(G[1:r,])) %*% chol(G[1:r,]) is no longer invertible (since it is
by
construction no longer of full rank), so one might need to use a generalized inverse, but whether this is actually an issue or not depends on whether one needs that inverse. Best, Wolfgang
-----Original Message----- From: R-sig-meta-analysis [mailto:
r-sig-meta-analysis-bounces at r-project.org] On
Behalf Of Viechtbauer, Wolfgang (NP) via R-sig-meta-analysis Sent: Sunday, 05 February, 2023 12:26 To: R Special Interest Group for Meta-Analysis; Yefeng Yang Cc: Viechtbauer, Wolfgang (NP) Subject: Re: [R-meta] A potential addition to metafor random-effect
structures
Hi all, I took a brief look at
and the description of the 'reduced rank' structure, but I can't
determine what
it is exactly doing. The description is focused on 'abundance data' (and
the two
references as well), so to what extent this is a method specific to this
kind of
data or whether the underlying principle can be abstracted away from
this
specific application and could be relevant for a multivariate
meta-analysis would
require digging into these references. This aside, I can see some value in using a lower-dimensional
approximation to
the var-cov matrix of a given random effect, but this is the first time
I
have
come across this idea. This might be a dissertation-level research
topic.
Best, Wolfgang
-----Original Message----- From: R-sig-meta-analysis [mailto:
r-sig-meta-analysis-bounces at r-project.org] On
Behalf Of Reza Norouzian via R-sig-meta-analysis Sent: Saturday, 04 February, 2023 19:25 To: Yefeng Yang Cc: Reza Norouzian; R meta Subject: Re: [R-meta] A potential addition to metafor random-effect
structures
Hi Yefeng, Thank you for your input. You're right. We could eliminate parameters in a UN structure (although I have not come across studies that have evaluated how this strategy compares to simplifying the overall structure to some known ones say HCS etc.). However, the point is that we don't want to do that. Rather, we want to still be able to, at least, get an approximation to the UN structure consistent with our research goals. The estimation process for the new "reduced rank" structure implemented in the glmmTMB package is different from the general approach to estimating the parameters in the variance-covariance matrices of multivariate-multilevel models (I'm linking a couple of references for the details). I'm also sharing some simulated data below where we have 20 studies and a bit of a busy categorical moderator (busy_cat) with 11 levels. Each pair of these levels co-ocurr in a good number of studies. As a result, a UN structure can, in theory, be used with this data. For a moment, pretend this is a regular multilevel model. The model troubles the 'lmer()' (which by default uses a UN-type structure) giving a warning saying: "Model failed to converge: degenerate Hessian with 1 negative eigenvalues" lmer(yi~1 + (0 + busy_cat | study), data = dat, control = lmerControl(check.nobs.vs.nRE = "ignore")) By contrast, the new structure in glmmTMB seems to approximate the variance-covariance matrix with relative ease: glmmTMB(yi~1 + rr(0 + busy_cat | study, d=9), data = dat) where *d* defines the rank of the reduced rank matrix and may be determined by consulting the indices of model fit (e.g., AICc) Returning to rma.mv, this set-up is possible but is, in principle, difficult to fit: rma.mv(yi, vi, random = ~ busy_cat | study, data = dat, struct = "UN") Of course, with larger models, fitting the UN becomes even more
challenging.
I'm not sure how much and/or what kind of assessment(s) of this new structure currently exists in the methodological literature. But on its surface, it seems that this might potentially offer some solution to the problem described above. A couple of references: https://doi.org/10.1016/j.tree.2015.09.007 https://doi.org/10.1111/2041-210X.13303 Reza #==== dat <- structure(list(study = c(1L, 1L, 1L, 1L, 1L, 1L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 3L, 3L, 3L, 3L, 3L, 3L, 3L, 4L, 4L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 5L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 6L, 7L, 7L, 7L, 7L, 7L, 7L, 7L, 7L, 7L, 7L, 7L, 8L, 8L, 8L, 9L, 9L, 9L, 9L, 9L, 9L, 9L, 9L, 9L, 9L, 9L, 10L, 10L, 10L, 10L, 10L, 10L, 10L, 10L, 10L, 10L, 10L, 11L, 11L, 11L, 11L, 11L, 11L, 11L, 11L, 11L, 11L, 11L, 12L, 12L, 12L, 12L, 13L, 13L, 13L, 13L, 13L, 14L, 14L, 14L, 14L, 14L, 14L, 14L, 14L, 15L, 15L, 15L, 15L, 15L, 15L, 15L, 15L, 15L, 15L, 15L, 16L, 16L, 16L, 16L, 16L, 16L, 16L, 16L, 16L, 16L, 16L, 17L, 17L, 17L, 17L, 17L, 17L, 17L, 17L, 17L, 18L, 18L, 18L, 18L, 18L, 18L, 18L, 18L, 18L, 18L, 18L, 19L, 19L, 19L, 19L, 19L, 19L, 19L, 19L, 19L, 19L, 19L, 20L, 20L, 20L, 20L, 20L, 20L, 20L, 20L, 20L, 20L, 20L), busy_cat =
c("A",
"B", "C", "D", "E", "F", "A", "B", "C", "D", "E", "F", "G", "H", "I", "J", "A", "B", "C", "D", "E", "F", "G", "A", "B", "A", "B", "C", "D", "E", "F", "G", "H", "I", "J", "K", "A", "B", "C", "D", "E", "F", "G", "H", "I", "J", "K", "A", "B", "C", "D", "E", "F", "G", "H", "I", "J", "K", "A", "B", "C", "A", "B", "C", "D", "E", "F", "G", "H", "I", "J", "K", "A", "B", "C", "D", "E", "F", "G", "H", "I", "J", "K", "A", "B", "C", "D", "E", "F", "G", "H", "I", "J", "K", "A", "B", "C", "D", "A", "B", "C", "D", "E", "A", "B", "C", "D", "E", "F", "G", "H", "A", "B", "C", "D", "E", "F", "G", "H", "I", "J", "K", "A", "B", "C", "D", "E", "F", "G", "H", "I", "J", "K", "A", "B", "C", "D", "E", "F", "G", "H", "I", "A", "B", "C", "D", "E", "F", "G", "H", "I", "J", "K", "A", "B", "C", "D", "E", "F", "G", "H", "I", "J", "K", "A", "B", "C", "D", "E", "F", "G", "H", "I", "J", "K"), yi = c(1.061232, 1.041498, 1.01212, 1.030785, 1.044044, 1.001106, 2.185649, 0.722472, 1.777883, -3.707484, -1.122784, -0.670218, 0.847478, -0.817499, -0.989279, -0.109316, 4.626342, -3.924966, -3.738935, 3.431953, -4.553619, 4.566486, -2.15679, 1.614955, -0.336294, 0.666886, 0.944594, 0.609327, 0.645289, 0.699607, 0.538795, 0.784128, 0.655853, 0.508843, 0.757394, 0.758911, 1.082832, 1.196894, 1.175294, 1.295843, 1.227679, 1.105074, 1.290589, 1.218274, 1.187507, 1.357256, 1.351238, 0.638988, 0.971827, 2.419408, 0.585523, 0.406051, 1.826986, 1.985741, 0.193095, 2.011539, -0.257707, -0.746552, 1.005255, 0.716152, 1.019061, 0.733833, 3.008501, 2.6867, -2.853071, 8.800149, 0.036545, 3.184804, 9.174492, 3.015245, 2.990873, 0.051937, 0.488265, 0.184795, 0.483777, 0.543814, 0.552761, 0.491369, 0.322948, 0.49982, 0.400504, 0.242472, 0.220942, 0.943544, 0.719185, 0.69481, 1.011191, 0.702694, 0.557346, 0.881592, 0.7111, 0.740943, 0.867288, 0.640434, 0.300603, -0.684002, -0.094158, 0.284941, 0.243352, 0.241784, 0.241475, 0.239465, 0.238137, 0.640465, 0.747771, 0.588343, 0.46087, 0.70215, 0.513497, 0.48016, 0.528353, -1.769666, -1.648741, 3.575931, 6.691991, 1.649686, 4.918442, 0.08828, 4.674842, 1.10579, 2.513216, -0.135326, 2.912121, 1.96222, 0.18392, 4.094903, 0.929668, -0.133818, -0.33387, 0.829354, 1.494821, -0.178365, 0.93133, 3.147792, 3.704992, 3.951833, 3.660468, 3.569095, 3.661126, 3.885397, 1.328482, 2.78613, 0.791192, 1.044481, 1.087288, 0.346358, 0.540191, 1.305145, 0.813937, 0.773615, 1.137649, 0.393769, 1.251533, 0.661325, 0.677076, 0.67366, 0.685895, 0.63965, 0.672944, 0.744525, 0.645676, 0.733516, 0.60801, 0.680176, 0.766859, 0.780003, 0.781483, 0.778856, 0.788837, 0.794191, 0.796418, 0.77453, 0.778382, 0.778289, 0.772888), vi = c(0.7474, 0.593578, 0.554073, 0.857758, 0.857874, 0.748354, 0.978036, 0.944594, 0.85457, 0.926669, 0.896048, 0.994489, 0.817942, 0.677917, 0.07352, 0.370259, 0.609184, 0.423549, 0.007955, 0.590989, 0.798034, 0.597558, 0.422634, 0.446426, 0.157957, 0.385686, 0.23976, 0.902443, 0.710356, 0.009568, 0.782508, 0.12907, 0.506288, 0.355324, 0.743611, 0.508864, 0.69998, 0.247733, 0.180387, 0.721025, 0.404901, 0.857251, 0.221638, 0.503539, 0.478148, 0.137207, 0.609075, 0.534511, 0.748871, 0.809561, 0.991642, 0.415739, 0.741997, 0.527954, 0.507353, 0.552872, 0.205752, 0.39845, 0.368142, 0.522139, 0.900889, 0.506931, 0.773571, 0.494464, 0.778418, 0.3104, 0.035137, 0.439232, 0.302984, 0.466369, 0.280695, 0.606802, 0.287334, 0.020198, 0.352628, 0.544663, 0.392097, 0.991331, 0.926401, 0.578156, 0.022029, 0.19, 0.909979, 0.934414, 0.122896, 0.363876, 0.966356, 0.450393, 0.05392, 0.211827, 0.831676, 0.772909, 0.460403, 0.833485, 0.423425, 0.458889, 0.706572, 0.140761, 0.812811, 0.378294, 0.258562, 0.821408, 0.289384, 0.005483, 0.485172, 0.533032, 0.260431, 0.915904, 0.155149, 0.229136, 0.50377, 0.131536, 0.765896, 0.77428, 0.541643, 0.179194, 0.919003, 0.654007, 0.915693, 0.82844, 0.753146, 0.656434, 0.24255, 0.396368, 0.135542, 0.615797, 0.606934, 0.214534, 0.510376, 0.555259, 0.114764, 0.112979, 0.09028, 0.811538, 0.063753, 0.558368, 0.186755, 0.851563, 0.677727, 0.418678, 0.917958, 0.182144, 0.426196, 0.206626, 0.710773, 0.619399, 0.793707, 0.213669, 0.289772, 0.151224, 0.04104, 0.919974, 0.193374, 0.675358, 0.700523, 0.141227, 0.912359, 0.570215, 0.432748, 0.299634, 0.580279, 0.305556, 0.174329, 0.796159, 0.382653, 0.469863, 0.926704, 0.648651, 0.088831, 0.339708, 0.624416, 0.73655, 0.490699, 0.550977, 0.505356), row_id = 1:175), class = "data.frame", row.names = c(NA, -175L)) #==== On Fri, Feb 3, 2023 at 9:05 PM Yefeng Yang <yefeng.yang1 at unsw.edu.au>
wrote:
Dear Reza, If I understand correctly, you are talking about simplifying the
configuration
of the random effects structure in case your designed model is over-
parameterized
(this is often the case for small meta-analyses). As far as I know,
metafor is
quite flexible in imposing constraints on the heterogeneity
variance-covariance
matrix. Briefly, arguments like sigma2, tau2, rho not only can be
estimated from
the model but also can be set manually. BTW, Wolfgang created many
shorthands of
simplified "UN". Those shorthands can meet most of the conditions (at
least for
my own cases). Unless you want to test whether a specific parameter is "significant" (if this is the case, one can use the likelihood ratio
test -
under
the null hypothesis, the statistics follow the chi-square
distribution). If I
provid any misleading answers, other experts like Wolfgang, James, and
Mike
may
want to correct me.
Best, Yefeng
________________________________ From: R-sig-meta-analysis <r-sig-meta-analysis-bounces at r-project.org
on
behalf
of Reza Norouzian <rnorouzian at gmail.com>
Sent: Saturday, 4 February 2023 5:00 To: R meta <r-sig-meta-analysis at r-project.org> Subject: [R-meta] A potential addition to metafor random-effect
structures
[You don't often get email from rnorouzian at gmail.com. Learn why this
is
important at https://aka.ms/LearnAboutSenderIdentification ]
Hi All, From time to time, I encounter situations where the number of levels for a categorical variable and/or a combination of such variables are large enough in each study (i.e., creating high-dimensional joint effect distributions) that I need to forgo adopting a "UN" structure associated with those levels in favor of a more restricted structure (e.g., "HCS"). Depending on my research goals, however, this strategy may be
suboptimal.
Recently, I noticed that the glmmTMB package has added a new random-effects structure for these cases called the "reduced rank" structure possible by imposing some restrictions on the matrix of random-effects to ensure the relevant parameters' identifiability. I'm not sure how much and/or what kind of assessment(s) of this new structure currently exists in the methodological literature. But on its surface, it seems that this might potentially offer some solution to the problem described above. Will be glad to hear your thoughts/comments on the potential of this new structure for multivariate-multilevel meta-regression models perhaps implemented in rma.mv(). Kind regards, Reza
_______________________________________________ R-sig-meta-analysis mailing list @ R-sig-meta-analysis at 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]]
_______________________________________________ R-sig-meta-analysis mailing list @ R-sig-meta-analysis at r-project.org To manage your subscription to this mailing list, go to: https://stat.ethz.ch/mailman/listinfo/r-sig-meta-analysis