# Random Effects Model - forest, funnel, summary (interpret negative symbol as positive)
meta <- escalc(measure="SMCR",m1i = mydata$i.wb.m.t1, m2i=mydata$i.wb.m.t2, sd1i= mydata$i.wb.sd.t1, sd2i = mydata$i.wb.sd.t2, ni= mydata$i.sample.size.t1, ri=ri, data=mydata)
# "SMCR" for the standardized mean change using raw score standardization (Becker, 1988)
mc <- rma.uni(yi,vi, data=meta, slab = meta$study.id)
summary(mc)
##
## Random-Effects Model (k = 12; tau^2 estimator: REML)
##
## logLik deviance AIC BIC AICc
## -12.9434 25.8869 29.8869 30.6826 31.3869
##
## tau^2 (estimated amount of total heterogeneity): 0.5194 (SE = 0.2253)
## tau (square root of estimated tau^2 value): 0.7207
## I^2 (total heterogeneity / total variability): 99.80%
## H^2 (total variability / sampling variability): 502.71
##
## Test for Heterogeneity:
## Q(df = 11) = 617.5099, p-val < .0001
##
## Model Results:
##
## estimate se zval pval ci.lb ci.ub
## -0.5281 0.2099 -2.5160 0.0119 -0.9395 -0.1167 *
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
confint(mc)
##
## estimate ci.lb ci.ub
## tau^2 0.5194 0.2627 1.9334
## tau 0.7207 0.5126 1.3905
## I^2(%) 99.8011 99.6075 99.9465
## H^2 502.7144 254.7680 1868.4379
# PPIs improve wellbeing by about .53 on average (pvalue = .01), 95% CI = -94, -.12
#Creating vectors to keep figures consistent (this is for forest)
xl=c(-20, 20)
al=c(-10,10)
cx=1.5
forest(mc, alim = xl, xlim=al,header=TRUE)
text(-4.5, -1, pos=1, cex=.75, bquote(paste(
"(estimate = -.53",
", se = ", .(mc$se), ", ",
.(fmtp2(mc$pval)), "; ",
I^2, " = ", .(fmtx(mc$I2, digits=1)), "%, ",
tau^2, " = ", .(fmtx(mc$tau2, digits=2)), ")")))

#Creating vectors to keep figures consistent (this is for funnel)
xl=c(-4, 4)
al=c(-10,10)
cx=0.75
funnel(mc, xlim = xl, cex.lab=cx,cex.axis=cx,back='White',offset=1.5,shade="lightgrey",hlines="lightgrey", cex=cx, slab=mydata$study.id, label=FALSE, legend="topleft")

