[R-meta] Assessing selection bias / multivariate meta-analysis
Pia-Magdalena Schmidt
p|@-m@gd@|en@@@chm|dt @end|ng |rom un|-bonn@de
Tue Dec 3 21:19:57 CET 2024
Dear James,
thank you so much for your detailed reply.
In fact, the negative effect size does not represent an improvement but a
deterioration. However, I would still need to change the code for fitting
the model as the problem of one-sided p-value and direction as described in
your email remains. Unfortunately, I'm not sure which of your suggestions (a
or b) I should go for. I have tried both (results below). Do you have
further recommendations?
Best,
Pia
a)
1.
# 3PSM-bootstrap
set.seed(20240916)
system.time(
mod_3PSM_boot <- selection_model(
data = DatMA,
yi = ES_all_spv$yi*-1,
sei = sei,
cluster = id_database,
selection_type = "step",
steps = .025,
CI_type = "percentile",
bootstrap = "multinomial",
R = 1999
)
)
print(mod_3PSM_boot, transf_gamma = TRUE, transf_zeta = TRUE)
param Est SE percentile_lower percentile_upper
beta 0.702 0.159 0.3011 1.059
tau2 0.325 0.186 0.0893 1.756
lambda1 0.165 0.138 0.0244 0.772
2.
# 3PSM-bootstrap
set.seed(20240916)
system.time(
mod_3PSM_boot <- selection_model(
data = DatMA,
yi = ES_LOR_spv$yi*-1,
sei = sei,
cluster = id_database,
selection_type = "step",
steps = .025,
CI_type = "percentile",
bootstrap = "multinomial",
R = 1999
)
)
print(mod_3PSM_boot, transf_gamma = TRUE, transf_zeta = TRUE)
param Est SE percentile_lower percentile_upper
beta 0.529 0.163 -2.11e-01 0.876
tau2 0.210 0.119 3.07e-02 0.731
lambda1 0.157 0.144 6.00e-17 0.697
b)
1.
system.time(
mod_3PSM_boot <- selection_model(
data = DatMA,
yi = ES_all_spv$yi,
sei = sei,
cluster = id_database,
selection_type = "step",
steps = .975,
CI_type = "percentile",
bootstrap = "multinomial",
R = 1999
)
)
print(mod_3PSM_boot, transf_gamma = TRUE, transf_zeta = TRUE)
param Est SE percentile_lower percentile_upper
beta -0.702 0.159 -1.0594 -0.301
tau2 0.325 0.186 0.0893 1.756
lambda1 6.073 5.092 1.2949 40.908
2.
# 3PSM-bootstrap
set.seed(20240916)
## Startpunkt des Zufallszahlengenerators = wird als Seed bezeichnet.
system.time(
mod_3PSM_boot <- selection_model(
data = DatMA,
yi = ES_LOR_spv$yi,
sei = sei,
cluster = id_database,
selection_type = "step",
steps = .975,
CI_type = "percentile",
bootstrap = "multinomial",
R = 1999
)
)
print(mod_3PSM_boot, transf_gamma = TRUE, transf_zeta = TRUE)
param Est SE percentile_lower percentile_upper
beta -0.529 0.163 -0.8756 2.11e-01
tau2 0.210 0.119 0.0307 7.31e-01
lambda1 6.363 5.848 1.4342 3.18e+15
On So, 1 Dez 2024 10:18:05 -0600
James Pustejovsky <jepusto using gmail.com> wrote:
> Hi Pia,
>
> You can ignore the warning messages that you're getting.
>(We haven't yet
> worked out how to suppress them in the bootstrapping
>code.)
>
> Your code looks good to me except for one subtle but
>potentially
> consequential issue. Based on the multivariate summary
>meta-analyses, it
> looks like you have a strongly negative average effect
>size. If the effects
> are coded so that negative values represent improvement,
>then this needs to
> be taken into account when fitting the selection models.
>The models
> implemented in metaselection are based on one-side
>p-values, where the null
> is mu <= 0 and the alternative is mu > 0 (i.e., positive
>effects are
> improvements). We have not yet implemented an option to
>change the
> direction of the null and alternative hypotheses
>(although it's high on my
> to-do list). In the mean time, there are a few
>alternative ways that you
> could modify the existing code to fit a more appropriate
>model. Either:
> a) Recode effect sizes so that positive values
>correspond to improvement
> (i.e., take yi = -yi) and interpret the beta estimate
>accordingly.
> or
> b) Change the threshold of the step function to steps =
>.975 and interpret
> the lambda (selection parameter) estimate as the
>relative probability that
> a statistically significant (and negative) effect size
>is reported compared
> to the probability that a non-significant or
>counter-therapeutic effect
> size is reported. For instance, if lambda = 4 then this
>means that a
> statistically significant, therapeutic effect is four
>times more likely to
> be reported than a non-significant or non-therapeutic
>effect.
> Again, you only need do (a) or (b) but not both. It's
>possible that making
> this change will shift the results in a meaningful way
>because the revised
> model will be capturing a more plausible form of
>selective reporting.
>
> James
>
> On Thu, Nov 28, 2024 at 11:24 AM Pia-Magdalena Schmidt <
> pia-magdalena.schmidt using uni-bonn.de> wrote:
>
>> Dear James & Wolfgang, dear all,
>>
>>
>> I have tried the metaselection package for some of my
>>meta-analyses and
>> would be very grateful if you could have a look at the
>>results to see if
>> it
>> looks plausible to you?
>>
>>
>> James, you implemented several models, did I understand
>>correctly that you
>> recommend using the modified 3PSM with bootstrapping?
>>
>> Furthermore, I received warnings after fitting the
>>mod_3PSM_boot (see
>> below). May I ignore them?
>>
>>
>> Results from fitted models with rma.mv, cluster robust &
>>mod_3PSM
>> bootstrapping for two analyses:
>>
>>
>> 1.
>> 1.1 Results rma.mv
>> res_spv <- rma.mv(yi=ES_all_spv$yi, V, random = ~ 1 |
>> id_database/effect_id,
>> data = dat)
>>
>>
>> Multivariate Meta-Analysis Model (k = 45; method: REML)
>> logLik Deviance AIC BIC AICc
>> -52.7333 105.4667 111.4667 116.8192 112.0667
>>
>> Variance Components:
>> estim sqrt nlvls fixed factor
>> sigma^2.1 0.0619 0.2488 30 no id_database
>> sigma^2.2 0.2018 0.4493 45 no id_database/effect_id
>>
>> Test for Heterogeneity:
>> Q(df = 44) = 187.0060, p-val < .0001
>>
>> Model Results:
>> estimate se zval pval ci.lb ci.ub
>> -1.0496 0.1073 -9.7845 <.0001 -1.2599 -0.8394 ***
>>
>>
>> 1.2 Results rma.mv cluster robust
>> res_spv.robust <- robust(res_spv, cluster = id_database,
>>clubSandwich =
>> TRUE)
>>
>>
>> Multivariate Meta-Analysis Model (k = 45; method: REML)
>> logLik Deviance AIC BIC AICc
>> -52.7333 105.4667 111.4667 116.8192 112.0667
>>
>> Variance Components:
>> estim sqrt nlvls fixed factor
>> sigma^2.1 0.0619 0.2488 30 no id_database
>> sigma^2.2 0.2018 0.4493 45 no id_database/effect_id
>>
>> Test for Heterogeneity:
>> Q(df = 44) = 187.0060, p-val < .0001
>>
>> Number of estimates: 45
>> Number of clusters: 30
>> Estimates per cluster: 1-3 (mean: 1.50, median: 1)
>>
>> Model Results:
>> estimate se¹ tval¹ df¹ pval¹ ci.lb¹ ci.ub¹
>> -1.0496 0.1071 -9.7997 24.28 <.0001 -1.2706 -0.8287
>>***
>> 1) results based on cluster-robust inference (var-cov
>>estimator: CR2,
>> approx t-test and confidence interval, df:
>>Satterthwaite approx)
>>
>>
>> 1.3 Results mod_3PSM bootstrapping
>> set.seed(20240916)
>> system.time(
>> mod_3PSM_boot <- selection_model(
>> data = dat,
>> yi = ES_all_spv$yi,
>> sei = sei,
>> cluster = id_database,
>> selection_type = "step",
>> steps = .025,
>> CI_type = "percentile",
>> bootstrap = "multinomial",
>> R = 1999
>> )
>> )
>> print(mod_3PSM_boot, transf_gamma = TRUE, transf_zeta =
>>TRUE)
>>
>>
>> user system elapsed
>> 269.861 0.799 270.760
>>
>> Warning messages:
>> 1: In optimx.run(par, optcfg$ufn, optcfg$ugr,
>>optcfg$uhess, lower, :
>> Hessian is reported non-symmetric with asymmetry ratio
>> 1.24952908365386e-11
>> 2: Hessian forced symmetric
>> 3: In optimx.run(par, optcfg$ufn, optcfg$ugr,
>>optcfg$uhess, lower, :
>> Hessian is reported non-symmetric with asymmetry ratio
>> 5.66259792559018e-11
>> 4: Hessian forced symmetric
>> 5: In optimx.check(par, optcfg$ufn, optcfg$ugr,
>>optcfg$uhess, lower, :
>> Parameters or bounds appear to have different
>>scalings.
>> This can cause poor performance in optimization.
>> It is important for derivative free methods like
>>BOBYQA, UOBYQA, NEWUOA.
>> 6: In optimx.check(par, optcfg$ufn, optcfg$ugr,
>>optcfg$uhess, lower, :
>> Parameters or bounds appear to have different
>>scalings.
>> This can cause poor performance in optimization.
>> It is important for derivative free methods like
>>BOBYQA, UOBYQA, NEWUOA.
>>
>>
>> > print(mod_3PSM_boot, transf_gamma = TRUE, transf_zeta
>>= TRUE)
>> param Est SE percentile_lower percentile_upper
>> beta -1.06e+00 0.1703 -1.63e+00 -8.09e-01
>> tau2 3.07e-01 0.1781 9.40e-02 2.20e+00
>> lambda1 8.07e+13 0.0604 1.13e+12 4.87e+14
>>
>>
>> 2.
>> 2.1 Results rma.mv
>> res_LOR <- rma.mv(yi=ES_LOR_spv$yi, V, random = ~ 1 |
>> id_database/effect_id,
>> data = dat)
>>
>>
>> Multivariate Meta-Analysis Model (k = 19; method: REML)
>> logLik Deviance AIC BIC AICc
>> -14.8934 29.7867 35.7867 38.4579 37.5010
>>
>> Variance Components:
>> estim sqrt nlvls fixed factor
>> sigma^2.1 0.0000 0.0000 13 no id_database
>> sigma^2.2 0.1889 0.4347 19 no id_database/effect_id
>>
>> Test for Heterogeneity:
>> Q(df = 18) = 88.4226, p-val < .0001
>>
>> Model Results:
>> estimate se zval pval ci.lb ci.ub
>> -0.8414 0.1257 -6.6954 <.0001 -1.0877 -0.5951 ***
>>
>>
>> 2.2 Results rma.mv cluster robust
>> res_LOR.robust <- robust(res_LOR, cluster = id_database,
>>clubSandwich =
>> TRUE)
>>
>>
>> Multivariate Meta-Analysis Model (k = 19; method: REML)
>>
>> logLik Deviance AIC BIC AICc
>> -14.8934 29.7867 35.7867 38.4579 37.5010
>>
>> Variance Components:
>> estim sqrt nlvls fixed factor
>> sigma^2.1 0.0000 0.0000 13 no id_database
>> sigma^2.2 0.1889 0.4347 19 no id_database/effect_id
>>
>> Test for Heterogeneity:
>> Q(df = 18) = 88.4226, p-val < .0001
>>
>> Number of estimates: 19
>> Number of clusters: 13
>> Estimates per cluster: 1-3 (mean: 1.46, median: 1)
>>
>> Model Results:
>> estimate se¹ tval¹ df¹ pval¹ ci.lb¹ ci.ub¹
>> -0.8414 0.1178 -7.1438 9.08 <.0001 -1.1075 -0.5754 ***
>> 1) results based on cluster-robust inference (var-cov
>>estimator: CR2,
>> approx t-test and confidence interval, df:
>>Satterthwaite approx)
>>
>>
>> 2.3 Results mod_3PSM bootstrapping
>> set.seed(20240916)
>> system.time(
>> mod_3PSM_boot <- selection_model(
>> data = dat,
>> yi = ES_LOR_spv$yi,
>> sei = sei,
>> cluster = id_database,
>> selection_type = "step",
>> steps = .025,
>> CI_type = "percentile",
>> bootstrap = "multinomial",
>> R = 1999
>> )
>> )
>> print(mod_3PSM_boot, transf_gamma = TRUE, transf_zeta =
>>TRUE)
>>
>> user system elapsed
>> 262.519 0.324 262.812
>> Warning messages:
>> 1: In optimx.run(par, optcfg$ufn, optcfg$ugr,
>>optcfg$uhess, lower, :
>> Hessian is reported non-symmetric with asymmetry ratio
>> 9.13585399793946e-13
>> 2: Hessian forced symmetric
>> 3: In optimx.run(par, optcfg$ufn, optcfg$ugr,
>>optcfg$uhess, lower, :
>> Hessian is reported non-symmetric with asymmetry ratio
>> 1.06505934682683e-11
>> 4: Hessian forced symmetric
>> >
>> > print(mod_3PSM_boot, transf_gamma = TRUE, transf_zeta
>>= TRUE)
>> param Est SE percentile_lower percentile_upper
>> beta -8.13e-01 0.1227 -1.14e+00 -6.07e-01
>> tau2 1.56e-01 0.0892 1.87e-02 3.67e-01
>> lambda1 1.12e+13 0.0566 1.36e+11 2.31e+14
>>
>>
>> Best,
>> Pia
>>
>> On Do, 21 Nov 2024 17:18:28 +0000
>> Viechtbauer, Wolfgang (NP) via R-sig-meta-analysis
>> <r-sig-meta-analysis using r-project.org> wrote:
>> > Thanks for providing additional info about
>> >metaselection. I think this package will be a very
>> >important tool in the meta-analytic toolbox!
>> >
>> > I was inspired by your recent blog post:
>> >
>> >
>>https://jepusto.com/posts/beta-density-selection-models/
>> >
>> > and added the truncated beta selection model also to
>> >selmodel(). I took a very quick peak at your code and I
>> >think you are using analytic gradients, which helps to
>> >speed up / stabilize the model fitting. But it's nice
>>to
>> >have two parallel implementations to cross-check the
>> >results.
>> >
>> > Best,
>> > Wolfgang
>> >
>> >> -----Original Message-----
>> >> From: James Pustejovsky <jepusto using gmail.com>
>> >> Sent: Thursday, November 21, 2024 15:10
>> >> To: R Special Interest Group for Meta-Analysis
>> >><r-sig-meta-analysis using r-
>> >> project.org>
>> >> Cc: Viechtbauer, Wolfgang (NP)
>> >><wolfgang.viechtbauer using maastrichtuniversity.nl>
>> >> Subject: Re: [R-meta] Assessing selection bias /
>> >>multivariate meta-analysis
>> >>
>> >> Was going to chime in about the metaselection package
>> >> (https://github.com/jepusto/metaselection)---it's
>>still
>> >>under development but
>> >> the core functionality and documentation is in place.
>> >>The package implements the
>> >> bootstrapped selection model as demonstrated in my
>>blog
>> >>post
>> >>
>>(https://jepusto.com/posts/cluster-bootstrap-selection-model/),
>> >>but with a much
>> >> easier interface and faster calculation; it also
>> >>implements selection models
>> >> with cluster-robust standard errors, though these
>>seem
>> >>to be not as accurate as
>> >> bootstrapping. Folks are welcome to give the package
>>a
>> >>try and to reach out
>> >> with questions or potential bugs if you run into
>> >>anything. We are working on a
>> >> paper describing the methods implemented in the
>>package
>> >>and reporting pretty
>> >> extensive simulations about their performance.
>> >>
>> >> My student Man Chen (now on the faculty at UT Austin)
>> >>has studied a whole bunch
>> >> of the available methods for selective reporting bias
>> >>correction, looking
>> >> specifically at how they perform in meta-analyses
>>with
>> >>dependent effect sizes,
>> >> and proposing adaptations of some of the methods to
>> >>better acknowledge
>> >> dependency. Our working paper (on this is
>> >> here: https://osf.io/preprints/metaarxiv/jq52s
>> >>
>> >> Pia asked about a few other possible techniques:
>> >> - The Egger's test / PET-PEESE approach with
>> >>cluster-robust variance estimation
>> >> is reasonable but, as Wolfgang noted, it is not
>> >>specifically diagnostic about
>> >> missing studies vs. missing effects. If the effect
>>sizes
>> >>nested within a given
>> >> study tend to have similar standard errors, then it
>>will
>> >>mostly be picking up on
>> >> association between study sample size and study-level
>> >>average effect size. And
>> >> of course, it also has the limitation that this
>> >>small-study association can be
>> >> caused by things other than selective reporting.
>> >> - Mathur & Vanderweele's sensitivity analysis is
>>quite
>> >>useful, though it does
>> >> not provide an estimate of the severity of selective
>> >>reporting (instead, it
>> >> provides information about the degree of potential
>>bias
>> >>assuming a specific
>> >> level of selection).
>> >> - For 3PSM, the cluster-bootstrap technique
>>implemented
>> >>in the metaselection
>> >> package is a way to deal with dependent effects, so
>>it
>> >>is no longer necessary to
>> >> use ad hoc approaches like ignoring dependence,
>> >>aggregating to the study level,
>> >> or selecting a single effect per study.
>> >>
>> >> James
>> >>
>> >> On Thu, Nov 21, 2024 at 6:37 AM Viechtbauer, Wolfgang
>> >>(NP) via R-sig-meta-
>> >> analysis <mailto:r-sig-meta-analysis using r-project.org>
>> >>wrote:
>> >> And I just stumbled across this:
>> >>
>> >> https://github.com/jepusto/metaselection
>> >>
>> >> James, don't hide all your good work from us!
>> >>
>> >> 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: Thursday, November 21, 2024 13:21
>> >> > To: R Special Interest Group for Meta-Analysis
>> >><r-sig-meta-analysis using r-
>> >> > http://project.org>
>> >> > Cc: Viechtbauer, Wolfgang (NP)
>> >> <mailto:wolfgang.viechtbauer using maastrichtuniversity.nl>
>> >> > Subject: Re: [R-meta] Assessing selection bias /
>> >>multivariate meta-analysis
>> >> >
>> >> > Dear Pia,
>> >> >
>> >> > Generally, I don't think there really is any method
>> >>that is going to be a
>> >> great
>> >> > choice here. The 'Egger sandwich' (i.e., an Egger
>>type
>> >>regression model using
>> >> > cluster-robust inference methods) is a decent
>>option,
>> >>since it logically
>> >> > generalizes the standard Egger regression method to
>> >>this context, but it is
>> >> > unclear what kind of bias/selection effect this may
>> >>pick up (missing studies,
>> >> > missing estimates within studies, a combination
>> >>thereof).
>> >> >
>> >> > Yes, for the 3PSM, you would have to either ignore
>>the
>> >>dependencies or select
>> >> > one estimate per study (and maybe repeat the latter
>>a
>> >>large number of times
>> >> for
>> >> > different subsets).
>> >> >
>> >> > I assume you are familiar with these papers. If
>>not,
>> >>they are directly
>> >> relevant:
>> >> >
>> >> > Rodgers, M. A., & Pustejovsky, J. E. (2021).
>> >>Evaluating meta-analytic methods
>> >> to
>> >> > detect selective reporting in the presence of
>> >>dependent effect sizes.
>> >> > Psychological Methods, 26(2), 141-160.
>> >>https://doi.org/10.1037/met0000300
>> >> >
>> >> > Fernández-Castilla, B., Declercq, L., Jamshidi, L.,
>> >>Beretvas, S. N., Onghena,
>> >> > P., & Van den Noortgate, W. (2021). Detecting
>> >>selection bias in meta-analyses
>> >> > with multiple outcomes: A simulation study. The
>> >>Journal of Experimental
>> >> > Education, 89(1), 125-144.
>> >>https://doi.org/10.1080/00220973.2019.1582470
>> >> >
>> >> > Nakagawa, S., Lagisz, M., Jennions, M. D.,
>>Koricheva,
>> >>J., Noble, D. W. A.,
>> >> > Parker, T. H., Sánchez-Tójar, A., Yang, Y., &
>>O'Dea,
>> >>R. E. (2022). Methods for
>> >> > testing publication bias in ecological and
>> >>evolutionary meta-analyses. Methods
>> >> > in Ecology and Evolution, 13(1), 4-21.
>> >>https://doi.org/10.1111/2041-210X.13724
>> >> >
>> >> > I think James is working on some methods related to
>> >>this topic:
>> >> >
>> >> >
>> >>https://jepusto.com/posts/cluster-bootstrap-selection-model/
>> >> >
>> >> > Best,
>> >> > Wolfgang
>> >> >
>> >> > > -----Original Message-----
>> >> > > From: R-sig-meta-analysis
>> >><mailto:r-sig-meta-analysis-bounces using r-project.org>
>> >> On
>> >> > Behalf
>> >> > > Of Pia-Magdalena Schmidt via R-sig-meta-analysis
>> >> > > Sent: Wednesday, November 20, 2024 21:58
>> >> > > To: mailto:r-sig-meta-analysis using r-project.org
>> >> > > Cc: Pia-Magdalena Schmidt
>> >><mailto:pia-magdalena.schmidt using uni-bonn.de>
>> >> > > Subject: [R-meta] Assessing selection bias /
>> >>multivariate meta-analysis
>> >> > >
>> >> > > Dear all,
>> >> > > Although this topic has been discussed several
>>times
>> >>and I read the archives
>> >> > > and referenced papers, I’m still not sure how to
>> >>assess and possibly correct
>> >> > > for selection bias in multivariate meta-analyses.
>> >> > >
>> >> > > I used the metafor package and ran meta-analyses
>> >>with SMCC as effect size
>> >> > > (all studies used within-designs) and fitted
>> >>http://rma.mv models as several
>> >> > > studies report more than one effect size.
>> >>Furthermore, I used cluster-robust
>> >> > > methods to examine the robustness of the models.
>> >> > > For a subset of my data, I used meta-regressions
>> >>with one continuous
>> >> > > moderator.
>> >> > > All effect sizes are from published journal
>> >>articles. The range of included
>> >> > > studies is between 30 and 6 with a number of
>>effect
>> >>sizes between 45 and 10.
>> >> > >
>> >> > > Since I want to take the dependencies into
>>account,
>> >>I would not use funnel
>> >> > > plots or trim and fill. I wonder if using Egger's
>> >>regression test adjusted
>> >> > > for http://rma.mv as well as PET-PEESE and
>>perhaps
>> >>the sensitivity analysis
>> >> > > suggested by Mathur & Vanderweele (2020) as well
>>as
>> >>3PSM would be a
>> >> > > reasonable way to go? Although the latter would
>>only
>> >>use one effect size per
>> >> > > study or an aggregated effect size, right?
>> >> > >
>> >> > > I would be very grateful for any recommendations!
>> >> > > Best,
>> >> > > Pia
>> >> > >
>> >> > > Below is an excerpt from my code:
>> >> > > ES_all <- escalc(measure="SMCC", m1i= m1i, sd1i=
>> >>sd1i, ni = ni, m2i= m2i,
>> >> > > sd2i= sd2i, pi= pi, ri = ri, data= dat)
>> >> > > V <- vcalc(vi=ES_all$vi, cluster=id_database, obs
>>=
>> >>effect_id, rho =0.605,
>> >> > > data=dat)
>> >> > > res <- http://rma.mv(yi=ES_all$yi, V, random = ~
>>1 |
>> >>id_database/effect_id,
>> >> data =
>> >> > > dat)
>> >> > > res.robust <- robust(res, cluster = id_database,
>> >>clubSandwich = TRUE)
>> >> > >
>> >> > > # subset
>> >> > > res_LOR <- http://rma.mv(yi=ES_LOR$yi, V, random
>>= ~
>> >>1 |
>> >> id_database/effect_id,
>> >> > > mods = ~ dose, data = dat)
>> > _______________________________________________
>> > 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
>>
>>
More information about the R-sig-meta-analysis
mailing list