RBayesflow

A Guided Path from Bayesian Textbook to Working Analysis

September 26, 2026 · John Peach

Introduction

Suppose you’ve taken an introductory statistics course in R, fit a few linear regressions in RStudio, and read a chapter or two about Bayesian statistics. You install Stan - Software for Bayesian Data Analysis, copy a model from a tutorial, and it runs. A table of numbers appears: estimates, standard errors, intervals. It looks exactly like the output from a regression you trust.

The problem is that this table appears whether or not the computation worked. Bayesian models are fit by a random-sampling algorithm, and when the algorithm fails it usually fails without generating error messages. The numbers still print, and nothing turns red. A beginner has no way to know the estimates are meaningless.

I recently found a textbook that guides users through the process of Bayesian analysis and thought it would be very useful to build it into a software system I called RBayesflow using the R statistics language.

The R community has built good tools for every step of a Bayesian analysis. brms fits models, bayesplot draws diagnostic plots, loo compares models, and posterior computes convergence checks. None of them tells you which to use first, what to check before you believe anything, or when to stop. Experienced analysts know the right order, but students may not.

RBayesflow is a set of R scripts, Quarto templates, and a small bookkeeping object that puts those tools in order. It follows the Bayesian workflow described by Andrew Gelman, Aki Vehtari, and colleagues, first as a 2020 paper and now as a 2026 book. It adds no new statistics — every calculation is done by an established library. What it adds is sequence, a record of every decision, and one rule it enforces: you cannot read results from a model whose diagnostics failed until you’ve looked at the failure and said so.

This article covers:

1. What Bayesian analysis is, and when you shouldn’t use it;
2. What the Bayesian workflow is and why it helps;
3. The design ideas behind RBayesflow;
4. How the code is organized, and the help available at each skill level;
5. A worked example from a real clinical trial;
6. Short descriptions of fifteen case studies from the book, each reproduced with RBayesflow;
7. The limitations you should know about before relying on it.

A glossary at the end defines every technical term used.

Bayesian Analysis in Brief

What it is

Bayesian analysis treats unknown quantities as things you can have beliefs about, and uses probability to describe those beliefs. You start with what you believe before seeing the data, then update it with the data. The rule for the update is Bayes’ theorem:

where

• $\theta$ (theta) is the set of unknown quantities you want to learn about, such as the effect of a drug or the average reaction time of a sleepy person.
• $y$ is the data you observed.
• $p(\theta)$ is the prior: how plausible each value of $\theta$ seemed before you saw the data.
• $p(y \mid \theta)$ is the likelihood: how probable your observed data would be if $\theta$ took a particular value.
• $p(\theta \mid y)$ is the posterior: how plausible each value of $\theta$ is after seeing the data.
• $\propto$ means “is proportional to.” The two sides differ only by a constant that makes the probabilities add up to one.

The intuition is a weighted vote. Each possible value of $\theta$ starts with a weight from the prior. The data then multiply each weight by how well that value explains what you saw. Values that explain the data well gain weight; values that don’t, lose it.

For almost every realistic model, the posterior can’t be written down as a formula. Instead, a program called Stan draws thousands of random samples from it using Markov chain Monte Carlo (MCMC), an algorithm that wanders through the space of possible parameter values, spending more time where the posterior is high. Summaries such as means and intervals are then computed from the samples. If you’re not careful, the summaries may look good but may be hiding problems with the analysis. Stan normally runs four independent chains, each starting from a different random point. If the chains end up describing the same distribution, that’s good evidence the sampler worked. If they disagree, at least one of them is wrong.

Where it helps

Bayesian methods are worth the extra effort in several situations. Grouped data — students in schools, patients in hospitals, repeated measurements on the same person — are hierarchical, and Bayesian hierarchical models share information between groups in a principled way called partial pooling No pooling estimates each group alone, so small groups get noisy answers. Complete pooling treats all groups as one, which hides real differences. Partial pooling is the compromise, and the data decide how much to share. . When you have real prior knowledge (reaction times aren’t negative and rarely exceed a few seconds), a prior can encode it. When you need the full uncertainty rather than a best guess, the posterior gives you that. And when a model is built from a physical or biological mechanism where the equations don’t fit standard regression software but can be written in Stan, Bayes is often the only practical path.

