[R-meta] Covariance matrix specification for multiple levels + specifying random slopes model

[Reza Norouzian]

Dear Mick,

As mentioned in your previous post, it might be helpful to provide a bit of
data structure (e.g., `subset(data_quad, cohort %in%
unique(cohort)[1:2])[c("cohort","timepoint","outcome",...)]` or similar)
and/or some additional context. Otherwise, any advice will primarily be
based on your syntax.

Regarding the V matrix, Wolfgang has a flexible `vcalc()` function that he
can hopefully help you with.

Please see my responses inline regarding your model-related queries:

>>>>> At the moment we landed on a potential model specification looking
like this:

rma.mv(yi, V,
                random = list( ~ timepoint | interaction(cohort, outcome),
~ outcome | cohort),
                mods = ~ outcome : log(timepoint),
                data = data_quad,
                struct = c("CAR", "HCS”))

Yes, if you want to borrow strength from other (more frequently
co-occurring) levels of outcome in cohorts, it is useful to use a single
model and let the outcome levels (co)vary according to some covariational
structure within cohorts (ideally "UN" but if not possible, perhaps a
simpler structure like "HCS" etc.).

However, regardless, the above model would likely require a bit of
modification to improve interpretation. Specifically, since you are using
"CAR" to measure the effect of timepoint in cohorts, it would be wiser to
distinguish between within-cohort effect of timepoint and between-cohort
effect of timepoint by first aggregating the timepoint at the cohort level
to enable capturing between-cohort changes in timepoint:

data_quad$timepoint_Bet  <- with(data_quad, ave(timepoint, cohort))

And subtracting those aggregates from the timepoints to enable capturing
within-cohort changes in timepoint:

data_quad$timepoint_Wit <- with(data_quad, timepoint - timepoint_Bet)

And then refitting your initial model as:

rma.mv(yi, V,
       random = list( ~  timepoint_Wit | interaction(cohort, outcome), ~
outcome | cohort),
       mods = ~ log(timepoint_Bet) + log(timepoint_Wit):outcome,
       data = data_quad,
       struct = c("CAR", "HCS"))

Depending on your research question/study goals, you may then want to
interpret the fixed effect of timepoint_Wit to your audience to show the
effect of timepoint within the cohorts (typically of immediate substantive
interest) and/or the cohort-aggregated effect of timepoint (timepoint_Bet)
which typically is not of immediate substantive interest (but to remain in
the model to help reduce the bias in the other model coefficients).

>>>>> Another approach we considered was using a random slopes and
intercept model, with a random slope for outcome (3 levels). Would this
just be a case of changing the struct to “GEN” for the 2nd parameter? But
then that is outcomes nested within cohorts. So it just needs to be 1 |
outcome then?

The "GEN" structure for your *outcome* variable is essentially a more
general form of "UN". For example, compare the equivalency between the
following two models:

(res1 <- rma.mv(yi~ outcome + 0, vi, data =  dat.berkey1998,
                random = ~ outcome | trial, struct = "UN"))

(res2 <- rma.mv(yi~ outcome + 0, vi, data =  dat.berkey1998,
                random = ~ outcome + 0 | trial, struct = "GEN"))

fitstats(res1, res2)

So, fitting a "GEN" structure just like a "UN" structure might be a bit
difficult in practice, which I'm assuming was why you opted for a simpler
"HCS" structure. So, your current HCS structure which nests the outcomes in
cohorts should work as long as you have checked its performance.

That said, you could *in theory* add a GEN structure to capture the
variation in the intercepts and slopes of cohort-specific regression lines
fit onto the true effects in cohorts, but that likely will be an overfit,
plus rma.mv() won't support three *~inner | outer* random terms in a single
model:

## DON'T RUN:
rma.mv(yi, V,
       random = list( ~  timepoint_Wit | interaction(cohort, outcome), ~
timepoint_Wit | interaction(cohort, outcome), ~ outcome | cohort),
       mods = ~  log(timepoint_Bet) + log(timepoint_Wit):outcome,
       data = data_quad,
       struct = c("GEN","CAR", "HCS"))

Please think of these discussions as simply suggestions primarily based on
your original syntax.

Kind regards,
Reza


On Wed, Feb 22, 2023 at 7:00 PM Mick Girdwood via R-sig-meta-analysis <
r-sig-meta-analysis using r-project.org> wrote:

> Hi everyone,
>
> Thank you in advance for your help, the advice from this board is always
> friendly and indispensable.
>
> I am conducting a longitudinal meta-analysis, (thank you Reza for your
> advice relating to a previous question to do with this
> https://stat.ethz.ch/pipermail/r-sig-meta-analysis/2023-January/004349.html
> ). I have two things I was hoping you could help me with. Similar to this
> previous question here
> https://stat.ethz.ch/pipermail/r-sig-meta-analysis/2021-April/002794.html
>
> At the moment I had planned to use the 'longitudinal meta-analysis' per
> measure (outcome) of interest. However it occurred to our team that we have
> many similar measures of the same outcome (in this case quadriceps muscle
> strength). We wondered whether it would be better to use all this
> information in one analysis (instead of running separately for each, and
> use each measure to borrow strength). The reason being that the different
> measures are correlated with each other (as well as the timepoints). At the
> moment we landed on a potential model specification looking like this:
>
> rma.mv(yi, V,
>                 random = list( ~ timepoint | interaction(cohort, outcome),
> ~ outcome | cohort),
>                 mods = ~ outcome : log(timepoint),
>                 data = data_quad,
>                 struct = c("CAR", "HCS”))
>
> We have 3 different levels for outcome and timepoint is a continuous
> measure (hence CAR)
>
> 1. The thing I’m struggling with is how to get an appropriate V covariance
> matrix. Previously before including outcome I had used
> impute_covariance_matrix() with ‘cohort’ as the cluster. But now there is
> another level of correlated effects to consider - outcome. How can I
> account for this? Or should I just fix phi as the correlation parameter for
> the 2nd random component relating to outcome?
>
> 2. Another approach we considered was using a random slopes and intercept
> model, with a random slope for outcome (3 levels). Would this just be a
> case of changing the struct to “GEN” for the 2nd parameter? But then that
> is outcomes nested within cohorts. So it just needs to be 1 | outcome then?
>
> Thank you so much for your help, I am indebted!
>
> Mick Girdwood
> La Trobe University | Australia
>
>
>
>
>
>
>         [[alternative HTML version deleted]]
>
> _______________________________________________
> 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]]