--- title: "Location-Scale Meta-Analysis" author: "František Bartoš" date: "28th of April 2026" output: rmarkdown::html_vignette: self_contained: yes bibliography: ../inst/REFERENCES.bib csl: ../inst/apa.csl link-citations: true vignette: > %\VignetteIndexEntry{Location-Scale Meta-Analysis} %\VignetteEngine{knitr::rmarkdown} %\VignetteEncoding{UTF-8} %\VignetteEngine{knitr::rmarkdown_notangle} --- ```{r child = "_vignette-nowrap.md", echo = FALSE, eval = TRUE} ``` ```{r setup, include = FALSE} source("vignette-cache.R", local = knitr::knit_global()) cached_fits <- vignette_cache( name = "v12-metafor-parity-location-scale", objects = c("fit1_brma", "fit2_brma", "fit3_brma"), packages = c("metafor", "metadat") ) knitr::opts_chunk$set( collapse = TRUE, comment = "#>", eval = vignette_cache_eval(cached_fits), message = FALSE, warning = FALSE, fig.width = 7, fig.height = 3.5, dev = "png", fig.retina = 3 ) if (.Platform$OS.type == "windows") { knitr::opts_chunk$set(dev.args = list(type = "cairo")) } ``` ```{r load-models, include = FALSE} library("metafor") library("RoBMA") data("dat.bangertdrowns2004", package = "metadat") dat <- dat.bangertdrowns2004 dat$ni100 <- dat$ni / 100 dat$meta <- as.factor(dat$meta) dat$imag <- as.factor(dat$imag) vignette_cache_load(cached_fits) ``` ```{r precompute-models, include = FALSE, eval = FALSE} library("metafor") library("RoBMA") data("dat.bangertdrowns2004", package = "metadat") dat <- dat.bangertdrowns2004 dat$ni100 <- dat$ni / 100 dat$meta <- as.factor(dat$meta) dat$imag <- as.factor(dat$imag) fit1_brma <- brma( yi = yi, vi = vi, measure = "SMD", data = dat, seed = 1 ) fit2_brma <- brma( yi = yi, vi = vi, measure = "SMD", mods = ~ ni100, scale = ~ ni100, data = dat, seed = 1 ) fit3_brma <- brma( yi = yi, vi = vi, measure = "SMD", mods = ~ ni100 + meta, scale = ~ ni100 + imag, data = dat, seed = 1 ) vignette_cache_save(cached_fits) ``` This vignette illustrates how to fit location-scale meta-analytic models with the `RoBMA` R package [@RoBMA]. We walk through the writing-to-learn dataset (`dat.bangertdrowns2004` in the `metadat` package), showing the `brma()` call with the `scale = ~ ...` argument alongside its `metafor::rma()` counterpart [@metafor], and build from a simple random-effects model up to a model with different predictors in the location and scale parts. A location-scale model parameterizes the between-study heterogeneity $\tau$ as a regression on study-level covariates, allowing $\tau$ to vary across subgroups or as a function of continuous moderators. The location side is the usual meta-regression on the effect. All `brma()` calls keep the default prior distributions; see the [*Prior Distributions*](v01-prior-distributions.html) vignette for the prior definitions and customization options. For the non-location-scale workflow showing more features (also applicable to location-scale models) see the [*Bayesian Meta-Analysis*](v02-bayesian-meta-analysis.html) vignette. ## Random-Effects Model The writing-to-learn dataset contains 48 standardized mean differences from studies on the effectiveness of school-based writing interventions on academic achievement. Effect sizes `yi` and sampling variances `vi` are already provided. ```{r data} data("dat.bangertdrowns2004", package = "metadat") dat <- dat.bangertdrowns2004 dat$ni100 <- dat$ni / 100 dat$meta <- as.factor(dat$meta) dat$imag <- as.factor(dat$imag) ``` We start from a standard random-effects model with no scale predictor. ```{r fit1-metafor} fit1_metafor <- metafor::rma(yi, vi, data = dat) fit1_metafor ``` ```{r fit1-brma, eval = FALSE} fit1_brma <- brma( yi = yi, vi = vi, measure = "SMD", data = dat, seed = 1 ) ``` ```{r fit1-brma-summary} summary(fit1_brma, include_mcmc_diagnostics = FALSE) ``` The pooled effect and the single heterogeneity parameter `tau` track the REML estimates from the `metafor` package. ## Adding a Scale Predictor Sample size is a natural candidate for explaining heterogeneity: smaller studies often produce more variable estimates than larger ones. We add the rescaled total sample size (`ni100 = ni / 100`) on both the location and the scale side. ```{r fit2-metafor} fit2_metafor <- suppressWarnings(metafor::rma( yi, vi, mods = ~ ni100, scale = ~ ni100, data = dat )) fit2_metafor ``` ```{r fit2-brma, eval = FALSE} fit2_brma <- brma( yi = yi, vi = vi, measure = "SMD", mods = ~ ni100, scale = ~ ni100, data = dat, seed = 1 ) ``` ```{r fit2-brma-summary} summary(fit2_brma, include_mcmc_diagnostics = FALSE) ``` The location coefficient on `ni100` matches the REML estimate from the `metafor` package. The scale block reports `exp(intercept)`, which is the baseline between-study heterogeneity $\tau$, and a multiplicative log-scale coefficient on `ni100`. Note that `metafor::rma()` regresses $\log(\tau^2)$ while `brma()` regresses $\log(\tau)$, so each `metafor` package scale coefficient is exactly twice the matching `brma()` coefficient. `regplot()` with `pi = TRUE` overlays the prediction interval on the bubble plot. Because $\tau$ varies with `ni100`, the `brma()` prediction band narrows toward larger sample sizes, where the model expects less between-study heterogeneity. `metafor::regplot()` does not currently draw a prediction interval for location-scale models, so only the regression line and confidence band appear in the left panel. ```{r regplot-pi, fig.width = 7, fig.height = 3.2} par(mfrow = c(1, 2), mar = c(4, 4, 2, 1), cex.axis = 0.8, cex.lab = 0.8, cex.main = 0.75) metafor::regplot(fit2_metafor, mod = "ni100", pi = TRUE, main = "metafor", xlim = c(0, 6), ylim = c(-1, 2)) regplot(fit2_brma, mod = "ni100", pi = TRUE, main = "RoBMA", xlim = c(0, 6), ylim = c(-1, 2)) ``` ## Different Variables in Location and Scale Variables in the location and scale parts can differ. The next fit adds metacognitive reflection (`meta`) on the location side and imaginative-component coding (`imag`) on the scale side. ```{r fit3-metafor} fit3_metafor <- suppressWarnings(metafor::rma( yi, vi, mods = ~ ni100 + meta, scale = ~ ni100 + imag, data = dat )) fit3_metafor ``` ```{r fit3-brma, eval = FALSE} fit3_brma <- brma( yi = yi, vi = vi, measure = "SMD", mods = ~ ni100 + meta, scale = ~ ni100 + imag, data = dat, seed = 1 ) ``` ```{r fit3-brma-summary} summary(fit3_brma, include_mcmc_diagnostics = FALSE) ``` Each scale coefficient is again roughly half of its `metafor` package counterpart, with `brma()` shrinking the larger `metafor` package estimates more strongly because of the weakly informative default prior distribution. ## Scale Predictor Prior Distributions Scale regression coefficients are not on the effect-size scale: they describe multiplicative changes in $\tau$ on the log scale and are therefore scale-free. The default prior distribution is `Normal(0, 0.5)`, which corresponds to a weakly informative prior distribution on the multiplicative effect. A coefficient of $\pm 0.5$ corresponds to roughly $\exp(\pm 0.5) \approx 0.6$ or $1.65$, i.e., a $40\%$ decrease or $65\%$ increase in $\tau$ per unit change in the predictor; values past $\pm 1$ correspond to a halving or doubling, which the prior distribution treats as extreme. `plot()` with `prior = TRUE` overlays the prior distribution on the posterior so the shift is visible. `parameter_scale = "ni100"` selects the scale-side coefficient (the analogous selector on the location side is `parameter_mods = "ni100"`). ```{r scale-posterior, fig.width = 5, fig.height = 3.5} par(mar = c(4, 4, 1, 1)) plot(fit2_brma, parameter_scale = "ni100", prior = TRUE, xlim = c(-2, 2)) ``` ## Other Inference Helpers All inference helpers available for any other `brma()` fit work the same way for the location-scale fit specified via `scale = ~ ...`. This includes meta-regression bubble plots, marginal means, model comparison via Bayes factors, multilevel structures via `cluster`, and the full set of model-fit diagnostics. See the [*Bayesian Meta-Analysis*](v02-bayesian-meta-analysis.html) vignette for more details. ## References