Where it doesn’t

Simpler methods are often better. A t-test or linear model answers a straightforward comparison in seconds on a large sample. Firth penalized logistic regression handles rare binary outcomes in small samples, and was designed specifically for that case. And if you can’t check whether the sampler worked, a well-understood frequentist method is the safer choice.

RBayesflow builds this judgment into its first phase. Before any model is fit, a function called assess_offramps() looks at the outcome type, the sample size, and the stated goal, then lists alternatives. Table 1 shows its decision rules.

Table 1: Off-ramp rules applied by assess_offramps() before any fitting.

The alternatives are shown with equal weight, and none is labelled “recommended.” If you pick a non-Bayesian method, exit_workflow() writes a short machine-readable log (a YAML file) recording what you chose, the evidence it was based on, and which alternatives you considered. Leaving the Bayesian path is a legitimate, documented outcome.

The Bayesian Workflow

From paper to book

In 2020, Gelman, Vehtari, and nine co-authors posted a long paper called Bayesian Workflow (arXiv:2011.01808). Its argument was that textbooks teach how to fit a model, while real analysis is a loop: build a model, check it, find what’s wrong, and build a better one. The paper tried to write down the unwritten habits of experienced analysts.

In 2026 that paper grew into a book of the same name, published by Chapman & Hall, with Richard McElreath joining as a lead author. The book has three parts on principles, statistical workflow, and computational workflow, followed by sixteen case studies (Chapters 16–31). Each case study takes one real or simulated problem and works through a sequence of models, showing the mistakes as well as the fixes. The authors maintain a companion website with code for every case study, and a free electronic edition for non-commercial use. The book describes itself as “not a checklist or cookbook.” RBayesflow is closer to a checklist than the authors intend their book to be. That’s a deliberate trade-off for beginners, who benefit from structure before they have the judgment to depart from it.

The seven phases

RBayesflow condenses the workflow into seven phases, each with its own Quarto template:

1. Goal declaration and data inspection: state what you want to learn, plot the data, and decide whether you need a Bayesian model at all.
2. Prior specification and prior predictive simulation: choose priors, then simulate fake data from the priors alone to check that they allow sensible values and rule out absurd ones.
3. Model fitting: run Stan through brms or cmdstanr.
4. MCMC diagnostics: check whether the sampler worked.
5. Posterior predictive checks: simulate new data from the fitted model and compare it with the real data.
6. Model comparison: estimate how well each model would predict new data, using leave-one-out cross-validation.
7. Reporting: produce a report in which every number traces back to a saved object.

Phases 2 to 6 usually repeat. A posterior predictive check reveals a problem, you change the model, and you go around again.

Why a workflow helps

Before a Bayesian result means anything, two separate questions have to be answered, and beginners tend to merge them. Did the computation work? That’s what Phase 4 checks — it’s about the sampler, not the science. Is the model any good? That’s what Phases 5 and 6 check. A sampler can work perfectly on a model that describes the data badly.

A fixed order makes it hard to skip either question, and it produces a record. When you return to an analysis six months later, or hand it to someone else, the reasons for each choice are written down.

RBayesflow Design Principles

Every computation is delegated to an established library — brms, cmdstanr, posterior, bayesplot, tidybayes, loo, ggplot2. If you want to know how a number was computed, the answer is in a well-documented package. RBayesflow is the ordering, not the math.

If Phase 4 finds a problem, the workflow withholds the coefficient table. It shows which checks failed and tells you to call wf$diagnose(), which displays the diagnostic plots and asks you to confirm you’ve reviewed them. Only then does output appear, and in practice mode with a warning attached. The acknowledgment is written into the record with a timestamp.

You choose one of three modes at the start:

• learn shows plots comparing the prior with the posterior instead of coefficient tables, so you see how the data changed your beliefs;
• practice shows standard output once diagnostics pass;
• expert shows everything with no diagnostic gate.

Off-ramps are offered before any fitting, with equal weight between the Bayesian and non-Bayesian options. When the optional AI assistant proposes code, it appears as text to read, not something that executes automatically.

