Short: this can be done to display each forest plot, but saving them programmatically will be difficult.
library(metafor)
#> Loading required package: Matrix
#> Loading 'metafor' package (version 2.4-0). For an overview
#> and introduction to the package please type: help(metafor).
#> Loading required package: Matrix
#> Loading 'metafor' package (version 2.4-0). For an overview
#> and introduction to the package please type: help(metafor).
dat <- data.frame(
StudyID = c(
"Australia", "Australia",
"Australia", "Australia", "Australia", "Australia",
"Australia", "Australia", "Australia", "Australia",
"Australia", "Australia", "Australia", "Australia",
"Australia", "Australia", "Newzealand", "Newzealand",
"Newzealand", "Newzealand", "Newzealand", "Newzealand",
"Newzealand", "Newzealand", "Newzealand", "Newzealand",
"Newzealand", "Newzealand", "Newzealand", "Newzealand",
"Newzealand", "Newzealand"
), variable_control = c(
"Mathew",
"Mathew", "Mathew", "Mathew", "Simp", "Simp", "Simp", "Simp",
"hutten", "hutten", "hutten", "hutten", "AB",
"AB", "AB", "AB", "Mathew", "Mathew", "Mathew", "Mathew", "Simp",
"Simp", "Simp", "Simp", "hutten", "hutten",
"hutten", "hutten", "AB", "AB", "AB", "AB"
),
mean_control = c(
103.125,
103.125, 103.125, 103.125, 3.20026104175422, 3.20026104175422,
3.20026104175422, 3.20026104175422, 0.655931978647405, 0.655931978647405,
0.655931978647405, 0.655931978647405, 5.6147916190302, 5.6147916190302,
5.6147916190302, 5.6147916190302, 86.3478260869565, 86.3478260869565,
86.3478260869565, 86.3478260869565, 2.70600451720365, 2.70600451720365,
2.70600451720365, 2.70600451720365, 0.776713622717005, 0.776713622717005,
0.776713622717005, 0.776713622717005, 2.25282841071637, 2.25282841071637,
2.25282841071637, 2.25282841071637
), sd_control = c(
14.0136837888306,
14.0136837888306, 14.0136837888306, 14.0136837888306, 0.278705816120361,
0.278705816120361, 0.278705816120361, 0.278705816120361, 0.0899465623045049,
0.0899465623045049, 0.0899465623045049, 0.0899465623045049, 6.28461426690684,
6.28461426690684, 6.28461426690684, 6.28461426690684, 22.9900629657485,
22.9900629657485, 22.9900629657485, 22.9900629657485, 0.492994678382565,
0.492994678382565, 0.492994678382565, 0.492994678382565, 0.129978730202719,
0.129978730202719, 0.129978730202719, 0.129978730202719, 7.7819081791642,
7.7819081791642, 7.7819081791642, 7.7819081791642
), n_control = c(
16L,
16L, 16L, 16L, 16L, 16L, 16L, 16L, 120L, 120L, 120L, 120L, 15L,
15L, 15L, 15L, 46L, 46L, 46L, 46L, 46L, 46L, 46L, 46L, 1035L,
1035L, 1035L, 1035L, 45L, 45L, 45L, 45L
),
se_control = c(
3.50342094720765,
3.50342094720765, 3.50342094720765, 3.50342094720765, 0.0696764540300902,
0.0696764540300902, 0.0696764540300902, 0.0696764540300902, 0.00821096019070354, 0.00821096019070354, 0.00821096019070354, 0.00821096019070354,
1.62268042620452, 1.62268042620452, 1.62268042620452, 1.62268042620452,
3.38969985579422, 3.38969985579422, 3.38969985579422, 3.38969985579422,
0.0726880997546798, 0.0726880997546798, 0.0726880997546798, 0.0726880997546798,
0.00404019302943354, 0.00404019302943354, 0.00404019302943354,
0.00404019302943354, 1.16005837888485, 1.16005837888485, 1.16005837888485,
1.16005837888485
),
lowCI_control = c(
95.6576350141697, 95.6576350141697,
95.6576350141697, 95.6576350141697, 3.05174919547557, 3.05174919547557,
3.05174919547557, 3.05174919547557, 0.6396734573882, 0.6396734573882,
0.6396734573882, 0.6396734573882, 2.13448824216196, 2.13448824216196,
2.13448824216196, 2.13448824216196, 79.5206201201125, 79.5206201201125,
79.5206201201125, 79.5206201201125, 2.55960316915644, 2.55960316915644,
2.55960316915644, 2.55960316915644, 0.768785709935317, 0.768785709935317,
0.768785709935317, 0.768785709935317, -0.0851156305499581, -0.0851156305499581,
-0.0851156305499581, -0.0851156305499581
), hiCI_control = c(
110.59236498583,
110.59236498583, 110.59236498583, 110.59236498583, 3.34877288803287,
3.34877288803287, 3.34877288803287, 3.34877288803287, 0.67219049990661,
0.67219049990661, 0.67219049990661, 0.67219049990661, 9.09509499589843,
9.09509499589843, 9.09509499589843, 9.09509499589843, 93.1750320538006,
93.1750320538006, 93.1750320538006, 93.1750320538006, 2.85240586525086,
2.85240586525086, 2.85240586525086, 2.85240586525086, 0.784641535498692,
0.784641535498692, 0.784641535498692, 0.784641535498692, 4.59077245198269,
4.59077245198269, 4.59077245198269, 4.59077245198269
), Variable_test = c(
"Mathew",
"Simp", "hutten", "AB", "Mathew", "Simp", "hutten",
"AB", "Mathew", "Simp", "hutten", "AB", "Mathew", "Simp",
"hutten", "AB", "Mathew", "Simp", "hutten", "AB",
"Mathew", "Simp", "hutten", "AB", "Mathew", "Simp",
"hutten", "AB", "Mathew", "Simp", "hutten", "AB"
),
mean_test = c(
94.85, 3.11469350939242, 0.638430601983837,
6.91683472226622, 94.85, 3.11469350939242, 0.638430601983837,
6.91683472226622, 94.85, 3.11469350939242, 0.638430601983837,
6.91683472226622, 94.85, 3.11469350939242, 0.638430601983837,
6.91683472226622, 86.0222222222222, 2.6390595091832, 0.759731493048648,
0.73921233900569, 86.0222222222222, 2.6390595091832, 0.759731493048648,
0.73921233900569, 86.0222222222222, 2.6390595091832, 0.759731493048648,
0.73921233900569, 86.0222222222222, 2.6390595091832, 0.759731493048648,
0.73921233900569
),
sd_test = c(
17.2909316055982, 0.202930255027022,
0.0919447770686176, 8.14939209465999, 17.2909316055982, 0.202930255027022,
0.0919447770686176, 8.14939209465999, 17.2909316055982, 0.202930255027022,
0.0919447770686176, 8.14939209465999, 17.2909316055982, 0.202930255027022,
0.0919447770686176, 8.14939209465999, 21.4163158796882, 0.451677388306635,
0.138197141547272, 0.908654293596164, 21.4163158796882, 0.451677388306635,
0.138197141547272, 0.908654293596164, 21.4163158796882, 0.451677388306635,
0.138197141547272, 0.908654293596164, 21.4163158796882, 0.451677388306635,
0.138197141547272, 0.908654293596164
),
n_test = c(
20L, 20L,
190L, 21L, 20L, 20L, 190L, 21L, 20L, 20L, 190L, 21L, 20L,
20L, 190L, 21L, 45L, 45L, 990L, 45L, 45L, 45L, 990L, 45L,
45L, 45L, 990L, 45L, 45L, 45L, 990L, 45L
),
se_test = c(
3.86636984644171,
0.0453765844931791, 0.00667037520849417, 1.77834314960258,
3.86636984644171, 0.0453765844931791, 0.00667037520849417,
1.77834314960258, 3.86636984644171, 0.0453765844931791, 0.00667037520849417,
1.77834314960258, 3.86636984644171, 0.0453765844931791, 0.00667037520849417,
1.77834314960258, 3.1925558756394, 0.0673320896102137, 0.00439219348020005,
0.135454184568538, 3.1925558756394, 0.0673320896102137, 0.00439219348020005,
0.135454184568538, 3.1925558756394, 0.0673320896102137, 0.00439219348020005,
0.135454184568538, 3.1925558756394, 0.0673320896102137, 0.00439219348020005,
0.135454184568538
), lowCI_test = c(
86.7575949081585, 3.01971922654131,
0.625272652671813, 3.20727591549959, 86.7575949081585, 3.01971922654131,
0.625272652671813, 3.20727591549959, 86.7575949081585, 3.01971922654131,
0.625272652671813, 3.20727591549959, 86.7575949081585, 3.01971922654131,
0.625272652671813, 3.20727591549959, 79.5880486308587, 2.50336059906323,
0.751112403965226, 0.466222367603537, 79.5880486308587, 2.50336059906323,
0.751112403965226, 0.466222367603537, 79.5880486308587, 2.50336059906323,
0.751112403965226, 0.466222367603537, 79.5880486308587, 2.50336059906323,
0.751112403965226, 0.466222367603537
),
hiCI_test = c(
102.942405091841,
3.20966779224354, 0.651588551295861, 10.6263935290328, 102.942405091841,
3.20966779224354, 0.651588551295861, 10.6263935290328, 102.942405091841,
3.20966779224354, 0.651588551295861, 10.6263935290328, 102.942405091841,
3.20966779224354, 0.651588551295861, 10.6263935290328, 92.4563958135858,
2.77475841930317, 0.768350582132071, 1.01220231040784, 92.4563958135858,
2.77475841930317, 0.768350582132071, 1.01220231040784, 92.4563958135858,
2.77475841930317, 0.768350582132071, 1.01220231040784, 92.4563958135858,
2.77475841930317, 0.768350582132071, 1.01220231040784
)
)
run_escalc <- function(x){
escalc(n1i = n_control,
n2i = n_test,
m1i = mean_control,
m2i = mean_test,
sd1i = sd_control,
sd2i = sd_test,
data = subset(dat,variable_control == x),
measure = "SMD",
append = TRUE)
}
the_vars <- unique(dat[,"variable_control"])
# this returns plot.default objects that print to screen only
Map(run_escalc,the_vars) -> the_list
myfunction <- function(x) {
forest.rma(rma(yi, vi, data = x), slab = x[[1]], header = x[[2]][1])
}
Map(myfunction,the_list) -> the_plots




