How to use this document. Every code block below can be pasted into an R session and run in order. To keep the vignette light (and fast to check on CRAN) the chunks are not executed while knitting (
eval = FALSEin the setup chunk); instead the expected results are reported and discussed in the text and tables. To run the analysis live, changeeval = FALSEtoeval = TRUE, or simply copy each block into R.
This vignette is based on “Challenge I” of Mulder, J., &
Pfadt, J. M. (2026), “Going in the right direction: A tutorial to
directional hypothesis testing using the BFpack module in JASP,”
Advances in Methods and Practices in Psychological Science. It adapts
that example into a hands-on tutorial for the BFpack R
package (the JASP module wraps the same engine). The data and R syntax
are available on the online supplement, https://osf.io/7q5pf/.
Substantive theories rarely predict merely that “some effect is
non-zero.” They predict a pattern: which predictors matter, and
how their effects are ordered relative to one another. The conventional
workflow answers such questions indirectly — one significance test per
coefficient and per contrast, followed by a multiplicity correction —
and then asks the reader to reassemble the pieces into support for a
theory. This tutorial shows the direct alternative: translating
each competing theory into a single constrained hypothesis and comparing
them head-to-head with Bayes factors using
BFpack.
A large literature in social and political psychology argues that perceived competition and status threat shape attitudes toward immigrant groups. Within Ethnic Competition Theory (Scheepers et al., 2002), socioeconomic position — social class, education, and income — is positively associated with attitudes toward immigrants. Different theoretical accounts, however, predict different orderings of these three effects. We compare three theory-consistent variants plus a complement, using European Values Study (EVS) data from Germany.
We fit a linear regression of attitude on
class, education, and income,
with gender as a covariate, and test:
\[ \begin{aligned} \mathcal{H}_1 &:\ \beta_{\text{class}} > \beta_{\text{education}} > \beta_{\text{income}} > 0 \quad(\text{strict ordering: class matters most}),\\ \mathcal{H}_2 &:\ \beta_{\text{education}} > (\beta_{\text{class}}, \beta_{\text{income}}) > 0 \quad(\text{education largest; class and income unordered}),\\ \mathcal{H}_3 &:\ \beta_{\text{class}} = \beta_{\text{education}} = \beta_{\text{income}} > 0 \quad(\text{three equal positive effects}),\\ \mathcal{H}_4 &:\ \text{complement (none of } \mathcal{H}_1,\mathcal{H}_2,\mathcal{H}_3). \end{aligned} \]
A traditional point-null in which all three effects are zero is not of substantive interest here and is therefore excluded.
By the end you will be able to:
get_estimates();BFpack hypothesis string using >,
=, (,), and ;;BF() and interpret them; andThe data set regression_EVS_Germany.csv lives on the OSF
project https://osf.io/7q5pf/. The osfr package can
pull it directly from that project ID in a single call — no manual
download, no hard-coded file link:
install.packages("osfr") # once, if needed
library("osfr") # must be loaded before the download call below
EVS_Germany <- read.csv(osf_download(
osf_ls_files(osf_retrieve_node("7q5pf"),
pattern = "regression_EVS_Germany.csv"),
conflicts = "overwrite")$local_path)Prefer to avoid a dependency? Download
regression_EVS_Germany.csvmanually from https://osf.io/7q5pf/ and read it locally withread.csv("regression_EVS_Germany.csv"). If you know the file’s own OSF short link,read.csv("https://osf.io/<id>/download")is an equivalent one-liner.
We drop incomplete cases, standardize the continuous
variables so that the coefficients are on a common scale (essential for
the ordering tests — an ordering of coefficients is then an ordering of
effect sizes), and treat gender as a factor:
EVS_Germany <- EVS_Germany[complete.cases(EVS_Germany), ]
EVS_Germany$attitude <- c(scale(EVS_Germany$attitude))
EVS_Germany$education <- c(scale(EVS_Germany$education))
EVS_Germany$income <- c(scale(EVS_Germany$income))
EVS_Germany$class <- c(scale(EVS_Germany$class))
EVS_Germany$gender <- as.factor(EVS_Germany$gender)
fit1 <- lm(attitude ~ education + income + gender + class, data = EVS_Germany)This step is essential. The names used in a
hypothesis string are exactly the names printed by
get_estimates(). For this model the three predictors of
interest appear simply as class, education,
and income (the factor gender appears as
gender1).
Each theory becomes one constrained hypothesis; hypotheses are
separated by semicolons. (class, income) in \(\mathcal{H}_2\) means both are larger than
0 but unordered relative to each other. We do
not write \(\mathcal{H}_4\) explicitly:
complement = TRUE adds it automatically.
set.seed(1) # the equality/order tests use sampling; fix the seed.
BF_App1 <- BF(fit1,
hypothesis = "class > education > income > 0;
education > (class, income) > 0;
class = education = income > 0",
complement = TRUE)
print(BF_App1)By default the four hypotheses receive equal prior probabilities of
\(1/4\). To weight them differently,
pass prior.hyp.conf (one weight per hypothesis, in order,
including the complement) — this shifts the posterior probabilities but
leaves the Bayes factors, which measure the evidence in the
data, unchanged.
Table 1. Bayes factors (rows vs. columns) and posterior probabilities for the four hypotheses (equal prior probabilities, seed = 1).
| Hypothesis | vs \(\mathcal{H}_1\) | vs \(\mathcal{H}_2\) | vs \(\mathcal{H}_3\) | vs \(\mathcal{H}_4\) | \(P(\mathcal{H}_t\mid\text{Data})\) |
|---|---|---|---|---|---|
| \(\mathcal{H}_1\): class > education > income > 0 | 1.000 | 0.312 | 1.059 | 49.53 | 0.193 |
| \(\mathcal{H}_2\): education > (class, income) > 0 | 3.207 | 1.000 | 3.397 | 158.8 | 0.620 |
| \(\mathcal{H}_3\): class = education = income > 0 | 0.944 | 0.294 | 1.000 | 46.75 | 0.183 |
| \(\mathcal{H}_4\): complement | 0.020 | 0.006 | 0.021 | 1.000 | 0.004 |
All three theory-consistent hypotheses receive strong evidence over the complement, and \(\mathcal{H}_2\) is favored among \(\mathcal{H}_1\)–\(\mathcal{H}_3\) (though the evidence is mild, with Bayes factors around 3). A possible write-up:
“The model comparison favored \(\mathcal{H}_2\) with a posterior probability of .620, suggesting the largest positive effect for education, followed by class and income. Evidence distinguishing \(\mathcal{H}_2\) from \(\mathcal{H}_1\) and \(\mathcal{H}_3\) was mild (Bayes factors of approximately 3). All three theory-consistent hypotheses strongly outperformed the complement (Bayes factors of 49.5, 159, and 46.7), indicating that parameter configurations inconsistent with Ethnic Competition Theory are unlikely given the data.”
Note how a single comparison expresses the relative support for each whole theory — exactly what the substantive question asks, and exactly what a collection of separate tests cannot deliver.
The classical route tests the three coefficients and their three
pairwise contrasts one at a time, then corrects for multiplicity.
Because jointly testing coefficients and contrasts is awkward in common
GUI software, we implement it in R with multcomp and
car. The contrast matrix K_eq has one row per
null hypothesis; its columns follow the model’s design-matrix order:
(Intercept), education, income, gender1, class.
library("multcomp")
library("car")
K_eq <- rbind(
"class - education = 0" = c(0, -1, 0, 0, 1),
"class - income = 0" = c(0, 0, -1, 0, 1),
"education - income = 0" = c(0, 1, -1, 0, 0),
"class = 0" = c(0, 0, 0, 0, 1),
"education = 0" = c(0, 1, 0, 0, 0),
"income = 0" = c(0, 0, 1, 0, 0)
)
glht_eq <- multcomp::glht(fit1, linfct = K_eq) # two-sided by default
ci_eq <- confint(glht_eq) # 95% (unadjusted) CIs
sm_eq <- summary(glht_eq, test = adjusted("holm")) # Holm-adjusted p-values
res <- data.frame(
Null_hypothesis = names(sm_eq$test$coefficients),
Estimate = as.numeric(sm_eq$test$coefficients),
LB_95 = ci_eq$confint[, "lwr"],
UB_95 = ci_eq$confint[, "upr"],
t_value = as.numeric(sm_eq$test$tstat),
p_adjusted = as.numeric(sm_eq$test$pvalues),
row.names = NULL
)
print(cbind(Null_hypothesis = res[, 1], round(res[, -1], 3)))Table 2. Null hypotheses on individual effects and pairwise contrasts. Confidence bounds are unadjusted; \(p\)-values are Holm-adjusted.
| Null hypothesis | Estimate | 95%-LB | 95%-UB | \(t\) value | Adjusted \(p\) |
|---|---|---|---|---|---|
| \(\beta_{\text{class}} - \beta_{\text{education}} = 0\) | −0.055 | −0.178 | 0.067 | −1.124 | 0.522 |
| \(\beta_{\text{class}} - \beta_{\text{income}} = 0\) | 0.081 | −0.022 | 0.184 | 1.962 | 0.150 |
| \(\beta_{\text{education}} - \beta_{\text{income}} = 0\) | 0.137 | 0.031 | 0.242 | 3.230 | 0.005 |
| \(\beta_{\text{class}} = 0\) | 0.106 | 0.034 | 0.178 | 3.672 | 0.001 |
| \(\beta_{\text{education}} = 0\) | 0.161 | 0.090 | 0.233 | 5.599 | 0.000 |
| \(\beta_{\text{income}} = 0\) | 0.025 | −0.039 | 0.089 | 0.960 | 0.522 |
This yields only fragmented evidence about the hypotheses of
interest. At \(\alpha = .05\) the
effects of class and education are significant
and positive while income is not; yet the contrast \(\beta_{\text{class}} -
\beta_{\text{income}}\) is nonsignificant while \(\beta_{\text{education}} -
\beta_{\text{income}}\) is significant. These partly conflicting
results arise because statistical significance is not transitive
across separate tests: class and
education both differing from zero while
income does not need not imply that every contrast lines up
accordingly.
Consequently the separate tests give no coherent ordering of the three effects and no way to quantify support for the competing theories. Eyeballing the estimates and CIs suggests the data most favor \(\mathcal{H}_2\) (largest effect for education) — but how much more than \(\mathcal{H}_1\), \(\mathcal{H}_3\), or \(\mathcal{H}_4\) remains unclear. The joint Bayes-factor test in the previous section answers exactly that question with a single set of numbers.
Mulder, J. (2016). Bayes factors for testing order-constrained hypotheses on correlations. Journal of Mathematical Psychology.
Mulder, J., & Pfadt, J. M. (2026). Going in the right direction: A tutorial to directional hypothesis testing using the BFpack module in JASP. Advances in Methods and Practices in Psychological Science.
Mulder, J., et al. (2021). BFpack: Flexible Bayes factor testing of scientific expectations in R. Journal of Statistical Software, 100(18), 1–63.
Scheepers, P., Gijsberts, M., & Coenders, M. (2002). Ethnic exclusionism in European countries. European Sociological Review, 18(1), 17–34.