A workflow-state object, wf_state, stores the mode, the priors and the reasons for them, a fingerprint of the data, every diagnostic result, and an audit trail of every decision. The data fingerprint is a SHA-256 hash, a short code computed from the data. If a single value changes, the code changes, so you can tell whether a saved model was fit to the data you think it was.

How RBayesflow Is Organized

The project folder

RBayesflow is a project folder, not an R package. You download it once and run all your analyses inside it.

RBayesflow/
├── R/            workflow functions (state object, diagnostics, off-ramps, help)
├── templates/    one Quarto template per phase, plus the report template
├── data/         one subfolder per analysis
├── docs/         the written user guide
└── renv.lock     exact versions of every package

Each analysis lives in its own subfolder under data/. A tiny start-up file in that subfolder loads the shared workflow code automatically, so several analyses can sit side by side without copying code or overwriting each other’s files. The renv.lock file renv is an R package that records the exact version of every library a project uses and can reinstall those versions later. Without it, a package update could silently change your results. pins every package to a tested version, so the analysis you run next year uses the same software as the one you ran today.

Getting started

You need R 4.3 or newer, RStudio, and a C++ compiler (on Windows, RTools). Stan compiles each model into a small C++ program, which is why the compiler is required. Installing the compiler is the step most likely to go wrong, especially on Windows. RBayesflow checks it before anything else and stops with instructions if it’s broken, because nothing after it can succeed.

# One manual step: install cmdstanr from the Stan project's repository
install.packages("cmdstanr", repos = c("https://stan-dev.r-universe.dev", getOption("repos")))

# From the RBayesflow folder: check and install everything else
source("R/install.R")
rbf_install()          # prints a ✓ or ✗ for each of seven checks

# Create a new analysis
rbf_new("my_first_analysis")

Then open data/my_first_analysis/ as your working directory in RStudio and start the workflow:

wf <- init_workflow(mode = "learn", stage = "explore")
guide(wf)

Three kinds of user, three kinds of help

RBayesflow was designed for two kinds of user, with a third accommodated but not catered to.

Students are the primary target: someone who knows R and some frequentist statistics, has read about Bayesian methods, and needs help with concepts rather than syntax. They should use learn mode. Working scientists who are competent in their domain, occasional brms users, but not statisticians, mainly want a reproducible and defensible result. They should use practice mode. Expert Bayesian statisticians already have their own workflow; they’re not who RBayesflow was built for, but expert mode removes the diagnostic gate and the prompts so it stays out of their way.

Help comes at three levels. The user guide in docs/user-guide/ covers installation, starting an analysis, each of the seven phases, and a plotting reference explaining how to read every plot the workflow produces. At any point, guide(wf) reads the current state and prints three lines: which phase you’re in, what normal output looks like, and the exact next command. Above that is Posit Assistant, which can be connected to RStudio or Positron. After each phase, RBayesflow writes a file called wf_context.json describing the current state, and each template tells the assistant to read it. In learn mode the assistant responds with a question rather than an answer; in practice mode it explains diagnostics on request; in expert mode it stays quiet.

Table 2 shows how the pieces fit together.

Table 2: User modes and the help available in each.

The AI layer is optional. Everything in RBayesflow works without it, and it’s the least tested part of the system (see Limitations). I connected Posit Assistant to OpenRouter, a service that gives access to many language models including some free ones, and confirmed that it works, but I haven’t yet used it through a full case study. Posit Assistant is a Posit product rather than open-source software. The free alternative is simply to go without it: the documentation and guide() cover every step.

A Worked Example: Did the Drug Reduce Cannabis Use?

Chapter 18 of the book reanalyzes a clinical trial by Lintzeris and colleagues (2019). The trial tested nabiximols, a cannabis-derived mouth spray, as a treatment for cannabis-use disorder. It enrolled 128 people, half on the drug and half on a placebo. At the start and at weeks 4, 8, and 12, each person reported how many of the previous 28 days they had used cannabis. After dropping missing values, there are 385 observations.

The outcome is a count between 0 and 28.

Phase 1: Should this be Bayesian at all?

The analysis starts in learn mode:

wf <- init_workflow(mode = "learn", stage = "explore")

The off-ramp assessment recognizes a bounded count with repeated measurements on each person and prints its options:

