Volume 7 · ten days
A posterior says what you wanted p to say.
Bayesian estimation gives you the probability distribution of the effect given your data — which is the quantity every misreading of a p-value was reaching for. The formula syntax is the one you already know; only the engine changes. The second half of this volume is meta-analysis: effect-size coding, random-effects models, heterogeneity, moderators, bias diagnostics, and a forest plot that publishes.
10
days, posterior to forest plot
1
formula syntax, two engines
4
diagnostics before any posterior is read
Day 67
45 minutes · plain language, no integrals
Prior, likelihood, posterior — and why intervals stop being awkward
Bayesian inference combines what you believed before (the prior) with what the data say (the likelihood) to give a probability distribution over parameter values (the posterior). The practical consequence: you can say "there is a 95% probability the effect lies between 0.31 and 0.94" — the sentence everyone wants to say about a confidence interval and technically cannot.
| Question | Frequentist answer | Bayesian answer |
|---|---|---|
| What is the probability the effect is positive? | Undefined — the effect is fixed, not random. | A number: the proportion of the posterior above zero. |
| What does the interval mean? | 95% of such intervals, computed this way forever, would contain the true value. | Given the data and prior, 95% probability the value is in here. |
| Can I get evidence FOR no effect? | Not from a p-value. Requires equivalence testing. | Yes — a posterior concentrated near zero, or a Bayes factor. |
| What about prior knowledge? | No formal role. | Stated explicitly, and its influence is checkable. |
| Does it depend on my intentions? | Yes — optional stopping invalidates p-values. | No — the posterior does not depend on the sampling plan. |
Where the objection usually lands
"The prior is subjective." It is explicit, which is different. Every analysis contains assumptions — a link function, a distributional family, which covariates to include, when to stop collecting; Bayesian analysis writes one more of them down where a reviewer can see and challenge it. The defensible practice is weakly informative priors plus a sensitivity analysis showing the conclusion does not hinge on them (day 69), and at that point the objection has an answer rather than a shrug.
Do this now · 15 minutes
Run the grid example, then change the prior to
dbeta(theta, 20, 20) (a strong belief in 0.5) and re-run. Watch the posterior move, and note how much data it would take to overcome it. That intuition is what day 69 formalises.
Day 68
55 minutes · the formula you already know
brms: same syntax, different engine
brm() takes lme4 formulas. Everything you learned in volumes 4 and 5 transfers directly — fixed effects, random effects, families — and the output is a posterior instead of a point estimate with a standard error. The costs are compilation time (30–90 seconds the first run) and the discipline of checking convergence before reading anything.
| Output column | Read it as |
|---|---|
Estimate | Posterior mean (or median with robust = TRUE). The analogue of a coefficient. |
Est.Error | Posterior SD. Loosely the analogue of a standard error. |
l-95% CI / u-95% CI | The credible interval: 95% posterior probability the parameter lies here. |
Rhat | Convergence. Must be ≤ 1.01 for every parameter or the numbers are meaningless. |
Bulk_ESS / Tail_ESS | Effective sample size. Want > 1,000 for stable intervals; low tail ESS means unreliable interval endpoints. |
sigma | Residual SD, estimated with its own posterior rather than plugged in. |
Do this now · 25 minutes
Install brms, verify the compiler works, and refit your main volume-4 model with
brm(). Put the frequentist and Bayesian estimates side by side — with weak priors and decent n they will be nearly identical, which is the reassurance you need to trust the machinery.
Day 69
55 minutes · explicit, not arbitrary
Weakly informative priors, and the check that justifies them
brms picks flat or very wide default priors, which is safe but wasteful: you always know something, if only the plausible scale of your measure. A weakly informative prior rules out the absurd (a 400-point BDI-II effect) while leaving every plausible value well supported. The way to justify one is a prior predictive check — simulate data from the prior alone and see whether it looks like data from your field.
| Prior type | Example | When |
|---|---|---|
| Flat / improper | uniform(-Inf, Inf) | Almost never by choice. Can prevent convergence in hierarchical models. |
| Weakly informative | normal(0, 5) on a 0–63 scale | The default you should set. Rules out the absurd, nothing else. |
| Regularising | normal(0, 1) on standardised predictors | Many predictors; shrinks noise toward zero. Volume 10's lasso, in Bayesian form. |
| Informative | normal(0.3, 0.1) from a meta-analysis | Strong prior evidence exists. Cite it; run the sensitivity analysis. |
| On variance components | student_t(3, 0, 2.5), exponential() | Hierarchical models with few clusters, where flat priors misbehave. |
Do this now · 25 minutes
Set explicit weakly informative priors for your model, run a prior predictive check, and then the three-prior sensitivity table. Keep the table — it belongs in your supplementary materials.
Day 70
50 minutes · never read an unconverged posterior
R-hat, ESS, divergences: four checks, in order
MCMC gives you samples from the posterior only if the sampler explored it properly. Four diagnostics catch almost every failure, and they are cheap. Run them before you look at a single estimate — an unconverged model can print perfectly reasonable-looking numbers that are simply wrong.
| Symptom | Cause | Fix |
|---|---|---|
| R-hat > 1.01 | Chains have not mixed | More iterations; check for a misspecified or unidentified model |
| Divergent transitions | Curved posterior geometry, often a hierarchical funnel | control = list(adapt_delta = 0.99); tighter priors on variance components; non-centred parameterisation |
| Low tail ESS | Heavy tails, poor exploration | More iterations; reparameterise; check for outliers |
| Max treedepth warnings | Sampler hitting its step limit | max_treedepth = 15. Usually efficiency, not validity |
| Chains stuck at different values | Multimodality or non-identification | Model problem, not a sampler problem. Rethink the specification |
| Bimodal posterior | Genuinely two solutions, or label switching in mixtures | Investigate; do not summarise with a mean |
The habit
Wrap the checks in a function and call it on every fit:
check_brms <- function(m) list(rhat = max(rhat(m), na.rm = TRUE), div = sum(...), ess = min(neff_ratio(m), na.rm = TRUE)). Three numbers, one line, and you will never again present a supervisor with an unconverged model. Report all three in the paper: "all R-hat ≤ 1.01, no divergent transitions, minimum bulk ESS 9,812" is a complete convergence statement.Do this now · 20 minutes
Write your own
check_brms() helper (it is also day 94's first function) and run it on every model you have fitted. Then deliberately break one — fit with iter = 300 — and watch which diagnostics catch it.
Day 71
55 minutes · the sentences to write
Credible intervals, probability of direction, and ROPE
Bayesian reporting has its own conventions, and bayestestR produces all of them in one call. The core set: a point estimate (median), a 95% credible interval (HDI or quantile, say which), the probability of direction, and — when the question is whether an effect is negligible — the proportion of the posterior inside a region of practical equivalence you defined in advance.
| Quantity | Meaning and use |
|---|---|
| Median | The posterior's central value. More robust than the mean for skewed posteriors. |
| 95% HDI | The narrowest interval containing 95% of the posterior. Say HDI or quantile — they differ for skewed posteriors. |
| pd (probability of direction) | Posterior probability the effect has the sign you estimated. Ranges 50–100%. Correlates closely with 1 − p/2, but means something. |
| ROPE | Region of practical equivalence: the interval of effects too small to matter. Define it from your measure, before fitting. |
| % in ROPE | How much posterior mass is negligible. Under 2.5% = effect is practically meaningful; over 97.5% = practically equivalent to null. |
| Bayes factor | Relative evidence for one model over another. A different question — day 72. |
Two sentences you can write, and one you cannot
Write: "CBT reduced post-treatment BDI-II by 3.41 points, 95% HDI [−4.82, −2.00], pd = 100%, with no posterior mass inside the ±1-point region of practical equivalence." Write: "The waitlist effect was practically equivalent to zero: 82% of its posterior lay within ±1 point." Do not write: "the effect was significant (pd > 95%)" — importing a threshold ritual into a framework built to avoid it defeats the purpose. Report the distribution and let the reader see the uncertainty.
Do this now · 25 minutes
Define a ROPE for your outcome in its own units and justify it in one sentence (smallest change a clinician or participant would notice). Then run
describe_posterior() and write your main result as a full Bayesian sentence.
Day 72
50 minutes · evidence for the null
Bayes factors: a different question, and a real caveat
Estimation asks how big the effect is. A Bayes factor asks which of two models the data support better — and crucially it can support the null, which no p-value can. The caveat is unavoidable: a Bayes factor depends on the prior on the effect size under the alternative, much more strongly than a posterior does. Report a sensitivity analysis or do not report the Bayes factor.
| BF10 | Conventional wording | Note |
|---|---|---|
| > 100 | Extreme evidence for H1 | Thresholds are conventions, not facts |
| 30–100 | Very strong | |
| 10–30 | Strong | |
| 3–10 | Moderate | The usual publishable band |
| 1–3 | Anecdotal | Effectively inconclusive; say so |
| 1/3–1 | Anecdotal for H0 | |
| < 1/3 | Moderate to strong for H0 | The claim a p-value cannot make: evidence of absence |
When to use which
Use estimation (days 68–71) when your question is how large the effect is — which is nearly always the scientific question, and what a preregistered SESOI implies. Use a Bayes factor when the question is genuinely model comparison: does this effect exist at all, is a null-hypothesis replication informative, should I keep this parameter. Use equivalence testing (frequentist TOST,
TOSTER) if you want a null claim without priors. Reporting all three is not showing off; it is showing that the conclusion does not depend on the framework.Do this now · 20 minutes
Compute a Bayes factor for your main comparison across five prior widths and tabulate them. Then write both sentences — the estimation one and the Bayes-factor one — and decide which answers your research question. Say why in the paper.
Day 73
55 minutes · where brms beats lme4
Hierarchical models: convergence without compromise
This is the day that converts people. The mixed models that made lme4 print singular-fit warnings — a random slope with 15 clusters, a complex crossed structure, a GLMM with few observations per cell — often fit cleanly in brms, because a mild prior on the variance components regularises the estimate instead of pushing it to the boundary. You also get uncertainty on the variance components themselves, which frequentist mixed models do not give you.
The reframing of "singular fit"
In lme4 a variance estimated at zero is a boundary problem you must respond to by simplifying. In brms the posterior for that SD simply concentrates near zero with a credible interval — which is the honest answer ("these data cannot distinguish much slope variation") rather than a warning. You keep the model your design justifies and report the uncertainty. That is a better scientific outcome than deleting the random slope and pretending you never wanted it.
Do this now · 25 minutes
Refit your most troublesome lme4 model in brms with priors on the SD and correlation parameters. Compare: does the fixed effect change, and what does the interval on the slope SD tell you that lme4's point estimate did not?
Day 74
55 minutes · synthesis starts with coding
From published papers to a computable effect-size dataset
A meta-analysis is 80% data extraction and 20% modelling, and the modelling is easy. The skill is converting whatever each paper reported — means and SDs, a t-value, an F, an odds ratio, a correlation — into a common effect size with its sampling variance. escalc() does the conversions; your job is the coding sheet and the decisions it records.
| Coding decision | Record it as a column |
|---|---|
| Which outcome, when a study reports several | outcome_used plus a rule decided in advance (primary outcome; first-listed; all, with multilevel modelling) |
| Which timepoint | timepoint — post-treatment and follow-up are different questions |
| Multiple comparison groups in one study | comparison; do not treat them as independent effects |
| Direction of the scale | reverse_scored — sign errors are the most common meta-analytic bug |
| Risk of bias | Several columns from a published tool (e.g. randomisation, blinding, attrition) |
| Moderators you preregistered | One column each: dose, population, format, country |
| Who coded it | coder — and double-code 20%, reporting kappa (volume 6 day 56) |
Dependence: the structural issue to decide now
If one study contributes three effects (three outcomes, or two treatment arms against a shared control), those effects are correlated and treating them as independent understates uncertainty. Three defensible options: average them within study, pick one by a preregistered rule, or fit a multilevel meta-analysis with
rma.mv() and random = ~ 1 | study/effect — the modern default, and combine it with cluster-robust inference (clubSandwich). Decide before extraction, because it changes the coding sheet.Do this now · 25 minutes
Build a coding sheet with columns for effect-size inputs, your preregistered moderators, risk of bias, and the coder. Extract five studies, run
escalc(), and check every sign by re-reading what the original paper's direction actually was.
Day 75
55 minutes · the pooled estimate is not the point
Random-effects models, heterogeneity, and meta-regression
A fixed-effect meta-analysis assumes every study estimates the same true effect. It almost never holds in behavioural science, so random-effects is the default: it estimates a distribution of true effects, whose spread (tau) is often more interesting than the mean. A prediction interval — where the next study's true effect is likely to fall — is the most informative number you can report and the one most often omitted.
| Statistic | What it says |
|---|---|
| tau | SD of true effects across studies, in effect-size units. The directly interpretable measure. |
| I² | Share of observed variability due to real differences rather than sampling error. A relative measure — it can be high with trivial tau when studies are large. Never report it alone. |
| Q and its p-value | Test of homogeneity. Underpowered with few studies, over-powered with many. A screening device. |
| Prediction interval | Where the next study's true effect will plausibly fall. The most practically useful line of output. |
| Number of studies (k) | Everything is unstable below k ≈ 10. With k < 5, tau is barely estimable — say so. |
Two limits on meta-regression
It is observational at the study level: dose was not randomly assigned to studies, so a dose moderator is confounded with everything that co-varies with dose across literatures. And it is badly underpowered — the informal rule is at least 10 studies per moderator, and testing four moderators on 24 studies is exploratory whatever the p-values say. Preregister the moderators, report them all, and describe the analysis as what it is.
Do this now · 25 minutes
Fit the random-effects model on your extracted data, report tau, I² and the prediction interval, then test your one preregistered moderator. Write the sentence a clinician would need: the pooled effect, and how much it varies between settings.
Day 76
55 minutes · the volume's deliverable
Small-study effects, p-curve, and a forest plot that publishes
Published effects are a biased sample of conducted effects: significant results are more likely to be written up, submitted and accepted. No diagnostic proves publication bias, and all of them confound it with genuine heterogeneity, so the standard is to run several and describe the pattern honestly. Then build the two figures every meta-analysis needs.
Portfolio artefact: the synthesis
A completed meta-analysis is the single most reusable object in an early research career — it justifies your next study's design, supplies the prior for day 69, gives the expected effect size for day 78's power analysis, and is publishable on its own. Deposit the coding sheet and the analysis script on OSF (volume 9 day 96) and it also becomes a living resource other people cite.
Do this now · 25 minutes
Run the funnel plot, Egger's test, trim-and-fill and PET-PEESE on your dataset, and write the honest paragraph: what the asymmetry could mean, and what range of pooled effects the diagnostics jointly support. Then export the forest plot at 170 mm and 600 dpi.
Reference
Seven steps, in this order
Fitting the model is step four. Most of the value is in the steps around it.
1 · write the model
Same formula as lme4. Family chosen from the outcome, not from habit.
2 · set priors
Weakly informative, on the scale of your measure.
get_prior() shows what exists.3 · prior predictive check
sample_prior = "only" then pp_check(). Does the prior imply plausible data?4 · fit
4 chains, 4,000 iterations,
cores = 4, a fixed seed.5 · diagnose
R-hat ≤ 1.01, zero divergences, ESS > 1,000, trace plots. Before reading anything.
6 · posterior predictive check
pp_check() and grouped statistics. Does the model generate your data?7 · report
Median, 95% HDI, pd, ROPE, prior sensitivity table, convergence statement.
Reference
From search to forest plot
| Stage | What it produces |
|---|---|
| Preregistered protocol (PROSPERO) | Question, inclusion criteria, search strategy, planned moderators, dependence rule |
| Search and screening | PRISMA flow diagram — PRISMA2020 generates it from your counts |
| Double coding | Coding sheet plus an agreement statistic (volume 6 day 56) |
escalc() | yi and vi per effect — the analysis dataset |
rma() / rma.mv() | Pooled estimate, tau, I², prediction interval |
| Meta-regression | Moderator tests, R², residual heterogeneity |
| Bias diagnostics | Funnel, Egger, trim-and-fill, PET-PEESE, p-curve |
| Figures | Forest plot and funnel plot at journal specification |
| Deposit | Coding sheet, script and PRISMA diagram on OSF |
Reference
Install these once
Checkpoint
Six questions before volume 8
{{ quizCounter }}
{{ quizScore }}
{{ quizQ }}
{{ quizFb }}
Next: volume 8
Everything so far happened after data collection. Volume 8 moves before it: analytic power, simulation-based power for designs no formula covers, power for mixed models, sensitivity analysis, sequential designs, and multiverse analysis as a robustness result.
Start day 77 →