# categorical moderator
# frequency
m3 <- rma(yi, vi, data=meta, slab = meta$study.id, mods = ~ frequency)
summary(m3)
##
## Mixed-Effects Model (k = 12; tau^2 estimator: REML)
##
## logLik deviance AIC BIC AICc
## -8.8069 17.6139 29.6139 29.2893 113.6139
##
## tau^2 (estimated amount of residual heterogeneity): 0.6724 (SE = 0.3642)
## tau (square root of estimated tau^2 value): 0.8200
## I^2 (residual heterogeneity / unaccounted variability): 99.74%
## H^2 (unaccounted variability / sampling variability): 378.43
## R^2 (amount of heterogeneity accounted for): 0.00%
##
## Test for Residual Heterogeneity:
## QE(df = 7) = 213.5241, p-val < .0001
##
## Test of Moderators (coefficients 2:5):
## QM(df = 4) = 2.6930, p-val = 0.6104
##
## Model Results:
##
## estimate se zval pval ci.lb
## intrcpt -0.4061 0.4108 -0.9887 0.3228 -1.2112
## frequencyMore than once per week -1.0202 0.7276 -1.4023 0.1608 -2.4462
## frequencyOnce 0.1830 0.7115 0.2573 0.7970 -1.2115
## frequencyVaried -0.0239 0.9176 -0.0261 0.9792 -1.8224
## frequencyWeekly 0.0108 0.6273 0.0172 0.9863 -1.2188
## ci.ub
## intrcpt 0.3990
## frequencyMore than once per week 0.4057
## frequencyOnce 1.5776
## frequencyVaried 1.7745
## frequencyWeekly 1.2403
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
confint(m3)
##
## estimate ci.lb ci.ub
## tau^2 0.6724 0.2848 3.1303
## tau 0.8200 0.5336 1.7693
## I^2(%) 99.7357 99.3782 99.9431
## H^2 378.4262 160.8316 1758.0524
# intercept reflects the estimated average change for daily activity use (this is the reference level)
# the other two coefficients estimate how much the average change for remaining frequencies differ from the reference level (daily)
# not significant, so frequency doesn't shape the effectiveness of the PPI
#continent
m4 <- rma(yi, vi, data=meta, slab = meta$study.id, mods = ~ continent)
summary(m4)
##
## Mixed-Effects Model (k = 12; tau^2 estimator: REML)
##
## logLik deviance AIC BIC AICc
## -9.4584 18.9168 28.9168 29.3140 58.9168
##
## tau^2 (estimated amount of residual heterogeneity): 0.5547 (SE = 0.2819)
## tau (square root of estimated tau^2 value): 0.7448
## I^2 (residual heterogeneity / unaccounted variability): 99.84%
## H^2 (unaccounted variability / sampling variability): 607.68
## R^2 (amount of heterogeneity accounted for): 0.00%
##
## Test for Residual Heterogeneity:
## QE(df = 8) = 543.8498, p-val < .0001
##
## Test of Moderators (coefficients 2:4):
## QM(df = 3) = 3.1864, p-val = 0.3638
##
## Model Results:
##
## estimate se zval pval ci.lb ci.ub
## intrcpt -1.2110 0.4422 -2.7383 0.0062 -2.0778 -0.3442 **
## continentAustralia 0.8496 0.8691 0.9776 0.3283 -0.8538 2.5529
## continentEurope 0.9269 0.5245 1.7672 0.0772 -0.1011 1.9548 .
## continentSouth America 0.7289 0.8672 0.8406 0.4006 -0.9707 2.4285
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
confint(m4)
##
## estimate ci.lb ci.ub
## tau^2 0.5547 0.2458 2.4342
## tau 0.7448 0.4958 1.5602
## I^2(%) 99.8354 99.6295 99.9625
## H^2 607.6848 269.8888 2663.4393
#no difference in wellbeing improvement by continent
#sample type
m5 <- rma(yi, vi, data=meta, slab = meta$study.id, mods = ~ sample.type)
summary(m5)
##
## Mixed-Effects Model (k = 12; tau^2 estimator: REML)
##
## logLik deviance AIC BIC AICc
## -11.3777 22.7554 30.7554 31.5443 40.7554
##
## tau^2 (estimated amount of residual heterogeneity): 0.6459 (SE = 0.3093)
## tau (square root of estimated tau^2 value): 0.8037
## I^2 (residual heterogeneity / unaccounted variability): 99.81%
## H^2 (unaccounted variability / sampling variability): 527.83
## R^2 (amount of heterogeneity accounted for): 0.00%
##
## Test for Residual Heterogeneity:
## QE(df = 9) = 374.8240, p-val < .0001
##
## Test of Moderators (coefficients 2:3):
## QM(df = 2) = 0.4931, p-val = 0.7815
##
## Model Results:
##
## estimate se zval pval ci.lb ci.ub
## intrcpt -0.6839 0.3323 -2.0577 0.0396 -1.3353 -0.0325 *
## sample.typeOlder adults 0.4370 0.6587 0.6634 0.5071 -0.8541 1.7281
## sample.typeStudents 0.2250 0.5222 0.4309 0.6665 -0.7985 1.2486
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
confint(m5)
##
## estimate ci.lb ci.ub
## tau^2 0.6459 0.3002 2.6164
## tau 0.8037 0.5479 1.6175
## I^2(%) 99.8105 99.5932 99.9532
## H^2 527.8295 245.8272 2135.0167
describeBy(meta$yi, meta$sample.type)
##
## Descriptive statistics by group
## group: General/Community
## vars n mean sd median trimmed mad min max range skew kurtosis se
## X1 1 6 -0.74 1.2 -0.4 -0.74 0.29 -3.16 0.08 3.24 -1.26 -0.24 0.49
## ------------------------------------------------------------
## group: Older adults
## vars n mean sd median trimmed mad min max range skew kurtosis se
## X1 1 2 -0.25 0.13 -0.25 -0.25 0.14 -0.34 -0.15 0.19 0 -2.75 0.09
## ------------------------------------------------------------
## group: Students
## vars n mean sd median trimmed mad min max range skew kurtosis se
## X1 1 4 -0.46 0.08 -0.47 -0.46 0.08 -0.54 -0.36 0.18 0.15 -2.05 0.04
#wellbeing improvement strongest for general/community sample (compared to older adults and students)
#age
m6 <- rma(yi, vi, data=meta, slab = meta$study.id, mods = ~ age.mean)
summary(m6)
##
## Mixed-Effects Model (k = 12; tau^2 estimator: REML)
##
## logLik deviance AIC BIC AICc
## -12.2825 24.5649 30.5649 31.4727 34.5649
##
## tau^2 (estimated amount of residual heterogeneity): 0.5873 (SE = 0.2670)
## tau (square root of estimated tau^2 value): 0.7664
## I^2 (residual heterogeneity / unaccounted variability): 99.82%
## H^2 (unaccounted variability / sampling variability): 554.20
## R^2 (amount of heterogeneity accounted for): 0.00%
##
## Test for Residual Heterogeneity:
## QE(df = 10) = 617.1709, p-val < .0001
##
## Test of Moderators (coefficient 2):
## QM(df = 1) = 0.0710, p-val = 0.7899
##
## Model Results:
##
## estimate se zval pval ci.lb ci.ub
## intrcpt -0.6693 0.5623 -1.1902 0.2340 -1.7714 0.4329
## age.mean 0.0037 0.0138 0.2664 0.7899 -0.0233 0.0307
##
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
confint(m6)
##
## estimate ci.lb ci.ub
## tau^2 0.5873 0.2851 2.2725
## tau 0.7664 0.5339 1.5075
## I^2(%) 99.8196 99.6290 99.9533
## H^2 554.2021 269.5299 2141.6029
#no age differences
## RESOURCES ##
# https://training.cochrane.org/handbook/current/chapter-15
# estimate of the average true standardized mean difference is .53
# Q-test suggests that the underlying true standardized mean differences are heterogeneous
# can put weighted=false in the rma.uni code, which would give us an estimate reflecting unweighted average of
# the true standardized mean differences of these studies
# for continuous outcome measures, we can present summary results for studies using natural units of measurement
# or as minimal important differences when all studies use the same scale (e.g. SWL)
# our studies measure the same construct but with different scales, so need to find a way to interpret the
# standardized mean difference
# resources for interpreting figures
# https://s4be.cochrane.org/blog/2023/06/27/how-to-read-a-funnel-plot/
# https://www.google.com/url?sa=i&url=https%3A%2F%2Fs4be.cochrane.org%2Fblog%2F2016%2F07%2F11%2Ftutorial-read-forest-plot%2F&psig=AOvVaw3QnJDUY9cHQ1Gc4fa1jKN0&ust=1739462356615000&source=images&cd=vfe&opi=89978449&ved=0CBQQjRxqFwoTCMCbzrHAvosDFQAAAAAdAAAAABAE
# https://www.google.com/url?sa=i&url=https%3A%2F%2Fbookdown.org%2FMathiasHarrer%2FDoing_Meta_Analysis_in_R%2Fforest.html&psig=AOvVaw3QnJDUY9cHQ1Gc4fa1jKN0&ust=1739462356615000&source=images&cd=vfe&opi=89978449&ved=0CBQQjRxqFwoTCMCbzrHAvosDFQAAAAAdAAAAABAJ
### forest plot
# line is of a null effect (i.e., 0 - no change in wellbeing)
# black box is the point estimate of the effect size (smc)
# any black box that crosses the null effect line is n.s. (none in this case)
# diamond (at R.E.M. line) is the point estimate and CIs when you combine and
# average all individual studies together
# I2 tells us heterogeneity of the studies (how consistent they are)
# reflects variability of true effect sizes
# over 50% means studies may be inconsistent due to a reason other than chance
# tau2 is the estimated standard deviation of the distribution of effect sizes
# quantifies proportion of total variation attributable to heterogenity
### funnel plot
# 95% confidence interval is the triangle lines
# y axis is study precision (higher up the dot, higher precision)
# x axis is the estimated effect size (standardized mean change)
# 95% of studies would be expected to be within the triangle but because the
# rule of thumb is to only interpret a funnel plot's asymmetry with 10 or more studies