=== Off-Ramp Assessment ===
Outcome type: count (bounded 0-28, repeated measures)
n = 385 observations, 128 participants
Goal: Estimate treatment effect of nabiximols vs placebo

Available analysis paths (no preference implied):

  [1] Poisson GLM + quasi-likelihood
  [2] Negative binomial GLM
  [3] Full Bayesian (Stan via brms)

We choose the Bayesian path because the goal is to compare several models in detail, which needs the full posterior. The choice goes into the audit trail.

Phases 2 and 3: Two reasonable models

The original paper used a Gaussian (normal) model with a separate starting level for each person. That model ignores the 0-to-28 limit: it’s happy to predict $-3$ days or $31$ days.

The obvious fix is a binomial model, which treats each of the 28 days as a separate yes-or-no event with some probability $p$ of use:

Here $\text{cu}_n$ is the number of use days in observation $n$; $p_n$ is the chance of use on any one day; $\text{logit}(p) = \log\bigl(p/(1-p)\bigr)$ stretches a probability between 0 and 1 onto the whole number line so that a regression can predict it; $b_0$ is the baseline level; $b_\text{group}$, $b_\text{week}$, and $b_\text{gw}$ are the effects of treatment group, time, and their interaction (the interaction is the treatment effect of interest); and $u_{j[n]}$ is a personal offset for participant $j$, allowing some people to be heavier users than others. A binomial count can never go below 0 or above 28, so this model respects the limit exactly.

In brms, the model is one call:

fit_binomial <- brm(
  formula = cu | trials(set) ~ group * week + (1 | id),
  data    = cu_df, family = binomial(link = "logit"),
  prior   = priors_binomial,
  backend = "cmdstanr", seed = 1234
)

Phase 4: The diagnostic failure

The Gaussian model passes its diagnostics but the binomial model doesn’t:

!!! DIAGNOSTIC FAILURE !!!
Failed criteria:
  Rhat_max = 1.0512 (threshold < 1.01)
  bulk_ESS_min = 69 (threshold > 400)
  tail_ESS_min = 173 (threshold > 400)
-> Call wf$diagnose() to review the failing checks.

R-hat compares the four chains with each other; values near 1 mean they agree, and 1.05 means they don’t. Effective sample size (ESS) counts how many truly independent draws the samples are worth; 4,000 correlated draws here are worth only about 69. Consecutive MCMC draws resemble each other because each step starts from the last. It’s like 4,000 photos of a slowly moving car: many frames, but few genuinely different views.

At this point a beginner using brms directly would still see a coefficient table. RBayesflow shows nothing but the failure. After wf$diagnose() and an acknowledgment, learn mode still shows only the failure message, now marked as reviewed. Saying you’ve looked at a problem doesn’t make it go away.

Phase 5: Why the “right” model is wrong

The sampler’s trouble is a symptom. The cause shows up when we simulate data from the model and compare it with reality. About 14% of the real observations are zero: people who stopped using entirely during a four-week window. The binomial model predicts only 2–4% zeros.

If each day is an independent coin flip with, say, a 50% chance of use, the chance of 28 non-use days in a row is $0.5^{28}$, about one in 270 million. Real people don’t behave like coin flips. Someone who abstains on Monday is much more likely to abstain on Tuesday. The binomial model has no way to express that, so its predictions are too narrow and miss both ends of the range.

For each observation, the LOO-PIT diagnostic asks where that observation falls within the model’s prediction for it, made without using that observation. A well-calibrated model produces a flat line. Figure 1 shows the binomial model’s U shape: far too many observations land at the extreme edges of their predictions.

The fix is a beta-binomial model, which lets each person’s daily probability vary from one four-week window to the next. It adds a single extra parameter that controls how spread out the counts are. Its diagnostics pass, and its LOO-PIT curve is close to flat (Figure 2).

Phase 6: The surprising ranking

Leave-one-out cross-validation estimates how well each model would predict an observation it hadn’t seen. The score is the expected log predictive density (elpd); higher is better, and differences of more than a few units are meaningful.

Table 3: Leave-one-out model comparison for the clinical trial.

The binomial model, chosen because it respects the 0-to-28 limit, finishes last, about 300 units behind a Gaussian model that predicts impossible values. The Gaussian spreads its predictions widely, so the real zeros and 28s don’t surprise it much. The binomial is confident and wrong. The beta-binomial beats both by a wide margin.

