The intermediate representation

Every @brm block parses to a typed intermediate representation before anything is emitted for a backend. That representation is what the IR pane — the second pane of every generated model comparison on this site — shows, automatically and everywhere. This page explains how to read it.

From formula block to BRMI

@brm lowers a formula block to a BRMI: a NamedTuple of operations keyed by left-hand-side name. Each entry is either a statement (loc ~ 1 + age, y ~ Normal(loc, err)) or a data column with its element type and length. Priors are first-class operations too, addressed by effect(...) / sd(...) / cor(...) keys, so a model and its priors compose with plain Base.merge.

IR demo model

ir_demo_model = (@brm begin
    y ~ Normal(loc, err)
    loc ~ 1 + age + (1 | subj)
    err ~ Exponential(1)
end)((;
    age=[21.0, 38.0, 55.0, 29.0, 47.0, 61.0],
    subj=[1, 1, 2, 2, 3, 3],
    y=[0.2, 1.1, -0.4, 0.7, 1.4, 1.0],
))

BRMI:
  loc ~ 1 + age + (1 | subj)
  err ~ Exponential(1)
  y ~ Normal(loc, err)
  age: data (eltype=Float64, n=6)
  subj: data (eltype=Int64, n=6)

SBBRMI with data keys = [:age, :subj, :total_A_loc, :total_group_loc, :total_location_loc, :total_ng_loc, :total_nk_loc, :total_np_loc, :total_precision_loc, :y]
configured submodels:
_brm_total_scales_configured_1 = Base.merge(BayesianRegressionModels._brm_total_scales, quote
            tau ~ (ValueFamily(brm_vector_prior_faeb6f6956d0662a))(0.0, 1.0; n = 1)
        end)
emitted @slic body:
begin
    err ~ exponential(1.0 ./ 1)
    total_scale_loc ~ _brm_total_scales_configured_1(; n = total_nk_loc)
    total_loc::matrix[total_ng_loc, total_nk_loc] ~ brm_total(total_scale_loc, total_A_loc, total_location_loc, total_precision_loc)
    population_loc = brm_total_recover_rng(total_loc, total_scale_loc, total_A_loc, total_location_loc, total_precision_loc)
    deviation_loc = brm_total_deviations(total_loc, total_A_loc * population_loc)
    total_Z_loc = hcat(rep_vector(1.0, num_elements(total_group_loc)))
    X_loc = hcat(age)
    pop_loc ~ popefs(; X = X_loc)
    loc = rows_dot_product(total_loc[total_group_loc, :], total_Z_loc) + pop_loc
    y ~ normal(loc, err)
end