the_plots
#> $Mathew
#> $Mathew$xlim
#> [1] -36.28 65.39
#>
#> $Mathew$alim
#> [1] -10 30
#>
#> $Mathew$at
#> [1] -10 0 10 20 30
#>
#> $Mathew$ylim
#> [1] -1.5 11.0
#>
#> $Mathew$rows
#> [1] 1 2 3 4 5 6 7 8
#>
#> $Mathew$cex
#> [1] 1
#>
#> $Mathew$cex.lab
#> [1] 1
#>
#> $Mathew$cex.axis
#> [1] 1
#>
#>
#> $Simp
#> $Simp$xlim
#> [1] -47.99 63.44
#>
#> $Simp$alim
#> [1] -10 30
#>
#> $Simp$at
#> [1] -10 0 10 20 30
#>
#> $Simp$ylim
#> [1] -1.5 11.0
#>
#> $Simp$rows
#> [1] 1 2 3 4 5 6 7 8
#>
#> $Simp$cex
#> [1] 1
#>
#> $Simp$cex.lab
#> [1] 1
#>
#> $Simp$cex.axis
#> [1] 1
#>
#>
#> $hutten
#> $hutten$xlim
#> [1] -54.10 30.23
#>
#> $hutten$alim
#> [1] -25 5
#>
#> $hutten$at
#> [1] -25 -20 -15 -10 -5 0 5
#>
#> $hutten$ylim
#> [1] -1.5 11.0
#>
#> $hutten$rows
#> [1] 1 2 3 4 5 6 7 8
#>
#> $hutten$cex
#> [1] 1
#>
#> $hutten$cex.lab
#> [1] 1
#>
#> $hutten$cex.axis
#> [1] 1
#>
#>
#> $AB
#> $AB$xlim
#> [1] -21.86 17.49
#>
#> $AB$alim
#> [1] -10 5
#>
#> $AB$at
#> [1] -10 -5 0 5
#>
#> $AB$ylim
#> [1] -1.5 11.0
#>
#> $AB$rows
#> [1] 1 2 3 4 5 6 7 8
#>
#> $AB$cex
#> [1] 1
#>
#> $AB$cex.lab
#> [1] 1
#>
#> $AB$cex.axis
#> [1] 1
# there is no way to turn these into back into plot objects;
# see lines 654-699 of source code for forest.rma
# https://github.com/wviechtb/metafor/blob/master/R/forest.rma.r
# that function would have to be modified to save the
# plot objects
# These will work individually, however, if myfunction is
# applied not to the_list, but to Simp, Mathew, etc. extracted
# as follows
make_frame <- function(x) assign(x,data.frame(Map(run_escalc,x)[[x]]), envir = .GlobalEnv)
invisible(Map(make_frame,the_vars))
myfunction(Simp)