The lesson for a student is that picking a distribution that matches the type of data isn’t enough — it also has to match how spread out the data are. That insight is the book’s, not the software’s. What RBayesflow did was stop anyone from reading treatment effects out of the binomial fit before the problem was understood, and keep a record of each step along the way.

The Case Studies

I reproduced the fifteen case studies in Chapters 16 through 30 with RBayesflow. Each has a detailed write-up with code and figures in the RBayesflow repository. Chapter 31, on simulation-based calibration, was omitted. The book presents it as a check you run on your modelling code before a workflow begins, not a workflow phase in itself.

The case study reports are in Markdown format (.md) and may be read from RStudio with the rmarkdown package. Install the package and open the report, e.g. ch30_planetary_motion.md from Chapter 30 Differential equation for planetary motion in the RStudio editor. At the top of the editor you should see a button called “Preview”. Click that and it will open a new window with the document neatly formatted.

Chapter 16: Movie ratings. Two films both average 4 stars, one from 2 ratings and one from 100; which should you trust? Three models on simulated data show why the larger sample deserves more trust, and then why simple averages mislead once tough and lenient raters choose different films. The final model estimates a quality for each movie and a strictness for each rater at the same time.

Chapter 17: Sleep study. Reaction times of 18 people over eight days of sleep deprivation are used to teach prior specification. Three models of increasing complexity each introduce a new parameter, and each prior is justified from what’s known about reaction times rather than taken from a default. The chapter shows how to split one variance prior into between-person and within-person parts as the model grows.

Chapter 18: Clinical trial. The worked example above. The binomial model, chosen because it respects the outcome’s limits, fails its diagnostics and loses to a Gaussian model; a beta-binomial model beats both.

Chapter 19: Coronavirus testing. An early-2020 antibody survey found 50 positives among 3,330 people, but the test’s false-positive rate alone could explain a third of them. Three models, ending in a hierarchical model of test accuracy across 13 calibration studies, give nearly the same central estimate of prevalence, about 1%. The width of the uncertainty changes enormously, and it depends almost entirely on one prior.

Chapter 20: Classification competition. A $100,000 challenge asked entrants to sort 1,000 synthetic time series into trending and trendless groups, scoring at least 900 correct. A three-part mixture model of the estimated trends turns the question into a decision: the chance of reaching 900 is about one in 240,000, so the entry fee isn’t worth it.

Chapter 21: Stochastic learning in dogs. Bush and Mosteller’s 1955 data follow 30 dogs learning to avoid a shock over 25 trials. Four models range from a logistic regression to a hierarchical version of the Bush–Mosteller learning model. Posterior predictive checks and cross-validation show that the data can’t clearly separate the mechanistic model from the flexible one. The conclusion is that they look the same.

Chapter 22: Cat adoptions. Do black cats wait longer to be adopted? Data on over 22,000 cats from an Austin, TX shelter include many that left by other routes, so their adoption time was never observed. Five models show that ignoring those censored cats badly overstates the adoption rate and hides the difference between colors, and that testing each model on simulated data first catches the error.

Chapter 23: World Cup football. A model of the 64 matches of the 2014 World Cup has perfect sampling diagnostics and plausible-looking estimates, and is still wrong because of a missing factor of 2 in one line of code. Only posterior predictive checks, which ask the model to reproduce the data it was trained on, expose it. The chapter ends with discrete Poisson models that suit goal counts naturally.

Chapter 24: Roaches. Cockroach counts from 262 apartments before and after pest control are used to teach leave-one-out cross-validation. Poisson, negative binomial, varying-intercept, and zero-inflated models each fix a problem the previous one exposed, and cross-validation rejects the plain Poisson model by thousands of elpd units.

Chapter 25: Golf putting. Starting from the geometry of a ball and a hole, the chapter derives the probability of sinking a putt at each distance. Eleven models are tested against professional golf data from two sources a decade apart. The newer, more detailed data expose what the original physical model left out, and each addition is guided by what the previous model’s errors show.

