| Title: | Bayesian Synthetic Control Models |
| Version: | 1.0.1 |
| Description: | Implements the synthetic control method of Abadie, Diamond, and Hainmueller (2010) <doi:10.1198/jasa.2009.ap08746> within a Bayesian framework, enabling straightforward uncertainty quantification of treatment effects and other quantities of interest. Supports time-varying covariates with potentially time-varying effects, single or multiple treated units, and staggered treatment adoption. Provides methods for model assessment, comparison, and selection based on placebo studies, cross-validation, and posterior predictive checks. Posterior sampling is performed using Markov chain Monte Carlo via Stan. |
| License: | GPL (≥ 3) |
| Encoding: | UTF-8 |
| RoxygenNote: | 7.3.3 |
| Biarch: | true |
| Depends: | R (≥ 4.1.0) |
| Imports: | checkmate, cli, dplyr, ggplot2, loo, methods, posterior, progressr, projpred, quadprog, Rcpp (≥ 0.12.0), RcppParallel (≥ 5.0.1), rlang, rstan (≥ 2.32.7), rstantools (≥ 2.6.0), tidyr |
| LinkingTo: | BH (≥ 1.66.0), Rcpp (≥ 0.12.0), RcppEigen (≥ 0.3.3.3.0), RcppParallel (≥ 5.0.1), rstan (≥ 2.32.7), StanHeaders (≥ 2.32.10) |
| SystemRequirements: | GNU make |
| LazyData: | true |
| LazyDataCompression: | xz |
| URL: | https://github.com/helske/bscm |
| BugReports: | https://github.com/helske/bscm/issues |
| Suggests: | testthat (≥ 3.0.0) |
| Config/testthat/edition: | 3 |
| NeedsCompilation: | yes |
| Packaged: | 2026-07-31 12:11:26 UTC; jvhels |
| Author: | Jouni Helske |
| Maintainer: | Jouni Helske <jouni.helske@iki.fi> |
| Repository: | CRAN |
| Date/Publication: | 2026-08-09 06:40:02 UTC |
The 'bscm' package
Description
Implements the synthetic control method of Abadie, Diamond, and Hainmueller (2010) within a Bayesian framework, enabling straightforward uncertainty quantification of treatment effects and other quantities of interest. Posterior sampling is performed using Markov chain Monte Carlo via Stan.
bscm supports time-varying covariates with potentially time-varying
effects, single treated unit as well as multiple treated units possibly with
staggered treatment adoption.
Output of bscm() is a bscmfit object with methods for model assessment,
comparison, and selection based on placebo studies, cross-validation, and
posterior predictive checks.
Author(s)
Maintainer: Jouni Helske jouni.helske@iki.fi (ORCID)
References
Abadie A, Diamond A, and Hainmueller J (2010). Synthetic Control Methods for Comparative Case Studies: Estimating the Effect of California’s Tobacco Control Program. Journal of the American Statistical Association, 105(490), 493–505, https://doi.org/10.1198/jasa.2009.ap08746.
Stan Development Team (2025). RStan: the R interface to Stan. R package version 2.32.7. https://mc-stan.org.
See Also
Useful links:
Extract posterior draws of model parameters as a data frame
Description
Returns a data.frame representation of the posterior sample of the model
parameters. For samples from posterior predictive distribution,
see posterior_epred() and posterior_predict(). For donor
weights, use donor_weights().
Usage
## S3 method for class 'bscmfit'
as.data.frame(
x,
row.names = NULL,
optional = FALSE,
parameters = NULL,
include = TRUE,
...
)
Arguments
x |
[ |
row.names |
Ignored. |
optional |
Ignored. |
parameters |
[ |
include |
[ |
... |
Ignored. |
Value
A data.frame containing model parameters in a wide format.
See Also
Examples
head(as.data.frame(fit_single_treated, parameters = c("alpha", "sigma")))
Return posterior draws from BSCM fit as a draws object
Description
Return posterior draws from BSCM fit as a draws object
Usage
## S3 method for class 'bscmfit'
as_draws(x, parameters = NULL, inc_warmup = FALSE, include = TRUE, ...)
Arguments
x |
[ |
parameters |
[ |
inc_warmup |
[ |
include |
[ |
... |
Ignored. |
Value
An object of class draws_array containing the posterior draws of
the specified parameters.
See Also
posterior::as_draws() for converting the output to other formats,
and posterior::summarise_draws() for computing posterior summaries of the
draws.
Examples
head(as_draws(fit_single_treated, parameters = c("alpha", "sigma")))
Pre- and post-treatment average effects of a Bayesian synthetic control model
Description
Pre- and post-treatment average effects of a Bayesian synthetic control model
Usage
average_treatment_effect(x, ...)
## S3 method for class 'bscmfit'
average_treatment_effect(
x,
average = TRUE,
summary = TRUE,
probs = c(0.025, 0.975),
...
)
Arguments
x |
[ |
... |
Ignored. |
average |
[ |
summary |
[ |
probs |
[ |
Value
A tibble of posterior summaries (summary = TRUE) or
posterior draws (summary = FALSE) in long format.
Examples
average_treatment_effect(fit_single_treated)
Estimate Bayesian synthetic control model
Description
Function bscm estimates a Bayesian synthetic control models of varying
types.
Usage
bscm(
formula,
data,
treatment = "treatment",
time = "time",
unit = "id",
error = "iid",
priors = NULL,
prior_only = FALSE,
mcmc_diagnostics = TRUE,
save_data = TRUE,
compute_predictions = TRUE,
...
)
Arguments
formula |
[ |
data |
[ |
treatment |
[ |
time |
[ |
unit |
[ |
error |
[ |
priors |
[ |
prior_only |
[ |
mcmc_diagnostics |
[ |
save_data |
[ |
compute_predictions |
[ |
... |
Additional parameters passed on to |
Details
To define formula in case of no predictors, use
outcome ~ 1 or outcome ~ 0. In the former case, as well by the
default when using predictors (e.g., outcome ~ x + z), the model
includes intercept term for each treated unit. Intercept can be omitted by
using 0 in the RHS, e.g., outcome ~ 0 or outcome ~ 0 + x + z
(equivalently, you can use -1 in place of 0). Formula should not contain
the variable defining the treatment, which is defined separately using the
argument treatment. In case the variable is present also in the
formula, it is automatically removed.
To specify predictors with time-varying coefficients in formula, wrap them
in tv(), e.g., outcome ~ x + tv(~ z, df = 10, type = "rw1") defines a
model where x has a time-constant coefficient and z has a time-varying
coefficient following penalized cubic spline with 10 spline basis
functions and random walk prior on the spline coefficients.
Terms inside tv() are automatically also included in the
time-constant part of the model, since the time-varying coefficients
are defined to have zero mean in the pre-treatment period (in case of
multiple treated units, minimum period). Note that a common time-varying
intercept is omitted as it would cancel out in the linear predictor.
Both the time-constant and time-varying regression part is assumed to apply
for all treated and donor units with common coefficients. If you
want to apply a covariate only for treated units, just set the covariate
values to zero for donors. This approach can be also used to vary behaviour
of intercept: Defining formula such as y ~ 0 + intercept + x, where
intercept is a name of constant column in the data, will define common
intercept for all units, instead of unit-specific intercepts
(e.g., fixed effects).
When model contains covariates X, their effect is subtracted from
donors, i.e., for treated unit i,
y_i \sim N(\alpha_i + X_i\beta + Z^\ast\omega_i, \sigma_i^2),
where column j of Z^\ast is
z_j - (\alpha_j + X_j \beta), i.e., the original donor outcome
vector z_j minus the donor-specific intercept \alpha_j and
the effect of covariates X_j \beta for that donor. Note that
\beta is common across donors and the treated units.
Model can contain missing values in the outcome variable of the treated, but not of the donors, nor in the covariates. Missing outcomes are automatically imputed during MCMC sampling under a missing at random (MAR) assumption.
Value
An object of class bscmfit.
See Also
summary.bscmfit(), as_draws.bscmfit(), rstan::sampling().
Examples
# skip diagnostics and use small number of iterations for CRAN checks
fit <- bscm(
y ~ 1, single_treated, "treatment", "time", "id",
priors = list(omega = dirichlet_pr(0.5)),
chains = 1, cores = 1, refresh = 0, iter = 1000,
mcmc_diagnostics = FALSE
)
fit
Prior distributions in the Bayesian synthetic control model
Description
The functions normal_pr(), student_pr(), exponential_pr(), gamma_pr(),
beta_pr(), half_normal_pr(), dirichlet_pr(), and logistic_normal_pr()
specify prior distributions for parameters of the Bayesian synthetic
control model estimated by bscm().
Usage
normal_pr(location, scale)
student_pr(df, location, scale)
gamma_pr(shape, rate)
exponential_pr(rate)
beta_pr(shape1, shape2)
half_normal_pr(scale)
dirichlet_pr(concentration)
logistic_normal_pr(scale)
Arguments
location |
The location parameter for the prior. |
scale |
The scale parameter for the prior. |
df |
The degrees of freedom for the Student's t prior. |
shape |
The shape parameter for the gamma prior. |
rate |
The rate parameter for the gamma and exponential priors. |
shape1 |
The first shape parameter for the beta prior. |
shape2 |
The second shape parameter for the beta prior. |
concentration |
The concentration parameter for the Dirichlet prior. |
Details
For parameters that can take any real value (alpha, beta),
currently supported prior distributions are normal_pr(location, scale) and
student_pr(df, location, scale). Note however that for the intercept the
prior corresponds to a value of the intercept when donors and predictors
corresponding to time-constant coefficients are centered.
For parameters constrained to be positive (sigma, kappa),
supported distributions are exponential_pr(rate), gamma_pr(shape, rate),
and half_normal_pr(scale).
For autoregressive coefficients of the residuals (rho) the prior corresponds
to transformation 0.5 * (1 + rho), so that the support is (0, 1). Only
supported prior distribution for this is beta_pr(shape1, shape2).
For the donor weight vector omega, supported priors are
dirichlet_pr(concentration) and logistic_normal_pr(scale). The
(symmetric) Dirichlet prior with concentration \alpha_\omega defined as
\omega \sim \text{Dirichlet}(\alpha_\omega, \ldots,\alpha_\omega).
Values \alpha_\omega < 1 concentrate weight on few donors, while
\alpha_\omega = 1 is uniform over the simplex, and
\alpha_\omega > 1 pulls weights toward the center of the simplex.
Logistic normal with with scale \sigma_\alpha is defined as
\omega = \text{softmax}(\eta),
where \eta \sim N(0, \sigma_\alpha^2 I) constrained to sum to zero.
Larger scale induces more concentrated (sparser) weights.
Value
A bscm_prior object.
See Also
Examples
normal_pr(c(0, 0), c(1, 2))
student_pr(df = 4, location = 0, scale = 2)
exponential_pr(rate = c(1, 0.1))
gamma_pr(2, 1)
beta_pr(2, 2)
half_normal_pr(1)
dirichlet_pr(1)
logistic_normal_pr(2)
Check the validity of the posterior of bscmfit object
Description
This function is automatically called at the end of bscm() to check that
the output can be trusted in terms of convergence of the MCMC sampling.
Checks consists of the common diagnostics of Hamiltonian Monte Carlo variant
used by Stan, as well as the Rhat values and effective sample sizes of model
parameters. See rstan::check_hmc_diagnostics() and
posterior::default_convergence_measures() for details on the definitions
of these.
Usage
check_mcmc_diagnostics(x, ...)
## S3 method for class 'bscmfit'
check_mcmc_diagnostics(x, warn = TRUE, ...)
Arguments
x |
[ |
... |
Ignored. |
warn |
[ |
Value
Invisibly returns a list containing the results of the check.
References
https://mc-stan.org/learn-stan/diagnostics-warnings.html
Examples
check_mcmc_diagnostics(fit_single_treated, warn = FALSE)
Extract regression coefficients of a Bayesian synthetic control model
Description
Extract regression coefficients of a Bayesian synthetic control model
Usage
## S3 method for class 'bscmfit'
coef(object, type = NULL, summary = TRUE, probs = c(0.025, 0.975), ...)
Arguments
object |
[ |
type |
[ |
summary |
[ |
probs |
[ |
... |
Ignored. |
Value
A tibble of posterior summaries (summary = TRUE) or
list of tibbles of posterior draws (summary = FALSE).
Examples
coef(fit_single_treated)
coef(fit_single_treated, type = "beta")
Covariate adjustments for Bayesian synthetic control model
Description
Computes the posterior distribution of the covariate adjustment to the synthetic control for each predictor, defined as the regression coefficient multiplied by the difference between the predictor value of the treated unit and the corresponding weighted predictor value of the donor units.
Usage
covariate_adjustment(x, ...)
## S3 method for class 'bscmfit'
covariate_adjustment(
x,
plot = TRUE,
probs = c(0.025, 0.975),
alpha = 0.5,
scales = "free_y",
...
)
Arguments
x |
[ |
... |
Ignored. |
plot |
[ |
probs |
[ |
alpha |
[ |
scales |
[ |
Details
For predictor k, treated unit i, and time point t, the
covariate adjustment is
b_{tk}
\left(
X_{y,itk} - \sum_j \omega_{ji} X_{z,jtk}
\right),
where b_{tk} = \beta_k for predictors with a fixed coefficient and
b_{tk} = \beta_k + \gamma_{tk} for predictors with a time-varying
coefficient.
If there are multiple treated units, the covariate adjustments are averaged over treated units within each posterior draw before computing posterior summaries.
Value
A data frame containing the posterior mean and requested posterior
quantiles of the covariate adjustment for each predictor and time point.
If plot = TRUE, the corresponding plot is printed as a side effect.
Covariate imbalance of Bayesian synthetic control model
Description
For models with covariates, returns and optionally visualizes the covariate imbalances
\delta_{t} =
\sqrt{\frac{1}{K}\sum_{k=1}^K(x_{k,0,t} - \bar x_{k,0,t})^2},
t=1,\ldots,T, where
\bar x_{k,0,t} = \sum_{j=1}^J \omega_j x_{k, j, t} and x_{k,0,t}
is the value of kth covariate of a treated unit at time t, and
similarly for donors j=1,\ldots,J. This is computed separately for
each treated unit in case of multiple treated units.
Usage
covariate_imbalance(x, ...)
## S3 method for class 'bscmfit'
covariate_imbalance(x, plot = TRUE, probs = c(0.025, 0.975), ...)
Arguments
x |
[ |
... |
Optional arguments passed to |
plot |
[ |
probs |
[ |
Value
A data.frame of posterior summaries of synthetic covariate
distances.
Examples
covariate_imbalance(fit_single_treated, plot = TRUE, probs = c(0.05, 0.95))
Donor ranking
Description
Returns ranking of donors either using posterior means of donor weights from
bscmfit object, or ranking from a vsel object
returned by projpred::varsel() or projpred::cv_varsel() when applied to
a bscmfit object.
Usage
donor_ranking(x)
Arguments
x |
[ |
Value
Character vector of donor IDs in order of selection.
Extract donor weights of a Bayesian synthetic control model
Description
Extract donor weights of a Bayesian synthetic control model
Usage
donor_weights(x, ...)
## S3 method for class 'bscmfit'
donor_weights(x, summary = TRUE, probs = c(0.025, 0.5, 0.975), ...)
Arguments
x |
[ |
... |
Ignored. |
summary |
[ |
probs |
[ |
Value
A tibble of posterior summaries (summary = TRUE) or
posterior draws (summary = FALSE) in long format.
Examples
donor_weights(fit_single_treated) |> head(5)
Extract the number of effective donors in a Bayesian synthetic control model
Description
Effective donors is defined as 1 / \sum_{j=1}^J \omega_j^2,
where \omega_j is the donor weight of control unit j.
Usage
effective_donors(x, ...)
## S3 method for class 'bscmfit'
effective_donors(
x,
average = FALSE,
summary = TRUE,
probs = c(0.025, 0.975),
...
)
Arguments
x |
[ |
... |
Ignored. |
average |
[ |
summary |
[ |
probs |
[ |
Value
A data.frame of posterior summaries (summary = TRUE) or
posterior draws (summary = FALSE) in long format.
Examples
effective_donors(fit_single_treated)
Example bscmfit object
Description
The object fit_single_treated contains a Bayesian synthetic
control model estimated as
fit <- bscm( formula = y ~ x, data = single_treated, treatment = "treatment", chains = 2, cores = 1, refresh = 0, iter = 2000, warmup = 1000, control = list(adapt_delta = 0.8) )
Format
A bscmfit object.
Examples
fit_single_treated
plot(fit_single_treated)
Expected values of posterior predictive distribution of Bayesian synthetic control model
Description
Returns posterior draws (or summaries) of expected values y_mean of the
posterior predictive distribution of the model.
Usage
## S3 method for class 'bscmfit'
fitted(object, summary = TRUE, probs = c(0.025, 0.975), ...)
Arguments
object |
[ |
summary |
[ |
probs |
[ |
... |
Ignored. |
Value
A data.frame of posterior summaries (summary = TRUE) or
posterior draws (summary = FALSE) in long format.
Examples
fitted(fit_single_treated) |> head()
Get prior specifications used in the Bayesian synthetic control model
Description
Get prior specifications used in the Bayesian synthetic control model
Usage
get_priors(x, ...)
## S3 method for class 'bscmfit'
get_priors(x, ...)
Arguments
x |
[ |
... |
Ignored. |
Value
A named list of bscm_prior objects for the parameters of the model.
Create a projpred reference model from a BSCM fit
Description
Creates a refmodel object from bscmfit to be used with projpred
package. This enables the usage of projpred::varsel() and
projpred::cv_varsel() for donor selection. Note that for predictions with
newdata, you need to call proj_predict_bscm() which is a wrapper of
projpred::proj_predict() to which converts the data to the wide format
used by projpred.
Usage
## S3 method for class 'bscmfit'
get_refmodel(object, ...)
Arguments
object |
A |
... |
Additional arguments passed to |
Details
Projection are based on only on the pretreatment period.
This function is experimental and currently only a single treated unit is supported.
For projpred integration, donors are treated as separate "predictors" in the formula. Currently only supported model is one without extra predictors. There are also other restrictions, namely lack of support for K-fold.
Value
An object of class refmodel.
See Also
projpred::varsel(), projpred::cv_varsel()
Examples
fit <- bscm(
y ~ 1, treatment = "treatment", time = "time", unit = "id",
data = single_treated, refresh = 0, chains = 2, cores = 2
)
refmodel <- get_refmodel(fit)
vs <- projpred::varsel(refmodel)
plot(vs, stats = c("elpd", "rmse"), alpha = 0.05)
# predictions using the projection:
predictions <- projpred::proj_predict(vs)
# posterior mean of the predictions
colMeans(predictions)
Get input data to Stan from bscmfit object
Description
Reconstructs the data list passed to Stan in bscm(). Used by methods such
as loo() with reloo = TRUE.Requires the original data
(save_data = TRUE).
Usage
get_standata(x, ...)
Arguments
x |
[ |
... |
Ignored. |
Value
A named list suitable for passing to rstan::sampling().
Extract the stanfit object from the bscmfit object
Description
This function returns the output (stanfit object) from rstan::sampling().
Usage
get_stanfit(x, ...)
## Default S3 method:
get_stanfit(x, ...)
## S3 method for class 'bscmfit'
get_stanfit(x, ...)
Arguments
x |
[ |
... |
Ignored. |
Value
Object of class stanfit.
Leave-out donor sensitivity of a Bayesian synthetic control model
Description
leave_donor_out() re-estimates the original model after omitting donors
from the donor pool. By default, it removes one donor at a time. Setting
cumulative = TRUE instead removes donors cumulatively following the chosen
order, which can be used to assess how sensitive the results are to the most
influential donors.
Usage
leave_donor_out(x, ...)
## S3 method for class 'bscmfit'
leave_donor_out(
x,
order = "descending",
cumulative = FALSE,
probs = c(0.025, 0.5, 0.975),
...
)
Arguments
x |
[ |
... |
Additional arguments passed on to |
order |
[ |
cumulative |
[ |
probs |
[ |
Details
The donor order can always be supplied explicitly. For one-donor leave-out
runs (cumulative = FALSE), that order controls the order of the returned
runs. For cumulative leave-out runs (cumulative = TRUE), donors are removed
from the ordered pool one by one. If no order is supplied, one-donor
leave-out uses the original donor order, while cumulative leave-out orders
donors by posterior mean donor weight in descending order.
For models with multiple treated units, automatic donor ranking is based on the average posterior mean donor weight across treated units.
Value
An object of class bscm_ldo with data frames effect,
rmse, rmse_ratio, weights, and diagnostics, and a metadata list.
The data frames contain posterior summaries for the original fit
(step = 0) and for each leave-out run, identified by step,
n_removed, and last_removed. The metadata list contains the
omitted donors, whether donor omission was cumulative, summary
probabilities, and model metadata needed for plotting. The result can be
visualized with plot_weights() and plot_effects().
Leave-Future-Out Cross-Validation
Description
Estimates the leave-future-out (LFO) expected log predictive density (ELPD)
for bscmfit models. The LFO-CV is performed over the pre-treatment period:
at each step t (from L to T_pre_min - 1), the model is evaluated on
its ability to predict the next pre-treatment observation y[t+1] having
been fitted on y[1:t]. For models with multiple treated units, all units
must have at least L + 1 pre-treatment observations, and only the first
min(T_pre) time points are used so the number of treated units stays
constant throughout the evaluation period. In this case, a joint ELPD of
all treated units is computed at each step.
Usage
lfo(x, ...)
## S3 method for class 'bscmfit'
lfo(x, L, exact = FALSE, k_threshold = 0.7, ...)
Arguments
x |
[ |
... |
Additional arguments passed on to |
L |
[ |
exact |
[ |
k_threshold |
[ |
Value
An object of class bscm_lfo, a list with components:
-
ELPD: Total expected log predictive density. -
ELPD_SE: A crude approximation of standard error of ELPD, ignoring serial dependency of ELPD estimates. -
ELPDs: Vector of per-step ELPDs (lengthT_pre_min - L). Elementkis the ELPD for predicting observations atL + k. -
ks: Paretokvalues (lengthT_pre_min - L - 1,NULLfor exact LFO). -
refits: Time indices at which the model was re-estimated. -
L: The value ofLused. -
T_pre_min: The minimum pre-treatment period length across treated units. -
times: Vector of all unique time values from the original data. -
time: Name of the time variable. -
k_thres: The Pareto k threshold used.
References
Bürkner PC, Gabry J, and Vehtari A (2020). Approximate leave-future-out cross-validation for Bayesian time series models. Journal of Statistical Computation and Simulation, 90(14), 2499–2523, https://doi.org/10.1080/00949655.2020.1783262.
Examples
lfo_approx <- lfo(fit_single_treated, L = 10)
lfo_approx
lfo_exact <- lfo(fit_single_treated, L = 10, exact = TRUE)
lfo_exact
Posterior draws of pointwise log-likelihood
Description
Returns posterior draws of pointwise log-likelihoods of the treated units per pre-treatment time point
Usage
## S3 method for class 'bscmfit'
log_lik(object, ...)
Arguments
object |
[ |
... |
Ignored. |
Value
A matrix where rows correspond to posterior draws and columns to observed (non-missing) time periods. Output columns are ordered by unit then time: all observed pre-treatment time points for the first treated unit, then the second, and so on.
Examples
log_lik(fit_single_treated) |> head()
Approximate leave-one-out (LOO) cross-validation for Bayesian synthetic control models
Description
Approximate leave-one-out (LOO) cross-validation for Bayesian synthetic control models
Usage
## S3 method for class 'bscmfit'
loo(x, r_eff = TRUE, reloo = FALSE, k_threshold = 0.7, ...)
Arguments
x |
[ |
r_eff |
[ |
reloo |
[ |
k_threshold |
[ |
... |
Additional arguments to |
Value
An output from loo::loo().
References
Vehtari A, Gelman A, and Gabry J (2017). Practical Bayesian model evaluation using leave-one-out cross-validation and WAIC. Statistics and Computing. 27(5), 1413–1432, https://doi.org/10.1007/s11222-016-9696-4.
Examples
loo(fit_single_treated)
Bayesian R-squared value
Description
bayes_R2 computes the Bayesian R^2 measure of model fit for a
Bayesian synthetic control model, while loo_R2 computes a leave-one-out
adjusted version of the same quantity.
Usage
## S3 method for class 'bscmfit'
loo_R2(object, summary = TRUE, probs = c(0.025, 0.975), fixed_seed = TRUE, ...)
## S3 method for class 'bscmfit'
bayes_R2(object, summary = TRUE, probs = c(0.025, 0.975), ...)
Arguments
object |
[ |
summary |
[ |
probs |
[ |
fixed_seed |
[logical(1)] |
... |
Ignored. |
Value
data.frame of posterior summary of (LOO-adjusted) R-squared values.
References
Gelman A, Goodrich B, Gabry J, and Vehtari A (2019). R-squared for Bayesian regression models. The American Statistician. 73(3), 307–309, https://doi.org/10.1080/00031305.2018.1549100.
Examples
bayes_R2(fit_single_treated)
Simulated example data with multiple treated units
Description
This simulated data is generated based on a latent factor model with two
factors, 40 time points and 53 units. For the first three units,
a "treatment" \tau = 2 for the last 5 time points
t=0,\ldots,4 is added to an outcome y. There are two time-varying
covariates x and z, with latter one having time-varying effect on y.
To see exactly how the data was generated, see data-raw folder on the
Github repository of the package.
Format
A data frame with 1,550 rows and 5 columns:
time: Time index from ranging from -30 to 9.
id: Unit index variable 'ranging from 1 to 33.
y: Outcome variable.
x: Time-varying predictor.
z: Time-varying predictor.
treatment: Binary indicator variable where 1 corresponds to the treatment.
Examples
head(multiple_treated)
multiple_treated |>
dplyr::mutate(treated = any(treatment), .by = id) |>
ggplot2::ggplot(ggplot2::aes(time, y)) +
ggplot2::geom_line(
ggplot2::aes(group = id, colour = treated, alpha = treated)
) +
ggplot2::scale_alpha_manual(values = c(0.3, 1)) +
ggplot2::theme_bw()
Return the number of chains of bscmfit object
Description
Return the number of chains of bscmfit object
Usage
## S3 method for class 'bscmfit'
nchains(x)
Arguments
x |
[ |
Value
Number of Markov chains used in sampling as a single integer value.
Return the number of posterior draws of a bscmfit object
Description
Return the number of posterior draws of a bscmfit object
Usage
## S3 method for class 'bscmfit'
ndraws(x)
Arguments
x |
[ |
Value
Number of posterior draws as a single integer value.
Placebo effects of a Bayesian synthetic control model
Description
Estimates in-space and in-time placebo effects based on the model definition of previously estimated Bayesian synthetic control model.
Usage
placebo_effects(x, ...)
## S3 method for class 'bscmfit'
placebo_effects(x, type, L = NULL, probs = c(0.025, 0.975), ...)
Arguments
x |
[ |
... |
Additional arguments passed on to |
type |
[ |
L |
[ |
probs |
[ |
Details
For the in-space placebo (type = "donor"), original model is re-estimated
using each donor as the treated unit in turn (omitting the original, true
treated unit). The obtained effect estimates should be around zero, under
the assumption that the treatment affected only the true treated unit.
For the in-time placebo (type = "time") we still
estimate the treatment effect for the original treated, but move the start
of the treatment from L + 1 to T_{pre} + 1, where L is the
minimum number of pre-treatment time points to use, and T_{pre} is the
true last pre-treatment time point. In all cases, the obtained treatment
effects should fluctuate around zero for time points before the true
treatment time, under the assumption of no anticipation effects.
In both cases, the output includes also the results of the input model.
Value
A list with data frames effect, rmse, rmse_ratio and
diagnostics, and a metadata list. The data frames contain posterior
summaries for each run, identified by a placebo column. For
type = "donor", placebo is the treated unit name for the original fit
and the donor name for each placebo run. For type = "time", placebo
is the assumed treatment start time for each run. The metadata list
contains the placebo type, summary probabilities, model setup, and the
placebo labels used. The result can be visualized with plot_effects().
Plot Pareto k diagnostics from LFO cross-validation
Description
Plots Pareto k values against the corresponding time variable values, with a horizontal line at the refitting threshold. Only available (and relevant) for PSIS-LFO (not exact LFO).
Usage
## S3 method for class 'bscm_lfo'
plot(x, ...)
Arguments
x |
[ |
... |
Ignored. |
Value
A ggplot object, or NULL invisibly for exact LFO.
Visualize BSCM estimates
Description
A plot of the posterior mean and posterior interval of the treatment effect and synthetic control over time for single treated unit models. For models with multiple treated units, plots the average treatment effect over time since treatment across all treated units.
Usage
## S3 method for class 'bscmfit'
plot(x, probs = c(0.025, 0.975), ...)
Arguments
x |
[ |
probs |
[ |
... |
Ignored |
Value
A ggplot object
Examples
plot(fit_single_treated, probs = c(0.05, 0.95))
Visualize regression coefficients of Bayesian synthetic control model
Description
Plots posterior means and posterior intervals for the model coefficients.
For time-constant coefficients (type = "fixed"), a horizontal point-range
plot is produced. For time-varying coefficients (type = "varying"), a
ribbon-and-line plot over time is produced, faceted by parameter.
Usage
plot_coefs(x, ...)
## S3 method for class 'bscmfit'
plot_coefs(
x,
type = c("fixed", "varying"),
probs = c(0.025, 0.975),
combine = FALSE,
alpha = 0.5,
scales = "free_y",
...
)
Arguments
x |
[ |
... |
Ignored. |
type |
[ |
probs |
[ |
combine |
[ |
alpha |
[ |
scales |
[ |
Value
A ggplot object.
Examples
plot_coefs(fit_single_treated)
Visualize BSCM treatment effects
Description
plot_effects() plots the posterior mean and posterior interval of the
treatment effect over time. For the output of
leave_donor_out() or placebo_effects(), it additionally overlays the
posterior mean of the treatment effect from each leave-out or placebo fit.
Usage
plot_effects(x, ...)
## S3 method for class 'bscmfit'
plot_effects(x, probs = c(0.025, 0.975), unit = NULL, ...)
## S3 method for class 'bscm_ldo'
plot_effects(x, probs = NULL, ...)
## S3 method for class 'bscm_placebo_effects'
plot_effects(x, probs = NULL, ...)
Arguments
x |
[ |
... |
Ignored. |
probs |
[ |
unit |
[ |
Details
For models with multiple treated units, plot_effects.bscmfit() defaults
to plotting the average treatment effect over time since treatment. Supply a
unit identifier to plot the effect for a specific treated unit instead.
Value
A ggplot object.
Examples
plot_effects(fit_single_treated)
Visualize BSCM residuals
Description
plot_residuals() plots the posterior mean and posterior interval of the
pre-treatment residuals over time (type = "time"), or the posterior
distribution of the residual autocorrelations up to a given number of lags
(type = "autocorrelation"), again only for pre-treatment time points.
Usage
plot_residuals(x, ...)
## S3 method for class 'bscmfit'
plot_residuals(x, type = "time", probs = c(0.025, 0.975), max_lag = 5L, ...)
Arguments
x |
[ |
... |
Ignored. |
type |
[ |
probs |
[ |
max_lag |
[ |
Details
For models with multiple treated units, returns a named list of per-unit plots.
Value
A ggplot object, or a named list of ggplot objects when x is
a bscmfit with multiple treated units.
Examples
plot_residuals(fit_single_treated)
plot_residuals(fit_single_treated, type = "autocorrelation", max_lag = 10L)
Visualize donor weights
Description
plot_weights() visualizes posterior summaries of donor weights for a
fitted Bayesian synthetic control model. When applied to the output of
leave_donor_out(), it instead visualizes how donor
weights or their ranks change across the leave-out runs.
Usage
plot_weights(x, ...)
## S3 method for class 'bscmfit'
plot_weights(
x,
point_estimate = "median",
order = NULL,
coverage = c(0.5, 0.95),
linewidth = 1,
point_size = 2,
reverse = NULL,
...
)
## S3 method for class 'bscm_ldo'
plot_weights(
x,
type = "weight",
point_estimate = "median",
coverage = 0.95,
linewidth = 1,
point_size = 2,
reverse = TRUE,
...
)
Arguments
x |
[ |
... |
Ignored. |
point_estimate |
[´character(1) |
order |
[ |
coverage |
[ |
linewidth |
[ |
point_size |
[ |
reverse |
[ |
type |
[ |
Value
A ggplot object, or a named list of ggplot objects for models
with multiple treated units.
Examples
plot_weights(fit_single_treated)
Posterior Draws of the Expected Predictive Distribution
Description
Returns draws from the posterior distribution of the expected value of the Bayesian synthetic control model.
Usage
## S3 method for class 'bscmfit'
posterior_epred(object, ...)
Arguments
object |
[ |
... |
Ignored. |
Value
A matrix where rows correspond to posterior draws and columns to time periods.
Posterior Draws of the Linear Predictor
Description
Returns draws from the posterior distribution of the linear predictor of
the Bayesian synthetic control model. Since the model uses a Gaussian
likelihood with identity link, this is equivalent to
posterior_epred.bscmfit().
Usage
## S3 method for class 'bscmfit'
posterior_linpred(object, transform = FALSE, ...)
Arguments
object |
[ |
transform |
[ |
... |
Ignored. |
Value
A matrix where rows correspond to posterior draws and columns to time periods.
Posterior Predictive Distribution of Bayesian Synthetic Control Model
Description
Returns draws from the posterior predictive distribution of the Bayesian synthetic control. Note that this function does not simulate new realizations from this distribution, but rather returns the posterior draws computed during the model estimation.
Usage
## S3 method for class 'bscmfit'
posterior_predict(object, ...)
Arguments
object |
[ |
... |
Ignored. |
Value
A matrix of posterior predictive draws of the synthetic control, where rows correspond to posterior draws and columns to time periods.
Print method for LFO cross-validation output
Description
Prints the summary of the leave-future-out cross-validation.
Usage
## S3 method for class 'bscm_lfo'
print(x, ...)
Arguments
x |
[ |
... |
Ignored. |
Value
Returns x invisibly.
Print method for bscmfit objects
Description
Print method for bscmfit objects
Usage
## S3 method for class 'bscmfit'
print(x, ...)
## S3 method for class 'summary_bscmfit'
print(x, ...)
Arguments
x |
[ |
... |
Ignored. |
Value
Returns x (invisibly).
Returns x (invisibly).
Print MCMC Diagnostics
Description
Print MCMC Diagnostics
Usage
## S3 method for class 'bscmfit_diagnostics'
print(x, print_table = NULL, ...)
Arguments
x |
[ |
print_table |
[ |
... |
Ignored. |
Value
input x (invisibility).
Predictions from a projected BSCM model
Description
Convenience wrapper around projpred::proj_predict() that accepts data in
the original long format used by bscm(). If newdata is supplied in long
format, it is converted internally to the wide format for projpred.
Usage
proj_predict_bscm(object, newdata = NULL, ...)
Arguments
object |
A projection object or variable selection object returned by
|
newdata |
Optional new data used for predictions. The outcome, unit,
and time variables should have should have names matching the original
data. If |
... |
Additional arguments passed to |
Value
The output of projpred::proj_predict().
Autocorrelation function for residuals of BSCM
Description
Autocorrelation function for residuals of BSCM
Usage
residual_acf(x, max_lag = 5L, probs = c(0.025, 0.975))
Arguments
x |
Posterior draws from |
max_lag |
[ |
probs |
[ |
Value
A data frame.
Examples
residual_acf(fit_single_treated, max_lag = 10)
Posterior residuals of a Bayesian synthetic control model
Description
Returns posterior draws or summaries of residuals, defined as the difference between the observed outcome and the posterior expected value.
Usage
## S3 method for class 'bscmfit'
residuals(
object,
summary = TRUE,
probs = c(0.025, 0.975),
pretreatment_only = TRUE,
...
)
Arguments
object |
[ |
summary |
[ |
probs |
[ |
pretreatment_only |
[ |
... |
Ignored. |
Value
A data.frame of posterior summaries (summary = TRUE) or
posterior draws (summary = FALSE) in long format.
Examples
residuals(fit_single_treated) |> head()
Extract root mean squared errors of a Bayesian synthetic control model
Description
Returns posterior draws (or summaries) of root mean squared errors (RMSEs) for the pre-treatment and post-treatment periods.
Usage
rmse(x, ...)
## S3 method for class 'bscmfit'
rmse(x, average = FALSE, summary = TRUE, probs = c(0.025, 0.975), ...)
Arguments
x |
[ |
... |
Ignored. |
average |
[ |
summary |
[ |
probs |
[ |
Value
A tibble of posterior summaries (summary = TRUE) or
posterior draws (summary = FALSE).
Examples
rmse(fit_single_treated, probs = c(0.01, 0.5, 0.8))
Extract the ratio of post- to pre-treatment RMSE
Description
Returns posterior draws (or summaries) of the ratio of the post-treatment RMSE to the pre-treatment RMSE.
Usage
rmse_ratio(x, ...)
## S3 method for class 'bscmfit'
rmse_ratio(x, average = FALSE, summary = TRUE, probs = c(0.025, 0.975), ...)
Arguments
x |
[ |
... |
Ignored. |
average |
[ |
summary |
[ |
probs |
[ |
Value
A tibble of posterior summaries (summary = TRUE) or
posterior draws (summary = FALSE).
Examples
rmse_ratio(fit_single_treated, probs = c(0.01, 0.5, 0.8))
Extract residual standard deviations of a Bayesian synthetic control model
Description
Extract residual standard deviations of a Bayesian synthetic control model
Usage
## S3 method for class 'bscmfit'
sigma(object, summary = TRUE, probs = c(0.025, 0.975), ...)
Arguments
object |
[ |
summary |
[ |
probs |
[ |
... |
Ignored. |
Value
A data.frame of posterior summaries (summary = TRUE) or
posterior draws (summary = FALSE) in long format.
Examples
sigma(fit_single_treated)
Simulated example data with single treated unit
Description
This simulated data is generated based on a latent factor model with three
factors, 50 time points and 51 units. For the first unit (id = 1),
a "treatment" \tau = 1 + t for the last 10 time points
t=0,\ldots,9 is added to an outcome y.
To see exactly how the data was generated, see data-raw folder on the
Github repository of the package.
Format
A data frame with 2550 rows and 5 columns:
time: Time index from ranging from -40 to 9.
id: Unit index variable 'ranging from 1 to 31.
y: Outcome variable.
x: Time-varying predictor.
treatment: Binary indicator variable where 1 corresponds to the treatment.
Examples
head(single_treated)
single_treated |>
ggplot2::ggplot(ggplot2::aes(time, y, group = id)) +
ggplot2::geom_line(alpha = 0.3) +
ggplot2::geom_line(
data = single_treated |> dplyr::filter(id == 1), colour = "tomato"
) +
ggplot2::theme_bw()
Summarise posterior draws of an estimated Bayesian synthetic control model
Description
Generates posterior summary statistics for a Bayesian synthetic control
model estimated with bscm(). Note that many (more) summaries are available
via separate methods, e.g., donor_weights().
Usage
## S3 method for class 'bscmfit'
summary(object, probs = c(0.025, 0.975), ...)
Arguments
object |
[ |
probs |
[ |
... |
Ignored. |
Value
A data frame with various posterior summaries of the model.
Synthetic control series of a Bayesian synthetic control model
Description
Returns posterior draws (or summaries) of the pre- and post-treatment
trajectories of treated units, i.e. draws y_repfrom the posterior
predictive distribution of the model.
Usage
synthetic_control(x, ...)
## S3 method for class 'bscmfit'
synthetic_control(x, summary = TRUE, probs = c(0.025, 0.975), ...)
Arguments
x |
[ |
... |
Ignored. |
summary |
[ |
probs |
[ |
Value
A data.frame of posterior summaries (summary = TRUE) or
posterior draws (summary = FALSE) in long format.
Examples
synthetic_control(fit_single_treated) |> tail()
Treatment effect estimates of a Bayesian synthetic control model
Description
Treatment effect estimates of a Bayesian synthetic control model
Usage
treatment_effect(x, ...)
## S3 method for class 'bscmfit'
treatment_effect(
x,
average = TRUE,
summary = TRUE,
probs = c(0.025, 0.975),
...
)
Arguments
x |
[ |
... |
Ignored. |
average |
[ |
summary |
[ |
probs |
[ |
Value
A tibble of posterior summaries (summary = TRUE) or
posterior draws (summary = FALSE) in long format.
Examples
treatment_effect(fit_single_treated) |> tail()
Define the time-varying coefficients for BSCM
Description
Function tv() validates its arguments and returns the necessary
variables needed for the construct the time-varying coefficients in
bscm(). You should not call tv() separately, instead it should be part
of a formula, e.g., y ~ x + tv(~ z * w, df = 20).
Usage
tv(tv_formula, df = 10, type = "rw1", noncentered = TRUE)
Arguments
tv_formula |
[ |
df |
[ |
type |
[ |
noncentered |
[ |
Value
Object of class tv_term (a list).