# you can use lapply to display them all
lapply(the_list,myfunction)




#> $Mathew
#> $Mathew$xlim
#> [1] -36.28 65.39
#>
#> $Mathew$alim
#> [1] -10 30
#>
#> $Mathew$at
#> [1] -10 0 10 20 30
#>
#> $Mathew$ylim
#> [1] -1.5 11.0
#>
#> $Mathew$rows
#> [1] 1 2 3 4 5 6 7 8
#>
#> $Mathew$cex
#> [1] 1
#>
#> $Mathew$cex.lab
#> [1] 1
#>
#> $Mathew$cex.axis
#> [1] 1
#>
#>
#> $Simp
#> $Simp$xlim
#> [1] -47.99 63.44
#>
#> $Simp$alim
#> [1] -10 30
#>
#> $Simp$at
#> [1] -10 0 10 20 30
#>
#> $Simp$ylim
#> [1] -1.5 11.0
#>
#> $Simp$rows
#> [1] 1 2 3 4 5 6 7 8
#>
#> $Simp$cex
#> [1] 1
#>
#> $Simp$cex.lab
#> [1] 1
#>
#> $Simp$cex.axis
#> [1] 1
#>
#>
#> $hutten
#> $hutten$xlim
#> [1] -54.10 30.23
#>
#> $hutten$alim
#> [1] -25 5
#>
#> $hutten$at
#> [1] -25 -20 -15 -10 -5 0 5
#>
#> $hutten$ylim
#> [1] -1.5 11.0
#>
#> $hutten$rows
#> [1] 1 2 3 4 5 6 7 8
#>
#> $hutten$cex
#> [1] 1
#>
#> $hutten$cex.lab
#> [1] 1
#>
#> $hutten$cex.axis
#> [1] 1
#>
#>
#> $AB
#> $AB$xlim
#> [1] -21.86 17.49
#>
#> $AB$alim
#> [1] -10 5
#>
#> $AB$at
#> [1] -10 -5 0 5
#>
#> $AB$ylim
#> [1] -1.5 11.0
#>
#> $AB$rows
#> [1] 1 2 3 4 5 6 7 8
#>
#> $AB$cex
#> [1] 1
#>
#> $AB$cex.lab
#> [1] 1
#>
#> $AB$cex.axis
#> [1] 1