Translation cookbook

Linear regression

R:

fit <- lm(y ~ x1 + x2, data=df)

Python:

import lapylace as lp

fit = lp.stan_glm(
    "y ~ x1 + x2",
    data=df,
    family=lp.gaussian(),
    prior=lp.normal(0, 2.5),
    prior_intercept=lp.normal(0, 5),
    prior_aux=lp.exponential(1),
    chains=4,
    parallel_chains=4,
    iter_warmup=1000,
    iter_sampling=1000,
    seed=123,
)
fit.summary(["alpha", "beta", "sigma"])

Bayesian Stan version:

data { int<lower=0> N; vector[N] y; matrix[N, K] X; }
parameters { vector[K] beta; real<lower=0> sigma; }
model { y ~ normal(X * beta, sigma); beta ~ normal(0, 2.5); sigma ~ exponential(1); }

Logistic regression

fit = lp.stan_glm(
    "y ~ x1 + x2",
    data=df,
    family=lp.bernoulli(),
    prior=lp.normal(0, 2.5),
    prior_intercept=lp.normal(0, 5),
    chains=4,
    parallel_chains=4,
    iter_warmup=1000,
    iter_sampling=1000,
    seed=124,
)
fit.summary(["alpha", "beta"])

Stan/BlackJAX versions use the same log posterior, with unconstrained parameters and explicit priors.

Multilevel regression

Default to non-centered parameterizations:

alpha_j = mu_alpha + tau_alpha * z_j,   z_j ~ Normal(0, 1)
y_i ~ Normal(alpha[group_i] + x_i beta, sigma)

Use CmdStanPy for readable model-first exposition; use BlackJAX when the point is explicit sampler mechanics or JAX-native speed.