functions {
// value UDF brm_vector_prior_faeb6f6956d0662a_lpdf
real brm_vector_prior_faeb6f6956d0662a_lpdf(
    vector x,
    real arg_1,
    real arg_2
) {
    if((x[1] < 0.0)) {
        return negative_infinity();
    }
    return lognormal_lpdf(x[1] | arg_1, arg_2);
}
real brm_total_lpdf(
    matrix total,
    vector tau,
    matrix A,
    vector location,
    vector precision
) {
    int j = dims(total)[1];
    int k = dims(total)[2];
    int p = dims(A)[2];
    if (dims(tau)[1] != k) reject("brm_total_lpdf: dim mismatch — `tau` dim 1 (= ", dims(tau)[1], ") does not match `k` (= ", k, "), inferred from `total` dim 2. `k` sizes: `total` dim 2 (= ", dims(total)[2], "), `tau` dim 1 (= ", dims(tau)[1], "), `A` dim 1 (= ", dims(A)[1], ").");
    if (dims(A)[1] != k) reject("brm_total_lpdf: dim mismatch — `A` dim 1 (= ", dims(A)[1], ") does not match `k` (= ", k, "), inferred from `total` dim 2. `k` sizes: `total` dim 2 (= ", dims(total)[2], "), `tau` dim 1 (= ", dims(tau)[1], "), `A` dim 1 (= ", dims(A)[1], ").");
    if (dims(location)[1] != p) reject("brm_total_lpdf: dim mismatch — `location` dim 1 (= ", dims(location)[1], ") does not match `p` (= ", p, "), inferred from `A` dim 2. `p` sizes: `A` dim 2 (= ", dims(A)[2], "), `location` dim 1 (= ", dims(location)[1], "), `precision` dim 1 (= ", dims(precision)[1], ").");
    if (dims(precision)[1] != p) reject("brm_total_lpdf: dim mismatch — `precision` dim 1 (= ", dims(precision)[1], ") does not match `p` (= ", p, "), inferred from `A` dim 2. `p` sizes: `A` dim 2 (= ", dims(A)[2], "), `location` dim 1 (= ", dims(location)[1], "), `precision` dim 1 (= ", dims(precision)[1], ").");
    matrix[dims(precision)[1], dims(precision)[1]] Q = brm_total_precision(tau, A, precision, j);
    vector[dims(total)[2]] average = brm_total_mean(total);
    vector[dims(precision)[1]] beta = brm_total_conditional_mean(total, tau, A, location, precision, Q);
    vector[dims(total)[2]] residual = (average - (A * beta));
    real quadratic = 0.0;
    real lp = ((-0.5 * ((j * k) - p) * 1.8378770664093453) - (j * sum(log(tau))));
    for(c in 1:k) {
        quadratic += ((j * square(residual[c])) / square(tau[c]));
        for(g in 1:j) {
            quadratic += (square((total[g, c] - average[c])) / square(tau[c]));
        }
    }
    for(a in 1:p) {
        if((precision[a] > 0.0)) {
            lp += (0.5 * (log(precision[a]) - 1.8378770664093453));
            quadratic += (precision[a] * square((beta[a] - location[a])));
        }
    }
    return (lp - (0.5 * (log_determinant(Q) + quadratic)));
}
matrix brm_total_precision(
    vector tau,
    matrix A,
    vector precision,
    int n_groups
) {
    int k = dims(tau)[1];
    int p = dims(A)[2];
    if (dims(A)[1] != k) reject("brm_total_precision: dim mismatch — `A` dim 1 (= ", dims(A)[1], ") does not match `k` (= ", k, "), inferred from `tau` dim 1. `k` sizes: `tau` dim 1 (= ", dims(tau)[1], "), `A` dim 1 (= ", dims(A)[1], ").");
    if (dims(precision)[1] != p) reject("brm_total_precision: dim mismatch — `precision` dim 1 (= ", dims(precision)[1], ") does not match `p` (= ", p, "), inferred from `A` dim 2. `p` sizes: `A` dim 2 (= ", dims(A)[2], "), `precision` dim 1 (= ", dims(precision)[1], ").");
    matrix[dims(precision)[1], dims(precision)[1]] out = diag_matrix(precision);
    for(a in 1:p) {
        for(b in 1:p) {
            for(c in 1:k) {
                out[a, b] += ((n_groups * A[c, a] * A[c, b]) / square(tau[c]));
            }
        }
    }
    return out;
}
vector brm_total_mean(
    matrix total
) {
    int j = dims(total)[1];
    int k = dims(total)[2];
    vector[k] out = rep_vector(0.0, k);
    for(c in 1:k) {
        out[c] = (sum(total[:, c]) / j);
    }
    return out;
}
vector brm_total_conditional_mean(
    matrix total,
    vector tau,
    matrix A,
    vector location,
    vector precision,
    matrix Q
) {
    int j = dims(total)[1];
    int k = dims(total)[2];
    int p = dims(A)[2];
    if (dims(tau)[1] != k) reject("brm_total_conditional_mean: dim mismatch — `tau` dim 1 (= ", dims(tau)[1], ") does not match `k` (= ", k, "), inferred from `total` dim 2. `k` sizes: `total` dim 2 (= ", dims(total)[2], "), `tau` dim 1 (= ", dims(tau)[1], "), `A` dim 1 (= ", dims(A)[1], ").");
    if (dims(A)[1] != k) reject("brm_total_conditional_mean: dim mismatch — `A` dim 1 (= ", dims(A)[1], ") does not match `k` (= ", k, "), inferred from `total` dim 2. `k` sizes: `total` dim 2 (= ", dims(total)[2], "), `tau` dim 1 (= ", dims(tau)[1], "), `A` dim 1 (= ", dims(A)[1], ").");
    if (dims(location)[1] != p) reject("brm_total_conditional_mean: dim mismatch — `location` dim 1 (= ", dims(location)[1], ") does not match `p` (= ", p, "), inferred from `A` dim 2. `p` sizes: `A` dim 2 (= ", dims(A)[2], "), `location` dim 1 (= ", dims(location)[1], "), `precision` dim 1 (= ", dims(precision)[1], "), `Q` dim 1 (= ", dims(Q)[1], "), `Q` dim 2 (= ", dims(Q)[2], ").");
    if (dims(precision)[1] != p) reject("brm_total_conditional_mean: dim mismatch — `precision` dim 1 (= ", dims(precision)[1], ") does not match `p` (= ", p, "), inferred from `A` dim 2. `p` sizes: `A` dim 2 (= ", dims(A)[2], "), `location` dim 1 (= ", dims(location)[1], "), `precision` dim 1 (= ", dims(precision)[1], "), `Q` dim 1 (= ", dims(Q)[1], "), `Q` dim 2 (= ", dims(Q)[2], ").");
    if (dims(Q)[1] != p) reject("brm_total_conditional_mean: dim mismatch — `Q` dim 1 (= ", dims(Q)[1], ") does not match `p` (= ", p, "), inferred from `A` dim 2. `p` sizes: `A` dim 2 (= ", dims(A)[2], "), `location` dim 1 (= ", dims(location)[1], "), `precision` dim 1 (= ", dims(precision)[1], "), `Q` dim 1 (= ", dims(Q)[1], "), `Q` dim 2 (= ", dims(Q)[2], ").");
    if (dims(Q)[2] != p) reject("brm_total_conditional_mean: dim mismatch — `Q` dim 2 (= ", dims(Q)[2], ") does not match `p` (= ", p, "), inferred from `A` dim 2. `p` sizes: `A` dim 2 (= ", dims(A)[2], "), `location` dim 1 (= ", dims(location)[1], "), `precision` dim 1 (= ", dims(precision)[1], "), `Q` dim 1 (= ", dims(Q)[1], "), `Q` dim 2 (= ", dims(Q)[2], ").");
    vector[dims(total)[2]] average = brm_total_mean(total);
    vector[dims(location)[1]] natural = (precision .* location);
    for(a in 1:p) {
        for(c in 1:k) {
            natural[a] += ((j * A[c, a] * average[c]) / square(tau[c]));
        }
    }
    return mdivide_left_spd(Q, natural);
}
vector brm_total_recover_rng(
    matrix total,
    vector tau,
    matrix A,
    vector location,
    vector precision
) {
    int j = dims(total)[1];
    int k = dims(total)[2];
    int p = dims(A)[2];
    if (dims(tau)[1] != k) reject("brm_total_recover_rng: dim mismatch — `tau` dim 1 (= ", dims(tau)[1], ") does not match `k` (= ", k, "), inferred from `total` dim 2. `k` sizes: `total` dim 2 (= ", dims(total)[2], "), `tau` dim 1 (= ", dims(tau)[1], "), `A` dim 1 (= ", dims(A)[1], ").");
    if (dims(A)[1] != k) reject("brm_total_recover_rng: dim mismatch — `A` dim 1 (= ", dims(A)[1], ") does not match `k` (= ", k, "), inferred from `total` dim 2. `k` sizes: `total` dim 2 (= ", dims(total)[2], "), `tau` dim 1 (= ", dims(tau)[1], "), `A` dim 1 (= ", dims(A)[1], ").");
    if (dims(location)[1] != p) reject("brm_total_recover_rng: dim mismatch — `location` dim 1 (= ", dims(location)[1], ") does not match `p` (= ", p, "), inferred from `A` dim 2. `p` sizes: `A` dim 2 (= ", dims(A)[2], "), `location` dim 1 (= ", dims(location)[1], "), `precision` dim 1 (= ", dims(precision)[1], ").");
    if (dims(precision)[1] != p) reject("brm_total_recover_rng: dim mismatch — `precision` dim 1 (= ", dims(precision)[1], ") does not match `p` (= ", p, "), inferred from `A` dim 2. `p` sizes: `A` dim 2 (= ", dims(A)[2], "), `location` dim 1 (= ", dims(location)[1], "), `precision` dim 1 (= ", dims(precision)[1], ").");
    matrix[dims(precision)[1], dims(precision)[1]] Q = brm_total_precision(tau, A, precision, j);
    vector[dims(precision)[1]] beta = brm_total_conditional_mean(total, tau, A, location, precision, Q);
    return multi_normal_rng(beta, inverse_spd(Q));
}
matrix brm_total_deviations(
    matrix total,
    vector mu
) {
    int j = dims(total)[1];
    int k = dims(total)[2];
    if (dims(mu)[1] != k) reject("brm_total_deviations: dim mismatch — `mu` dim 1 (= ", dims(mu)[1], ") does not match `k` (= ", k, "), inferred from `total` dim 2. `k` sizes: `total` dim 2 (= ", dims(total)[2], "), `mu` dim 1 (= ", dims(mu)[1], ").");
    matrix[dims(total)[1], dims(total)[2]] out = total;
    for(c in 1:k) {
        out[:, c] = (total[:, c] - rep_vector(mu[c], j));
    }
    return out;
}
matrix hcat(vector x) {
    int n = dims(x)[1];
    return to_matrix(x, n, 1);
}
vector normal_lpdfs(
    vector obs,
    vector loc,
    real scale
) {
    return jbroadcasted_normal_lpdfs(obs, loc, scale);
}
vector jbroadcasted_normal_lpdfs(
    vector x1,
    vector x2,
    real x3
) {
    int n = dims(x1)[1];
    vector[n] rv;
    for(i in 1:n) {
        rv[i] = normal_lpdfs(broadcasted_getindex(x1, i), broadcasted_getindex(x2, i), x3);
    }
    return rv;
}
real normal_lpdfs(
    real args1,
    real args2,
    real args3
) {
    return normal_lpdf(args1 | args2, args3);
}
real broadcasted_getindex(vector x, int i) {
    return x[i];
}
vector normal_vector_rng(
    int anontok__1,
    vector a,
    real b
) {
    int n = anontok__1;
    if((n == 0)) {
        vector[n] rv;
        return rv;
    } else {
        return to_vector(normal_rng(a, b));
    }
}
}
data {
    int total_ng_loc;
    int total_nk_loc;
    int total_A_loc_m;
    int total_A_loc_n;
    matrix[total_A_loc_m, total_A_loc_n] total_A_loc;
    int total_location_loc_n;
    vector[total_location_loc_n] total_location_loc;
    int total_precision_loc_n;
    vector[total_precision_loc_n] total_precision_loc;
    int total_group_loc_n;
    array[total_group_loc_n] int total_group_loc;
    int age_n;
    vector[age_n] age;
    int y_n;
    vector[y_n] y;
}
transformed data {
    matrix[num_elements(total_group_loc), 1] total_Z_loc = hcat(rep_vector(1.0, num_elements(total_group_loc)));
    matrix[age_n, 1] X_loc = hcat(age);
    int pop_loc_n_covariates = 1;
}
parameters {
    real<lower=0.0> err;
    vector<lower=0.0>[1] total_scale_loc_tau;
    matrix[total_ng_loc, total_nk_loc] total_loc;
    vector[pop_loc_n_covariates] pop_loc_beta_pop;
}
transformed parameters {
    vector<lower=0.0>[1] total_scale_loc = total_scale_loc_tau;
    vector[age_n] pop_loc = (X_loc * pop_loc_beta_pop);
    vector[num_elements(total_group_loc)] loc = (rows_dot_product(total_loc[total_group_loc, :], total_Z_loc) + pop_loc);
}
model {
    err ~ exponential((1.0 ./ 1));
    total_scale_loc_tau ~ brm_vector_prior_faeb6f6956d0662a(0.0, 1.0);
    total_loc ~ brm_total(total_scale_loc, total_A_loc, total_location_loc, total_precision_loc);
    pop_loc_beta_pop ~ std_normal();
    y ~ normal(loc, err);
}
generated quantities {
    vector[total_precision_loc_n] population_loc = brm_total_recover_rng(
        total_loc,
        total_scale_loc,
        total_A_loc,
        total_location_loc,
        total_precision_loc
    );
    matrix[total_ng_loc, total_A_loc_m] deviation_loc = brm_total_deviations(total_loc, (total_A_loc * population_loc));
    vector[y_n] y_likelihood = normal_lpdfs(y, loc, err);
    vector[y_n] y_gen = normal_vector_rng(y_n, loc, err);
}

