Two independent defects in Horseshoe::GlobalUpdate() have been fixed. Both
affected every fit using prior_type = "horseshoe" (including HorseTrees()
and CausalHorseForest()). Results from earlier versions will not reproduce
under this release.
global_hp had no effect. The scale entered the auxiliary draw as
1.0 / alpha_global * alpha_global, which by operator precedence is the
constant 1.0. The user's global_hp was silently discarded, so every fit
behaved as though global_hp = 1 regardless of what was supplied. This was
reported by a referee of the accompanying manuscript.
The error scale and the forest scale were swapped. EtaPrior::GlobalUpdate()
declares (..., sigma, omega, ...) but was being called with (omega, sigma).
Both are const double&, so the swap compiled silently. The tau draw
therefore divided by sigma where it should have divided by omega. The
effect is small where omega = 1 (HorseTrees(), ShrinkageTrees()) and
substantial where omega = 1/2 (the causal forests), where it grew with the
number of covariates.
The old default of k = 0.1 was chosen while defect (1) was present, so it was
tuned for a model in which the global scale was pinned at 1. With the defect
corrected, that value shrinks far too aggressively: in our simulations it drove
pointwise coverage of the conditional treatment effect from roughly 0.95 down to
about 0.4. Re-running existing scripts against this release without updating
k is therefore worse than not upgrading at all.
The defaults are now calibrated for the corrected sampler:
| function | old | new |
|---|---|---|
| HorseTrees() | k = 0.1 | k = 1.0 |
| CausalHorseForest() | k = 0.1 | k = 1.5 |
| ShrinkageTrees() | local_hp, global_hp required | default to 1.0 / sqrt(number_of_trees) |
| CausalShrinkageForest() | four scales required | default to 1.5 / sqrt(number_of_trees_*) |
Values of k between roughly 0.5 and 1.5 (single-forest models) or 1 and 2
(causal models) worked well across a range of simulated settings, with smaller
values shrinking more aggressively and larger values being more conservative.
The defaults sit in the middle of each range. The
causal range is higher because the treatment forest enters as b * tau(x) with
b = +/- 1/2, halving its contribution to the response; measured on the scale
of that contribution the two recommendations nearly coincide. k remains
exposed on all four functions and is a natural target for cross-validation.
The local_hp / global_hp arguments are also no longer mandatory for
prior_type = "horseshoe"; either may be supplied alone, with the other taking
its default. Only their product identifies the prior, so setting them equal
loses nothing. prior_type = "horseshoe_fw" is unchanged and still requires
both to be given explicitly.
left_time and right_time columns and a
duplicated year_of_diagnosis.1; the frame is now 357 x 1004. Old fits leaked
the outcome, and positional column selection now picks different covariates.ovarian_truth: the prognostic surface, treatment effect, propensity
score and uncensored event time behind the semi-synthetic cohort. Its times
are in months; ovarian$OS_time is in days.year_of_diagnosis was documented as an instrumental variable; it is a
confounder. Wording only.predict() on a ShrinkageTrees object handed the stored training design
matrix straight to the C++ sampler. The C++ layer reads designs as a flat
row-major buffer, but R coerces a matrix to a vector column-major, so
the training data was transposed inside the sampler on every call whenever
n != p. The test design on the same call was flattened correctly, so the two
disagreed. It ran and returned plausible numbers.
predict() on a CausalShrinkageForest object was unaffected: it already
flattened both designs correctly.
Lower-dimensional posterior summarisation following Woody, Carvalho and Murray (2021). Every posterior draw of the fitted function, or of the treatment-effect and prognostic surfaces for causal fits, is projected onto a simpler summary model: a linear model (unpenalised, ridge, lasso, or elastic net), a spline-additive model, or a shallow CART tree. Returns a posterior over summaries with credible intervals, plus the posterior of the summary R-squared.
The penalised linear families are the exception: they project the posterior mean of the surface once and report point estimates, with no intervals. The spread of a penalised solution across draws has no coverage guarantee for the projection of the true surface, and the implied inclusion frequencies read as selection probabilities while being driven by the penalty.
glmnet and rpart are new suggested packages.
right_time onto left_time because the C++ layer needs a finite boundary,
but that collapse happened before the mean and scale were estimated. Those
come from survival::Surv(left_time, right_time, type = "interval2"), under
which left == right denotes an exact event, so every right-censored subject
was treated as observed. The centring mean and the sigma prior were biased
downward, by more the heavier the censoring. All four interval-censored paths
were affected: HorseTrees(), ShrinkageTrees(), CausalHorseForest() and
CausalShrinkageForest(). SurvivalBCF() and SurvivalShrinkageBCF() were
not. Interval-censored fits containing right-censored observations will
change.predict() silently dropped the Dirichlet splitting rule and re-ran under the
standard prior, affecting every SurvivalDART() and SurvivalShrinkageBCF()
fit. Predictions from those models will change.HorseTrees() base (0.95, not 1.05) and q
(0.90, not 1.00); CausalShrinkageForest() listed three of seven
prior_type values; plot(type = "path") shows the cross-validation curve,
not the summary R-squared.stats, utils and splines in Imports.There is no setting that reproduces pre-2.1.0 results exactly, because the
sampler itself has changed. The closest correspondence is that the old code was
equivalent to the corrected code with global_hp = 1, since the global scale
was inert.
summary(), plot(type = "ate"), and predict() for
CausalShrinkageForest (and CausalHorseForest) now default to a
Bayesian-bootstrap posterior for the average treatment effect: at each
MCMC iteration the per-observation CATEs are reweighted with
Dirichlet(1, ..., 1) weights before being summed, giving a draw from the
posterior of the population ATE (PATE). Credible intervals are
correspondingly wider than before because they now propagate uncertainty
in the covariate distribution, not only in tau(x).
bayesian_bootstrap = TRUE on summary.CausalShrinkageForest(),
plot.CausalShrinkageForest() (for type = "ate"), and
predict.CausalShrinkageForest(). Set to FALSE to recover the
previous equal-weight mixed ATE (MATE).predict() now also returns an ate summary for the new data and
retains the posterior CATE sample matrix as cate_samples.bayesian_bootstrap_ate() returns PATE and MATE
summaries (means, CIs, and full posterior draws) from either a fitted
CausalShrinkageForest or a CausalShrinkageForestPrediction.This is a breaking change for printed/plotted numerics: existing
scripts will report wider CIs than before. Use
bayesian_bootstrap = FALSE to reproduce previous output.
The ovarian dataset is now a single data frame (previously a list with
$clinical and $X elements). Clinical columns (OS_time, OS_event,
treatment, age, figo_stage, tumor_grade) and the 2000 gene
expression columns are combined into one data frame with 2006 columns.
This simplifies data access and aligns the format with the pdac dataset.
Code that previously used ovarian$clinical or ovarian$X must be
updated — see ?ovarian for the new structure.
plot(fit, type = "ate") and plot(fit, type = "cate") for
CausalShrinkageForest models incorrectly subtracting the control forest
predictions from the treatment forest predictions.Added the ovarian dataset: a processed TCGA-OV cohort (n = 357) for
high-dimensional survival prediction and causal inference.
See ?ovarian and examples/test-ovarian.R for a full worked example
covering survival prediction (SurvivalBART, SurvivalDART, HorseTrees) and
causal inference (SurvivalBCF, SurvivalShrinkageBCF, CausalHorseForest).
All causal model functions — CausalHorseForest(), CausalShrinkageForest(),
SurvivalBCF(), and SurvivalShrinkageBCF() — now accept a treatment_coding
argument controlling how the treatment indicator enters the BCF decomposition
y = f(x) + b * tau(x) + epsilon. Four options are available:
"centered" (default): b_i in {-1/2, 1/2}. This is the original behaviour."binary": b_i in {0, 1}. Standard binary coding."adaptive": b_i = z_i - e_hat(x_i), where e_hat(x_i) is the estimated
propensity score. This follows Hahn, Murray & Carvalho (2020) and is
implemented in the bcf package. Requires a propensity vector."invariant": Parameter-expanded (invariant) treatment coding. The coding
parameters b_0 and b_1 are assigned N(0, 1/2) priors and estimated within
the Gibbs sampler via conjugate normal updates, yielding a parameterisation
that is invariant to the coding of the treatment indicator (Hahn et al.,
2020, Section 5.2). The treatment effect is tau(x) = (b_1 - b_0) * tau_tilde(x).
Posterior draws of b_0 and b_1 are returned in the fitted object.The predict() method for CausalShrinkageForest objects automatically
carries forward the treatment coding used at training time. A propensity_test
argument is available for supplying test-set propensity scores (defaults to 0.5).
All survival-capable functions now support interval-censored data in
addition to right-censored data. Supply left_time and right_time
vectors (with outcome_type = "interval-censored") instead of y and
status. Three censoring types are distinguished:
left_time == right_time.left_time < right_time.right_time = Inf.This convention follows survival::Surv(type = "interval2"). Censored
event times are augmented within the
AFT Gibbs sampler. The following functions are affected:
HorseTrees(), ShrinkageTrees() (single-forest models)CausalHorseForest(), CausalShrinkageForest() (causal models)SurvivalBART(), SurvivalDART(), SurvivalBCF(),
SurvivalShrinkageBCF() (survival wrappers)All four primary model-fitting functions — ShrinkageTrees, HorseTrees,
CausalHorseForest, and CausalShrinkageForest — now accept an n_chains
argument (default 1). When n_chains > 1:
parallel::mclapply
(falls back to sequential execution on Windows).min(n_chains, parallel::detectCores()).N_post * n_chains total draws.ShrinkageTrees or CausalShrinkageForest
instance, so all existing print, summary, and predict methods work
without modification.SurvivalBART, SurvivalDART, SurvivalBCF, and
SurvivalShrinkageBCF inherit n_chains support through ....print and summary output adapts automatically: single-chain models show
Posterior draws, multi-chain models show Chains and Draws per chain,
with per-chain acceptance ratios listed separately.ShrinkageTrees and CausalShrinkageForest with constructors in constructors.R.print methods for both classes, displaying model specification, MCMC settings, acceptance ratio, and posterior mean sigma.summary methods for both classes, returning an inspectable object with posterior sigma (mean, SD, 95% CI), prediction summaries, variable importance (posterior inclusion probabilities), and — for causal models — ATE with credible interval (when store_posterior_sample = TRUE) and CATE heterogeneity.predict method for ShrinkageTrees, enabling posterior predictive inference on new data by re-running the sampler with stored training data and hyperparameters. Returns a ShrinkageTreesPrediction object with posterior mean and credible interval bounds. For survival models, the prediction object additionally stores predictions_sample (full posterior draws on the original scale) and sigma (posterior draws on the log-time scale), enabling posterior predictive survival curve plotting.predict method for CausalShrinkageForest, returning a CausalShrinkageForestPrediction object with three components: prognostic ($\mu(X)$), cate ($\tau(X)$), and total ($\mu(X) + \tau(X)$), each with posterior mean and credible interval bounds. For survival models with timescale = "time", predictions are back-transformed to the original time scale and the CATE becomes a multiplicative time ratio.print and summary methods for ShrinkageTreesPrediction and CausalShrinkageForestPrediction.as.mcmc.list() S3 method for ShrinkageTrees objects, converting
the sigma posterior (split by chain) into a coda::mcmc.list. This enables
all standard coda diagnostics: Gelman–Rubin R-hat, effective sample size,
Geweke test, Heidelberger–Welch test, autocorrelation plots, and more.summary() now automatically reports effective sample size (ESS) and —
for multi-chain fits — the Gelman–Rubin R-hat when the suggested package
coda is installed.coda to Suggests in DESCRIPTION.S3 plot() methods added for ShrinkageTrees, CausalShrinkageForest,
and ShrinkageTreesPrediction. Requires the suggested package ggplot2.
plot(fit, type = "trace") — sigma traceplot; one line per chain, useful for assessing mixing.plot(fit, type = "density") — overlaid posterior density of sigma, one curve per chain.plot(fit, type = "vi") — posterior credible intervals for variable inclusion probabilities (top n_vi predictors).plot(fit, type = "ate") — posterior density of the ATE with 95 % credible region (causal models only; requires store_posterior_sample = TRUE).plot(fit, type = "cate") — point estimates and 95 % credible intervals for the CATE of each training observation, sorted by posterior mean (causal models only; requires store_posterior_sample = TRUE).plot(fit, type = "vi", forest = "both") — side-by-side VI for the control and treatment forests (causal models only).plot(fit, type = "survival") — posterior survival curves
$S(t | x_i) = 1 - \Phi((\log t - \mu_i) / \sigma)$ derived from the AFT
log-normal model (survival outcomes only).
obs = NULL): computes
$\bar{S}(t) = n^{-1} \sum_i S(t | x_i)$ at each MCMC iteration with
pointwise credible bands.obs = c(1, 5, ...)): one curve per selected
training observation with its own credible band.level controls the credible band width (default 0.95).t_grid allows a custom time grid; auto-generated if NULL.km = TRUE overlays the Kaplan–Meier estimate as a dashed black
step function (population-averaged plot only; requires survival
package). Ignored with a message when obs is not NULL.plot(pred, type = "survival") — posterior predictive survival curves
for new (out-of-sample) data from predict(). Same obs, t_grid, and
level arguments as above. The KM overlay is not available for prediction
objects.HorseTrees, ShrinkageTrees,
SurvivalBART, SurvivalDART, SurvivalBCF, SurvivalShrinkageBCF,
CausalHorseForest, CausalShrinkageForest), all S3 methods (print,
summary, predict, plot), multi-chain MCMC, and a full TCGA PAAD
case study.SurvivalBART(), SurvivalDART(), SurvivalBCF(), and
SurvivalShrinkageBCF() now accept store_posterior_sample as an
explicit parameter (default TRUE), avoiding a "matched by multiple
actual arguments" error when passing it via ....CausalHorseForest and CausalShrinkageForest failing with
"argument 'y_train' is missing" when called directly (broken constructor
call introduced in 2.0.0 S3 refactor).plot(..., type = "vi") crashing with "argument must be coercible to non-negative integer" in all four model functions: covariate matrices
were being stored as flat numeric vectors instead of matrices, making
ncol() return NULL.HorseTrees and ShrinkageTrees where sigma_hat,
y_mean, and lambda were not initialised in the binary (probit) branch.SurvivalBCF wrapper for AFT-based Bayesian Causal Forests.SurvivalShrinkageBCF with Dirichlet structural sparsity.SurvivalDART for sparse high-dimensional survival modeling.SurvivalBART for survival modeling using standard BART.....std::map structure with a more efficient vector-based lookup, improving overall computational speed by approximately 30%.standard prior type option for ShrinkageTrees, corresponding to the conventional BART implementation without reversible jumps.CONTRIBUTING.md file with guidelines for contributors🎉 First CRAN release of ShrinkageTrees!
This package provides Bayesian regression tree models with shrinkage priors, supporting:
It includes four core functions:
HorseTrees(): fits a single regression tree with a standard Horseshoe prior.ShrinkageTrees(): fits a single tree with customizable shrinkage priors.CausalHorseForest(): fits a causal forest using the standard Horseshoe prior.CausalShrinkageForest(): fits a flexible causal forest with user-defined shrinkage priors and tuning options.The ...Trees functions use a single learner to estimate the outcome model directly. In contrast, the Causal...Forest variants fit separate models for the treated and control regression function. This enables estimation of conditional average treatment effects (CATEs).