Chapter 26: Shark movement. GPS tags record where white sharks off South Africa go, but not what they’re doing. Hidden Markov models infer an unobserved behavioral state from step length and turning angle, then ask how time of day, sex, and bait in the water change how often sharks switch states. The final model gives each shark its own tendencies.

Chapter 27: Birthdays. Why are fewer US babies born on Sundays? Twenty years of daily birth counts are pulled apart into a long-term trend, a seasonal cycle, a weekday pattern that grows over time, and special days such as holidays. Eleven models each add one piece, and the choice of prior for the special-day effects changes how many holidays the model detects.

Chapter 28: Student grades. Predicting Portuguese students’ grades from 26 candidate variables tempts a model to fit noise. Four priors are compared by what they imply about the proportion of variance explained ($R^2$) before any data are seen. A variable-selection method then finds the small set of predictors that does nearly as well as all 26.

Chapter 29: No vehicles in the park. Over 51,000 survey responses judge whether things like an ambulance or a kite break a “no vehicles in the park” rule. The model is a standard hierarchical logistic regression, yet the sampler struggles. Four versions show that the problem is how the parameters are written, not the model itself; rewriting fixes it at no statistical cost.

Chapter 30: Planetary motion. Forty noisy positions of a simulated planet are used to estimate the strength of gravity. Five of eight chains get stuck in the wrong part of a rugged posterior and look stable while reporting wrong answers. Starting the chains from points found by a fast approximate method, Pathfinder, fixes it in about a second.

Lessons learned

Reproducing fifteen case studies taught me a few things that aren’t in any single chapter.

Chapters 23 and 18 are the clearest reminders that passing diagnostics doesn’t mean the model is good. In Ch 23, the sampler worked perfectly on a model with a missing factor of 2 in one line of code. In Ch 18, the distribution chosen specifically to match the data type came last in the LOO ranking. The diagnostic gate catches computational failures, while posterior predictive checks and cross-validation catch the modelling ones.

In some models the quantities you care about are computed from other, hidden quantities that the sampler moves through. The computed quantities can look healthy while the hidden ones are in trouble. Run diagnostics on what the sampler actually samples, not the derived output.

Several of my early write-ups described figures from the plotting code rather than from the rendered image, and some of those descriptions were wrong. Render the figure and look at it before writing up the report.

Eight of the fifteen chapters need models that brms can’t express: hidden Markov models, differential equations, custom likelihoods. For those I called Stan directly through cmdstanr and recorded each fit in the workflow state by hand. The workflow still applied, but manually.

Chapters 29 and 30 fail not because the models are wrong but because of how the parameters are arranged, or where the chains start. Changing the arrangement or the starting points fixes both at no statistical cost.

Limitations

The diagnostic gate is a guardrail, not a lock. Calling summary(fit) directly on the fitted model, or editing the workflow state by hand, bypasses it. It stops accidents; it doesn’t stop someone who’s trying to skip steps on purpose.

The AI layer is the least tested part. Language-model responses vary from one run to the next and can’t be verified by automated tests. Posit Assistant connects and can read the context file, but none of the fifteen case studies relied on it.

RBayesflow is not an R package, so updating means downloading a new version of the folder and letting renv restore the pinned packages. That keeps analyses reproducible, but it’s less convenient than install.packages().

Stan needs a working C++ compiler. On Windows, that means RTools, and getting it set up correctly is the most common point of failure. The installation script checks it first and stops with instructions if it’s broken.

When a model can’t be written in brms, phases have to be recorded by hand. Eight of the fifteen case studies hit this — hidden Markov models, differential equations, and custom likelihoods are all outside what brms expresses. The workflow still applies, but less of it runs automatically.

Expert Bayesian analysts will probably find RBayesflow confining. It imposes an order that experienced practitioners break on purpose, which is expected.

Conclusion

RBayesflow is for students and working scientists who can use R but haven’t yet built the habits that experienced Bayesian analysts carry around. It makes you ask whether you need a Bayesian model before building one. It keeps results from a failed sampler hidden until you’ve looked at the failure. It separates “did the computation work?” from “is the model any good?” and makes you answer both questions, in order. Every decision goes into the record.

If you already have a Bayesian workflow you trust, you don’t need it. If you’re working through the Bayesian Workflow book, or starting your first real Bayesian analysis, it gives you the book’s structure in a form you can run, with fifteen worked case studies alongside.