#= line 0 =# Turing.@model(function brm_model(y, X_loc, group_effects_loc_1)
        err ~ Distributions.Exponential(1)
        beta_pop ~ Distributions.product_distribution(Base.fill(Distributions.Normal(), 2))
        eta_loc = X_loc * beta_pop
        group_effect_1 = Base.zeros(Base.length(y))
        group_1_1 ~ DynamicPPL.to_submodel(BRM.turing_group_effect(group_effects_loc_1, (nothing,), nothing))
        group_effect_1 = group_effect_1 + group_1_1.effect
        eta_loc = eta_loc + group_effect_1
        loc = eta_loc
        begin
            for i = Base.eachindex(y)
                y[i] ~ Distributions.Normal(loc[i], err)
            end
        end
        (; loc = loc, err = err, response = y)
    end)

Read the IR pane top to bottom: statements first, then the data each statement consumes. Grouping parentheses stay where Julia's precedence needs them ((1 | subj) binds looser than ~), redundant ones are dropped.

The lowering chain

One BRMI feeds every backend through two further stages:

• VBRMI — the prepared, modification-frozen model: prior addresses matched, group declarations resolved, fitted preprocessing/replay recorded. Its summary prints the flat parameter dimension plus the materialized columns and blocks.

• SBBRMI / TuringBRMI — the backend lowerings shown in the remaining panes: the StanBlocks model and Stan source on one side, the native Turing model on the other.

Nothing in the chain re-parses the formula: backends consume the same operations, so the IR pane and the emitted panes cannot drift apart — they are rendered from one object in one build step.

Introspection entry points

The API documents the query surface over these types, including prior lookup (term_priors, effect_priors) and the brm_descriptor views each backend exposes. Start from the IR pane when an emitted backend surprises you: if a statement or data column is wrong there, it is wrong in every pane below it.