[R-meta] R-sig-meta-analysis Digest, Vol 89, Issue 13
St Pourcain, Beate
Be@te@StPourc@|n @end|ng |rom mp|@n|
Tue Dec 17 10:37:01 CET 2024
Dear Wolfgang,
Many thanks for the new deltamethod in metafor. This helps a lot and makes many calculations easier (I used metafor_4.7-57 from your development site now). However, I do have troubles estimating more twin correlations for more complex models taking further fixed effects into account, using this function.
For example:
model.n.age.informant.par.ml <- rma.mv(
zi,
vzi,
random = list(~ MZDZ_factor | Study_ID, ~ MZDZ_factor | ESID_code), struct="UN",
data = data,
mods = ~ 0 + MZDZ_factor + Age_c + code_par,
method = "ML",
tdist = TRUE)
model.n.age.informant.par.ml
I would like to estimate
h2funct2<- function(b4,b3,b2,b1,age,rep) 2*(transf.ztor(b2+b3*age+b4*rep) - transf.ztor(b1+b3*age+b4*rep))
deltamethod(model.n.age.informant.par.ml, fun=function(b4,b3,b2,b1,age,rep) 2*(transf.ztor(b2+b3*age+b4*rep) - transf.ztor(b1+b3*age+b4*rep)))
I have troubles doing this directly in the new metafor version as the estimates do not really make sense (e.g. around -1 for twin correlations). It works better with the msm deltamethod. How, the msm deltamethod is now masked by the new metafor deltamethod (I need to load an older version of metafor to load the msm deltamethod).
Could you please advise what the best way would be to estimate these adjusted twin correlations?
Many thanks,
Beate
Beate St Pourcain, PhD
Senior Investigator & Group Leader
Room A207
Max Planck Institute for Psycholinguistics | Wundtlaan 1 | 6525 XD Nijmegen | The Netherlands
@bstpourcain
Tel: +31 24 3521964
Fax: +31 24 3521213
ORCID: https://orcid.org/0000-0002-4680-3517
Web: https://www.mpi.nl/departments/language-and-genetics/projects/population-variation-and-human-communication/
Further affiliations with:
MRC Integrative Epidemiology Unit | University of Bristol | UK
Donders Institute for Brain, Cognition and Behaviour | Radboud University | The Netherlands
My working hours may not be your working hours. Please do not feel obligated to reply outside of your normal working schedule.
-----Original Message-----
From: Viechtbauer, Wolfgang (NP) <wolfgang.viechtbauer using maastrichtuniversity.nl>
Sent: Friday, December 6, 2024 10:28 AM
To: St Pourcain, Beate <Beate.StPourcain using mpi.nl>; R Special Interest Group for Meta-Analysis <r-sig-meta-analysis using r-project.org>
Subject: RE: R-sig-meta-analysis Digest, Vol 89, Issue 13
I assume you are referring to model 'model0r.mlcon' below and coef(model0r.mlcon)[1] is the DZ coefficient (i.e., the estimated average r-to-z transformed ICC) and coef(model0r.mlcon)[2] the MZ coefficient and vcov(model0r.mlcon) gives the corresponding var-cov matrix. Then
res <- deltamethod(model0r.mlcon, fun=function(z1,z2) c(tanh(z1), tanh(z2))) res
will give you the back-transformed estimates (i.e., estimated ICCs) with corresponding SEs and CIs and
res$vcov
is the full var-cov matrix of these back-transformed estimates. And
deltamethod(model0r.mlcon, fun=function(z1,z2) 2*(tanh(z2)-tanh(z1)))
will give you the estimated H^2 based on Falconer's formula with corresponding SE/CI.
Best,
Wolfgang
> -----Original Message-----
> From: St Pourcain, Beate <Beate.StPourcain using mpi.nl>
> Sent: Tuesday, December 3, 2024 19:47
> To: Viechtbauer, Wolfgang (NP)
> <wolfgang.viechtbauer using maastrichtuniversity.nl>; R Special Interest
> Group for Meta-Analysis <r-sig-meta-analysis using r-project.org>
> Subject: RE: R-sig-meta-analysis Digest, Vol 89, Issue 13
>
> Dear Wolfgang,
>
> Many thanks for your advice and sorry about the delay.
> I added my answers below! Most things are sorted and clear! However, I
> did not follow your reply about the back-transformation of the
> directly estimated co/variance for fixed effects, i.e. from the rZ to
> the r the level, using the deltamethod. I can see how you derive the
> SE for the fixed effects at the back- transformed level with
> deltamethod, but not the variance and covariance. Could you please help here?
>
> Many thanks,
> Beate
>
> -----Original Message-----
> From: Viechtbauer, Wolfgang (NP)
> <wolfgang.viechtbauer using maastrichtuniversity.nl>
> Sent: Monday, November 11, 2024 4:13 PM
> To: St Pourcain, Beate <Beate.StPourcain using mpi.nl>; R Special Interest
> Group for Meta-Analysis <r-sig-meta-analysis using r-project.org>
> Subject: RE: R-sig-meta-analysis Digest, Vol 89, Issue 13
>
> Dear Beate,
>
> Some brief comments from my side to your questions below.
>
> Best,
> Wolfgang
>
> > -----Original Message-----
> > From: St Pourcain, Beate <Beate.StPourcain using mpi.nl>
> > Sent: Tuesday, November 5, 2024 19:16
> > To: Viechtbauer, Wolfgang (NP)
> > <wolfgang.viechtbauer using maastrichtuniversity.nl>; R Special Interest
> > Group for Meta-Analysis <r-sig-meta-analysis using r-project.org>
> > Subject: RE: R-sig-meta-analysis Digest, Vol 89, Issue 13
> >
> > Dear Wolfgang,
> > We have taken the project a few steps further now following
> vanhouwelingen2002.
> > Thanks for all your advice! However, bivariate models have of course
> > a slightly different structure compared to nested models that we
> > usually use and we have transformations and multiple random effects.
> > Thus, we had a few questions. We have a model like this:
> >
> > Study_ID: Unique ID for each independent cohort
> > ESID_code: Unique ID for study-specific assessment coded as nested
> > within each study (following Austerberry)
> > MZDZ: MZ coded 1, DZ coded 0
> > zi and vzi: Z transformed correlations and their variance
> >
> > model0r.mlcon <- rma.mv(
> > zi,
> > vzi,
> > random = list(~ MZDZ_factor | Study_ID, ~ MZDZ_factor |
> > ESID_code), struct="UN",
> > data = data,
> > mods = ~ 0 + MZDZ_factor,
> > method = "ML",
> > tdist = TRUE,
> > cvvc="varcov",
> > control=list(nearpd=TRUE))
> >
> > 1) We changed now from REML to ML to allow for model building with fixed
> > effects. Is there any major drawback?
>
> Unless your sample sizes are small, the difference between ML and REML
> should be negligible. For model selection including fixed effects, ML
> is the sensible choice.
>
> R: Perfect. Thank you!
>
> > 2) What would be the best way to estimate I2 for each MZ and DZ group using
> Z
> > and the raw correlations? (We have an idea how to do this with
> > nested models but not for bivariate models).
>
> Maybe this is useful:
>
> https://www.metafor-
> project.org/doku.php/tips:i2_multilevel_multivariate#multivariate_mode
> ls
>
> R: Very helpful, thank you!
>
> > Also does it make sense to estimate heterogeneity in rZ as:
> >
> > tau_het<-round(model0r.mlcon$tau2[1] + model0r.mlcon$tau2[2] -
> > 2*model0r.mlcon$rho*sqrt(model0r.mlcon$tau2[1]*model0r.mlcon$tau2[2]
> > ),
> > 3)
> >
> > gamma_het<-round(model0r.mlcon$gamma2[1] + model0r.mlcon$gamma2[2] -
> > 2*model0r.mlcon$phi*sqrt(model0r.mlcon$gamma2[1]*model0r.ml$gamma2[2
> > ])
> > , 3)
> >
> > res_het<-tau_study + gamma_type
>
> I can't answer that because you would first have to define what you
> mean by 'residual heterogeneity'.
>
> R: We aimed to estimate the heterogeneity in rZ-scores, analogous to
> the bivariate example approach (but now with two random effects)
> https://www.metafor-project.org/doku.php/analyses:vanhouwelingen2002
>
> > 3) We can estimate the SE for the raw scores
> > rDZse<- deltamethod(~ (exp(2*x1) - 1) / (exp(2*x1) + 1),
> > coef(model0r.mlcon),
> > vcov(model0r.mlcon))
> > rMZse<- deltamethod(~ (exp(2*x2) - 1) / (exp(2*x2) + 1),
> > coef(model0r.mlcon),
> > vcov(model0r.mlcon))
> > Is there an option to backtransform vcov(model0r.mlcon) from the
> > modelled rZ to raw r scores? Especially we are interested in the
> > covariance for raw r to derive a correlation between raw rMZ and
> > rDZ,
> including SE.
>
> The deltamethod() function from metafor can give you the full
> back-transformed var-cov matrix.
>
> R: Sorry for being slow here, I do not know how to proceed. The
> variance covariance matrix is directly estimated and there is not
> transformation involved. Could you please advise
>
> > 4) For adding further covariates as fixed effects, we were wondering how to
> > best predict manually rMZ and rDZ for different levels of the
> > covariate e.g. age
> >
> > model.age.r.mlcon <- rma.mv(
> > zi,
> > vzi,
> > random = list(~ MZDZ_factor | Study_ID, ~ MZDZ_factor |
> > ESID_code), struct="UN",
> > data = data,
> > mods = ~ 0 + MZDZ_factor + MZDZ_factor:I(Age -24), #age centered at 24 years
> > method = "ML",
> > tdist = TRUE)
> > model.age.r.mlcon snippet
> >
> > estimate se tval df pval
> > MZDZ_factor0 0.7264 0.1201 7.7141 78 <.0001
> > MZDZ_factor1 1.3193 0.2464 6.1661 78 <.0001
> > MZDZ_factor0:I(Age - 24) 0.0160 0.0040 4.0432 78 0.0001
> > MZDZ_factor1:I(Age - 24) 0.0276 0.0038 7.2905 78 <.0001
> >
> > We could now predict the twin correlations at different ages using
> > the predict function (and did this). However, this does not help
> > with the deltamethod. Is there a more elegant way to estimate the SE?
> > e.g.
> > age.20 <- -4
> > rDZ_20 <-(exp(2*(model.age.r.mlcon$b[1] +
> > model.age.r.mlcon$b[3]*age.20)) - 1) /
> > (exp(2*(model.age.r.mlcon$b[1]
> > + model.age.r.mlcon$b[3]*age.20)) + 1)
> >
> > rMZ_20 <-(exp(2*(model.age.r.mlcon$b[2] +
> > model.age.r.mlcon$b[4]*age.20)) - 1) /
> > (exp(2*(model.age.r.mlcon$b[2]
> > + model.age.r.mlcon$b[4]*age.20)) + 1)
> >
> > rDZ_20.se <-deltamethod(~(exp(2*(x1 + x3*age.20)) - 1) / (exp(2*(x1
> > +
> > x3*age.20)) + 1), coef(model.age.r.mlcon), vcov(model.age.r.mlcon))
> >
> > rMZ_20.se <-deltamethod(~(exp(2*(x2 + x4*age.20)) - 1) / (exp(2*(x2
> > +
> > x4*age.20)) + 1), coef(model.age.r.mlcon), vcov(model.age.r.mlcon))
>
> I don't understand the question. What do you mean by 'more elegant'?
> R: No worries, we will sort this out. I just meant a "shorter" way,
> but this is not a problem.
>
> > Many thanks again for all your help, Beate
> >
> > Beate St Pourcain, PhD
> > Senior Investigator & Group Leader
> > Room A207
> > Max Planck Institute for Psycholinguistics | Wundtlaan 1 | 6525 XD
> > Nijmegen | The Netherlands
> >
> > @bstpourcain
> > Tel: +31 24 3521964
> > Fax: +31 24 3521213
> > ORCID: https://orcid.org/0000-0002-4680-3517
> > Web:
> > https://www.mpi.nl/departments/language-and-genetics/projects/popula
> > ti
> > on-
> > variation-and-human-communication/
> > Further affiliations with:
> > MRC Integrative Epidemiology Unit | University of Bristol | UK
> > Donders Institute for Brain, Cognition and Behaviour | Radboud
> > University | The Netherlands
> >
> > My working hours may not be your working hours. Please do not feel
> > obligated to reply outside of your normal working schedule.
> >
> > -----Original Message-----
> > From: Viechtbauer, Wolfgang (NP)
> > <wolfgang.viechtbauer using maastrichtuniversity.nl>
> > Sent: Thursday, October 31, 2024 1:22 PM
> > To: R Special Interest Group for Meta-Analysis
> > <r-sig-meta-analysis using r- project.org>
> > Cc: St Pourcain, Beate <Beate.StPourcain using mpi.nl>
> > Subject: RE: R-sig-meta-analysis Digest, Vol 89, Issue 13
> >
> > Hi Beate,
> >
> > The simulations I mentioned are not published, but Fisher (1925)
> > directly states this as well, so you can stick to that reference.
> >
> > By the way, I just added a deltamethod() function to the development
> > version of the metafor package. So, coming back to this model:
> >
> > model <- rma.mv(
> > zi,
> > vzi,
> > random = list(~ MZDZ_factor | Study_ID, ~ MZDZ_factor | ESID), struct="UN",
> > data = data,
> > mods = ~ 0 + MZDZ_factor,
> > method = "REML",
> > tdist = TRUE)
> >
> > you can now just do:
> >
> > deltamethod(model, fun=function(b1,b2) 2*(transf.ztor(b2) -
> > transf.ztor(b1)))
> >
> > and you will directly get the estimate of h^2 and the corresponding SE and CI.
> >
> > Best,
> > Wolfgang
> >
> > > -----Original Message-----
> > > From: R-sig-meta-analysis
> > > <r-sig-meta-analysis-bounces using r-project.org>
> > > On Behalf Of St Pourcain, Beate via R-sig-meta-analysis
> > > Sent: Sunday, October 20, 2024 21:25
> > > To: r-sig-meta-analysis using r-project.org
> > > Cc: St Pourcain, Beate <Beate.StPourcain using mpi.nl>
> > > Subject: Re: [R-meta] R-sig-meta-analysis Digest, Vol 89, Issue 13
> > >
> > > Hi Michael,
> > >
> > > Thanks for pointing this out! The Fisher reference will certainly do.
> > > I had hoped to get the reference for "simulation studies I have
> > > done confirm this", but that's an added bonus.
> > >
> > > Have a nice evening,
> > > Beate
> > >
> > > Date: Sat, 19 Oct 2024 15:16:22 +0100
> > > From: Michael Dewey <lists using dewey.myzen.co.uk>
> > > To: R Special Interest Group for Meta-Analysis
> > > <r-sig-meta-analysis using r-project.org>, "Viechtbauer, Wolfgang (NP)"
> > > <wolfgang.viechtbauer using maastrichtuniversity.nl>
> > > Cc: "St Pourcain, Beate" <Beate.StPourcain using mpi.nl>
> > > Subject: Re: [R-meta] Meta-analysis of intra class correlation
> > > coefficients
> > >
> > > Dear Beate
> > >
> > > Somewhere buried in this thread Wolfgang said
> > >
> > > ==========================
> > > This goes back to Fisher (1925; Statistical methods for research workers).
> > >
> > > In your application (where you dealing with pairs), n is the
> > > number of pairs and m is 2. In that case, you can treat ICC(1)
> > > values like regular correlations. However, if you do apply the
> > > r-to-z transformation, then Fisher suggests to use 1/(n-3/2) as
> > > the variance (instead of 1/(n-3) as we typically use for r-to-z
> > > transformed Pearson product-moment correlation coefficients) and
> > > simulation studies I have done
> > confirm this.
> > > ======================
> > >
> > > Michael
> > >
> > > On 18/10/2024 19:43, St Pourcain, Beate via R-sig-meta-analysis wrote:
> > > > Dear Wolfgang,
> > > > No worries, I am aware of the difference and fully agree with
> > > > your comments. I
> > > was just surprised by the similarity in estimates and had hoped
> > > for an approximation that might provide a quick workaround in the
> > > current
> situation.
> > > Thanks again for all your help, we will take it from here.
> > > >
> > > > In case you would have (at some point) a reference for the
> > > > variance of Z
> > > scores for ICCs as
> > > >
> > > > 1/( n-3/2)
> > > >
> > > > that would be great, no rush!
> > > >
> > > > Enjoy your weekend,
> > > > Beate
More information about the R-sig-meta-analysis
mailing list