To try it:

git clone https://github.com/XerxesZorgon/RBayesflow.git

Then follow docs/user-guide/01-installation.md, run rbf_install(), and start with guide(wf).

BayesFlow

After publishing this article I found a Python package with a similar name, BayesFlow. It does simulation-based amortized Bayesian inference: neural networks are trained on simulations and reused for later inference. RBayesflow follows the Gelman et al. process of iterating through priors, fitting, diagnostics, and checks. The “R” in RBayesflow just means it’s for the R language.

RBayesflow Functions

Functions available in RBayesflow are grouped into core categories:

Setup and installation (R/install.R, new in v0.2.0)

• rbf_install() — checks R version, C++ toolchain, renv, CmdStan, and runs a Stan test; prints ✓/✗ for each of 7 steps. Has a dry_run argument that skips the actual install calls.
• rbf_new(name) — creates data/<name>/ with a .Rprofile, an empty wf_context.json, and a README.md.
• rbf_analysis_path() — a one-line helper that returns the current analysis subfolder if you’re inside one, or the project root otherwise. Mostly used internally by export_context().

Core workflow state (R/wf_state.R, R/init.R)

• init_workflow(mode, stage) — constructs the wf_state object and starts a session.
• new_wf_state(mode, stage, formula, data, family) — the S3 constructor itself. S3 is an object-oriented system for R - details are not important for RBayesflow users.
• print.wf_state(), summary.wf_state(), format.wf_state() — S3 methods that implement the display contract (the diagnostic gate lives here).

Diagnostics (R/diagnostics.R, R/diagnostic_registry.R)

• run_diagnostics(fit, wf_state) — runs Rhat/ESS/divergence/BFMI/treedepth checks plus the family-specific registry entry (Bernoulli, Poisson/NegBinomial, hierarchical, time-series).
• wf$diagnose() — shows the bayesplot diagnostic panels, prompts you to confirm you’ve reviewed them, and only then unlocks the coefficient output.
• detect_parameterization(fit) — returns "centered" or "non-centered".
• refit_noncentered(wf) — one-call switch to reparameterize a hierarchical model when divergences show up.

Off-ramps and exit (R/offramps.R, R/exit_workflow.R)

• assess_offramps(outcome_type, n, event_rate, goal) — returns the equal-weighted non-Bayesian alternatives.
• exit_workflow(wf, method, alternatives) — writes the YAML exit log and closes the stage.

Posit Assistant integration (R/context.R)

• export_context(wf, path) — serializes wf_state to JSON. As of v0.2.0 it defaults to writing into the current analysis subfolder rather than the project root.

The RBayesflows functions you are most likely to use in practice are:

• rbf_install()
• rbf_new()
• init_workflow()
• guide()
• wf$diagnose()
• assess_offramps() and exit_workflow()
• export_context()

To get help during a workflow use guide(wf) which reads wf$audit_trail and the diagnostic fields, figures out where you are in the seven-phase sequence, and prints three lines:

Current phase : <name and number>
What to expect: <what normal output looks like here>
Next step     : <exact function call or template to open>

It checks ten conditions in order:

Glossary

Audit trail. A timestamped list, kept inside the workflow state, of every decision made during an analysis.

Bayes’ theorem. The rule for updating beliefs with data: the posterior is proportional to the likelihood times the prior.

Beta-binomial distribution. A version of the binomial distribution in which the probability of success itself varies, allowing counts to be more spread out than a binomial permits.

brms. An R package that writes Stan models from R’s familiar regression formula syntax.

Calibration. A model is calibrated if its 90% prediction intervals contain the real value about 90% of the time.

Censoring. When an event such as an adoption hasn’t happened by the time observation stops, so you know only that it takes at least that long.

Chain. One run of the MCMC sampler. Stan normally runs four, starting from different random points.

cmdstanr. An R interface to Stan that runs models directly, used when brms can’t express a model.

Diagnostic gate. RBayesflow’s rule that results from a model with failed diagnostics are withheld until the user reviews and acknowledges the failure.

Divergence. A warning from Stan’s sampler that it hit a region of the posterior it couldn’t explore accurately. Even a few divergences can bias results.

