[R-meta] A potential addition to metafor random-effect structures

[Viechtbauer, Wolfgang (NP)]

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 using 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
>
>https://cran.r-project.org/web/packages/glmmTMB/vignettes/covstruct.html
>
>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 using 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 using 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 using r-project.org> on
>behalf
>>of Reza Norouzian <rnorouzian using gmail.com>
>>> Sent: Saturday, 4 February 2023 5:00
>>> To: R meta <r-sig-meta-analysis using r-project.org>
>>> Subject: [R-meta] A potential addition to metafor random-effect structures
>>>
>>> [You don't often get email from rnorouzian using 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