elpd (expected log predictive density). A score for how well a model predicts data it hasn’t seen. Higher is better.

Effective sample size (ESS). How many independent draws a set of correlated MCMC samples is worth. RBayesflow requires at least 400.

Hierarchical model. A model for grouped data in which each group has its own parameters, drawn from a shared distribution.

Likelihood. How probable the observed data would be for a given value of the unknown parameters.

LOO-CV (leave-one-out cross-validation). Estimating predictive accuracy by predicting each observation from a model fit without it. The loo package approximates this without refitting the model hundreds of times.

LOO-PIT. A calibration check showing where each observation falls within its leave-one-out prediction. A flat plot means good calibration.

MCMC (Markov chain Monte Carlo). A family of algorithms that draw samples from a probability distribution by taking a guided random walk through parameter values.

Off-ramp. A non-Bayesian alternative offered by RBayesflow before any model is fit.

Overdispersion. When data are more spread out than a model’s distribution allows.

Parameterization. The particular way a model’s parameters are written. Two parameterizations can describe the same model but be much easier or harder for the sampler to explore.

Partial pooling. In a hierarchical model, each group’s estimate is pulled toward the overall average, more strongly for groups with little data.

Posterior. The distribution of plausible parameter values after combining the prior with the data.

Posterior predictive check. Simulating new data from a fitted model and comparing it with the real data to see what the model gets wrong.

Prior. The distribution of plausible parameter values before seeing the data.

Prior predictive simulation. Simulating data from the priors alone, before fitting, to check that the priors allow sensible data.

R-hat. A statistic comparing chains with each other. Values close to 1 mean the chains agree. RBayesflow requires R-hat below 1.01.

Stan. A programming language and software for fitting Bayesian models, mainly by MCMC.

wf_state. RBayesflow’s workflow-state object, holding the mode, priors, diagnostics, and audit trail for one analysis.

Software

• R - R is a language and environment for statistical computing and graphics.
• RStudio - RStudio is an integrated development environment (IDE) for R. It includes a console, syntax-highlighting editor that supports direct code execution, as well as tools for plotting, history, debugging and workspace management.

References

Bürkner, P.-C. (2017). brms: An R package for Bayesian multilevel models using Stan. Journal of Statistical Software, 80(1), 1–28.

Firth, D. (1993). Bias reduction of maximum likelihood estimates. Biometrika, 80(1), 27–38.

Gabry, J., Simpson, D., Vehtari, A., Betancourt, M., & Gelman, A. (2019). Visualization in Bayesian workflow. Journal of the Royal Statistical Society: Series A, 182(2), 389–402.

Gelman, A., Vehtari, A., McElreath, R., Simpson, D., Margossian, C. C., Yao, Y., Kennedy, L., Gabry, J., Bürkner, P.-C., Modrák, M., & Leos Barajas, V. (2026). Bayesian Workflow. Chapman & Hall. Companion website: avehtari.github.io/Bayesian-Workflow. Free electronic edition for non-commercial use: users.aalto.fi/~ave/Bayesian-Workflow.pdf.

Gelman, A., Vehtari, A., Simpson, D., Margossian, C. C., Carpenter, B., Yao, Y., Kennedy, L., Gabry, J., Bürkner, P.-C., & Modrák, M. (2020). Bayesian workflow. arXiv:2011.01808.

Lintzeris, N., et al. (2019). Nabiximols for the treatment of cannabis dependence: A randomized clinical trial. JAMA Internal Medicine, 179(9), 1242–1253.

Peach, J. (2026). RBayesflow [Software]. github.com/XerxesZorgon/RBayesflow.

Stan Development Team. Stan modeling language users guide and reference manual. mc-stan.org.

Vehtari, A., Gelman, A., & Gabry, J. (2017). Practical Bayesian model evaluation using leave-one-out cross-validation and WAIC. Statistics and Computing, 27(5), 1413–1432.

Vehtari, A., Gelman, A., Simpson, D., Carpenter, B., & Bürkner, P.-C. (2021). Rank-normalization, folding, and localization: An improved $\widehat{R}$ for assessing convergence of MCMC. Bayesian Analysis, 16(2), 667–718.

Image Credits

• Hero: Bayes Theorem, Mathematical Mysteries.