Package 'stochtree'

Description

Flexible stochastic tree ensemble software. Robust implementations of Bayesian Additive Regression Trees (BART) (Chipman, George, McCulloch (2010) doi:10.1214/09-AOAS285) for supervised learning and Bayesian Causal Forests (BCF) (Hahn, Murray, Carvalho (2020) doi:10.1214/19-BA1195) for causal inference. Enables model serialization and parallel sampling and provides a low-level interface for custom stochastic forest samplers. Includes the grow-from-root algorithm for accelerated forest sampling (He and Hahn (2021) doi:10.1080/01621459.2021.1942012), a log-linear leaf model for forest-based heteroskedasticity (Murray (2020) doi:10.1080/01621459.2020.1813587), and the cloglog BART model of Alam and Linero (2025) doi:10.48550/arXiv.2502.00606 for ordinal outcomes.

Author(s)

Maintainer: Drew Herren [email protected] (ORCID)

Authors:

• Richard Hahn

• Jared Murray

• Carlos Carvalho

• Jingyu He

Other contributors:

• Pedro Lima [contributor]

• Entejar Alam [contributor]

• stochtree contributors [copyright holder]

• Eigen contributors (C++ source uses the Eigen library for matrix operations, see inst/COPYRIGHTS) [copyright holder]

• xgboost contributors (C++ tree code and related operations include or are inspired by code from the xgboost library, see inst/COPYRIGHTS) [copyright holder]

• treelite contributors (C++ tree code and related operations include or are inspired by code from the treelite library, see inst/COPYRIGHTS) [copyright holder]

• Microsoft Corporation (C++ I/O and various project structure code include or are inspired by code from the LightGBM library, which is a copyright of Microsoft, see inst/COPYRIGHTS) [copyright holder]

• Niels Lohmann (C++ source uses the JSON for Modern C++ library for JSON operations, see inst/COPYRIGHTS) [copyright holder]

• Daniel Lemire (C++ source uses the fast_double_parser library internally, see inst/COPYRIGHTS) [copyright holder]

• Victor Zverovich (C++ source uses the fmt library internally, see inst/COPYRIGHTS) [copyright holder]

See Also

Useful links:

• https://stochtree.ai/

• https://github.com/StochasticTree/stochtree

• Report bugs at https://github.com/StochasticTree/stochtree/issues

Description

Run the BART algorithm for supervised learning.

Usage

bart(
  X_train,
  y_train,
  leaf_basis_train = NULL,
  rfx_group_ids_train = NULL,
  rfx_basis_train = NULL,
  X_test = NULL,
  leaf_basis_test = NULL,
  rfx_group_ids_test = NULL,
  rfx_basis_test = NULL,
  observation_weights = NULL,
  num_gfr = 5,
  num_burnin = 0,
  num_mcmc = 100,
  previous_model_json = NULL,
  previous_model_warmstart_sample_num = NULL,
  general_params = list(),
  mean_forest_params = list(),
  variance_forest_params = list(),
  random_effects_params = list()
)

Arguments

Value

List of sampling outputs and a wrapper around the sampled forests (which can be used for in-memory prediction on new data, or serialized to JSON on disk).

Examples

n <- 100
p <- 5
X <- matrix(runif(n*p), ncol = p)
f_XW <- (
    ((0 <= X[,1]) & (0.25 > X[,1])) * (-7.5) +
    ((0.25 <= X[,1]) & (0.5 > X[,1])) * (-2.5) +
    ((0.5 <= X[,1]) & (0.75 > X[,1])) * (2.5) +
    ((0.75 <= X[,1]) & (1 > X[,1])) * (7.5)
)
noise_sd <- 1
y <- f_XW + rnorm(n, 0, noise_sd)
test_set_pct <- 0.2
n_test <- round(test_set_pct*n)
n_train <- n - n_test
test_inds <- sort(sample(1:n, n_test, replace = FALSE))
train_inds <- (1:n)[!((1:n) %in% test_inds)]
X_test <- X[test_inds,]
X_train <- X[train_inds,]
y_test <- y[test_inds]
y_train <- y[train_inds]

bart_model <- bart(X_train = X_train, y_train = y_train, X_test = X_test,
                   num_gfr = 10, num_burnin = 0, num_mcmc = 10)

Description

BART models contains external pointers to C++ objects, which means they cannot be correctly serialized to .Rds from an R session in their default state. This group of serialization functions allow us to convert between C++ data structures and a persistent JSON representation. The CppJson class wraps a performant C++ JSON API, and the functions saveBARTModelToJson and createBARTModelFromJson save to and load from this format. This representation, of course, also relies on external C++ pointers, so in order to save and reload BART models across sessions, we provide two other interfaces.

saveBARTModelToJsonString converts a BART model to an in-memory string containing the model's JSON representation and createBARTModelFromJsonString converts this representation back to a BART model object.

saveBARTModelToJsonFile and createBARTModelFromJsonFile save or reload a BART model directly to / from a .json file.

Finally, for cases in which multiple BART models have been sampled (for instance, multiple processes run via doParallel), we offer createBARTModelFromCombinedJson and createBARTModelFromCombinedJsonString for loading a new combined BART model from a list of BART JSON objects / strings.

Usage

saveBARTModelToJson(object)

saveBARTModelToJsonFile(object, filename)

saveBARTModelToJsonString(object)

createBARTModelFromJson(json_object)

createBARTModelFromJsonFile(json_filename)

createBARTModelFromJsonString(json_string)

createBARTModelFromCombinedJson(json_object_list)

createBARTModelFromCombinedJsonString(json_string_list)

Arguments

Value

saveBARTModelToJson return an object of type CppJson. saveBARTModelToJsonString returns a string dump of the BART model's JSON representation. saveBARTModelToJsonFile returns nothing, but writes to the provided filename.

createBARTModelFromJson, createBARTModelFromJsonFile, createBARTModelFromJsonString, createBARTModelFromCombinedJson, and createBARTModelFromCombinedJsonString all return objects of type bartmodel.

Examples

# Generate data
n <- 100
p <- 5
X <- matrix(runif(n*p), ncol = p)
y <- X[,1] + rnorm(n, 0, 1)

# Sample BART model
bart_model <- bart(X_train = X, y_train = y, num_gfr = 0,
                   num_burnin = 0, num_mcmc = 10)

# Save to in-memory JSON
bart_json <- saveBARTModelToJson(bart_model)
# Save to JSON string
bart_json_string <- saveBARTModelToJsonString(bart_model)
# Save to JSON file
tmpjson <- tempfile(fileext = ".json")
saveBARTModelToJsonFile(bart_model, file.path(tmpjson))

# Reload BART model from in-memory JSON object
bart_model_roundtrip <- createBARTModelFromJson(bart_json)
# Reload BART model from JSON string
bart_model_roundtrip <- createBARTModelFromJsonString(bart_json_string)
# Reload BART model from JSON file
bart_model_roundtrip <- createBARTModelFromJsonFile(file.path(tmpjson))
unlink(tmpjson)
# Reload BART model from list of JSON objects
bart_model_roundtrip <- createBARTModelFromCombinedJson(list(bart_json))
# Reload BART model from list of JSON strings
bart_model_roundtrip <- createBARTModelFromCombinedJsonString(list(bart_json_string))

Description

Run the Bayesian Causal Forest (BCF) algorithm for regularized causal effect estimation.

Usage

bcf(
  X_train,
  Z_train,
  y_train,
  propensity_train = NULL,
  rfx_group_ids_train = NULL,
  rfx_basis_train = NULL,
  X_test = NULL,
  Z_test = NULL,
  propensity_test = NULL,
  rfx_group_ids_test = NULL,
  rfx_basis_test = NULL,
  observation_weights = NULL,
  num_gfr = 5,
  num_burnin = 0,
  num_mcmc = 100,
  previous_model_json = NULL,
  previous_model_warmstart_sample_num = NULL,
  general_params = list(),
  prognostic_forest_params = list(),
  treatment_effect_forest_params = list(),
  variance_forest_params = list(),
  random_effects_params = list()
)

Arguments

Value

List of sampling outputs and a wrapper around the sampled forests (which can be used for in-memory prediction on new data, or serialized to JSON on disk).

Examples

n <- 500
p <- 5
X <- matrix(runif(n*p), ncol = p)
mu_x <- (
    ((0 <= X[,1]) & (0.25 > X[,1])) * (-7.5) +
    ((0.25 <= X[,1]) & (0.5 > X[,1])) * (-2.5) +
    ((0.5 <= X[,1]) & (0.75 > X[,1])) * (2.5) +
    ((0.75 <= X[,1]) & (1 > X[,1])) * (7.5)
)
pi_x <- (
    ((0 <= X[,1]) & (0.25 > X[,1])) * (0.2) +
    ((0.25 <= X[,1]) & (0.5 > X[,1])) * (0.4) +
    ((0.5 <= X[,1]) & (0.75 > X[,1])) * (0.6) +
    ((0.75 <= X[,1]) & (1 > X[,1])) * (0.8)
)
tau_x <- (
    ((0 <= X[,2]) & (0.25 > X[,2])) * (0.5) +
    ((0.25 <= X[,2]) & (0.5 > X[,2])) * (1.0) +
    ((0.5 <= X[,2]) & (0.75 > X[,2])) * (1.5) +
    ((0.75 <= X[,2]) & (1 > X[,2])) * (2.0)
)
Z <- rbinom(n, 1, pi_x)
noise_sd <- 1
y <- mu_x + tau_x*Z + rnorm(n, 0, noise_sd)
test_set_pct <- 0.2
n_test <- round(test_set_pct*n)
n_train <- n - n_test
test_inds <- sort(sample(1:n, n_test, replace = FALSE))
train_inds <- (1:n)[!((1:n) %in% test_inds)]
X_test <- X[test_inds,]
X_train <- X[train_inds,]
pi_test <- pi_x[test_inds]
pi_train <- pi_x[train_inds]
Z_test <- Z[test_inds]
Z_train <- Z[train_inds]
y_test <- y[test_inds]
y_train <- y[train_inds]
mu_test <- mu_x[test_inds]
mu_train <- mu_x[train_inds]
tau_test <- tau_x[test_inds]
tau_train <- tau_x[train_inds]
bcf_model <- bcf(X_train = X_train, Z_train = Z_train, y_train = y_train,
                 propensity_train = pi_train, X_test = X_test, Z_test = Z_test,
                 propensity_test = pi_test, num_gfr = 10,
                 num_burnin = 0, num_mcmc = 10)

Description

BCF models contains external pointers to C++ objects, which means they cannot be correctly serialized to .Rds from an R session in their default state. These functions allow us to convert between C++ data structures and a persistent JSON representation. The CppJson class wraps a performant C++ JSON API, and the functions saveBCFModelToJson and createBCFModelFromJson save to and load from this format. This representation, of course, also relies on external C++ pointers, so in order to save and reload BCF models across sessions, we provide two other interfaces.

saveBCFModelToJsonFile and createBCFModelFromJsonFile save or reload a BCF model's JSON representation directly to / from a .json file.

saveBCFModelToJsonString and createBCFModelFromJsonString handle in-memory strings containing JSON data, which can be written to disk or passed between processes.

Finally, for cases in which multiple BCF models have been sampled (for instance, sampled in multiple processes via doParallel), we offer createBCFModelFromCombinedJson and createBCFModelFromCombinedJsonString for loading a new combined BCF model from a list of BCF JSON objects or strings.

Usage

saveBCFModelToJson(object)

saveBCFModelToJsonFile(object, filename)

saveBCFModelToJsonString(object)

createBCFModelFromJson(json_object)

createBCFModelFromJsonFile(json_filename)

createBCFModelFromJsonString(json_string)

createBCFModelFromCombinedJson(json_object_list)

createBCFModelFromCombinedJsonString(json_string_list)

Arguments

Value

saveBCFModelToJson return an object of type CppJson. saveBCFModelToJsonFile returns nothing, but writes to the provided filename. saveBCFModelToJsonString returns a string dump of the BCF model's JSON representation.

createBCFModelFromJson, createBCFModelFromJsonFile, createBCFModelFromJsonString, createBCFModelFromCombinedJson, and createBCFModelFromCombinedJsonString all return objects of type bcfmodel.

Examples

# Generate data
n <- 100
p <- 5
X <- matrix(runif(n*p), ncol = p)
pi_X <- runif(n, 0.3, 0.7)
Z <- rbinom(n, p = pi_X, size = 1)
y <- X[,1] + Z + rnorm(n, 0, 1)

# Sample BCF model
bcf_model <- bcf(X_train = X, Z_train = Z, propensity_train = pi_X, y_train = y,
                 num_gfr = 0, num_burnin = 0, num_mcmc = 10)

# Save to in-memory JSON
bcf_json <- saveBCFModelToJson(bcf_model)
# Save to JSON string
bcf_json_string <- saveBCFModelToJsonString(bcf_model)
# Save to JSON file
tmpjson <- tempfile(fileext = ".json")
saveBCFModelToJsonFile(bcf_model, file.path(tmpjson))

# Reload BCF model from in-memory JSON object
bcf_model_roundtrip <- createBCFModelFromJson(bcf_json)
# Reload BCF model from JSON string
bcf_model_roundtrip <- createBCFModelFromJsonString(bcf_json_string)
# Reload BCF model from JSON file
bcf_model_roundtrip <- createBCFModelFromJsonFile(file.path(tmpjson))
unlink(tmpjson)
# Reload BCF model from list of JSON objects
bcf_model_roundtrip <- createBCFModelFromCombinedJson(list(bcf_json))
# Reload BCF model from list of JSON strings
bcf_model_roundtrip <- createBCFModelFromCombinedJsonString(list(bcf_json_string))

Description

Calibrate the scale parameter on an inverse gamma prior for the global error variance as in Chipman et al (2022)

This function is intended for advanced use cases in which users require detailed control of sampling algorithms and data structures. Minimal input validation and error checks are performed – users are responsible for providing the correct inputs. For tutorials on the "proper" usage of the stochtree's advanced workflow, we provide several vignettes at https://stochtree.ai/

Chipman, H., George, E., Hahn, R., McCulloch, R., Pratola, M. and Sparapani, R. (2022). Bayesian Additive Regression Trees, Computational Approaches. In Wiley StatsRef: Statistics Reference Online (eds N. Balakrishnan, T. Colton, B. Everitt, W. Piegorsch, F. Ruggeri and J.L. Teugels). doi:10.1002/9781118445112.stat08288

Usage

calibrateInverseGammaErrorVariance(
  y,
  X,
  W = NULL,
  nu = 3,
  quant = 0.9,
  standardize = TRUE
)

Arguments

Value

Value of lambda which determines the scale parameter of the global error variance prior (sigma^2 ~ IG(nu,nu*lambda))

Examples

n <- 100
p <- 5
X <- matrix(runif(n*p), ncol = p)
y <- 10*X[,1] - 20*X[,2] + rnorm(n)
nu <- 3
lambda <- calibrateInverseGammaErrorVariance(y, X, nu = nu)
sigma2hat <- mean(resid(lm(y~X))^2)
mean(var(y)/rgamma(100000, nu, rate = nu*lambda) < sigma2hat)

Description

Compute posterior credible intervals for specified terms from a fitted BART model. Supports intervals for mean functions, variance functions, random effects, and overall outcome predictions.

Usage

computeBARTPosteriorInterval(
  model_object,
  terms,
  level = 0.95,
  scale = "linear",
  X = NULL,
  leaf_basis = NULL,
  rfx_group_ids = NULL,
  rfx_basis = NULL
)

Arguments

Value

A list containing the lower and upper bounds of the credible interval for the specified term. If multiple terms are requested, a named list with intervals for each term is returned.

Examples

n <- 100
p <- 5
X <- matrix(rnorm(n * p), nrow = n, ncol = p)
y <- 2 * X[,1] + rnorm(n)
bart_model <- bart(y_train = y, X_train = X)
intervals <- computeBARTPosteriorInterval(
 model_object = bart_model,
 terms = c("mean_forest", "y_hat"),
 X = X,
 level = 0.90
)

Description

Compute posterior credible intervals for specified terms from a fitted BCF model. Supports intervals for prognostic forests, CATE forests, variance forests, random effects, and overall mean outcome predictions.

Usage

computeBCFPosteriorInterval(
  model_object,
  terms,
  level = 0.95,
  scale = "linear",
  X = NULL,
  Z = NULL,
  propensity = NULL,
  rfx_group_ids = NULL,
  rfx_basis = NULL
)

Arguments

Value

A list containing the lower and upper bounds of the credible interval for the specified term. If multiple terms are requested, a named list with intervals for each term is returned.

Examples

n <- 100
p <- 5
X <- matrix(rnorm(n * p), nrow = n, ncol = p)
pi_X <- pnorm(0.5 * X[,1])
Z <- rbinom(n, 1, pi_X)
mu_X <- X[,1]
tau_X <- 0.25 * X[,2]
y <- mu_X + tau_X * Z + rnorm(n)
bcf_model <- bcf(X_train = X, Z_train = Z, y_train = y,
                 propensity_train = pi_X)
intervals <- computeBCFPosteriorInterval(
 model_object = bcf_model,
 terms = c("prognostic_function", "cate"),
 X = X,
 Z = Z,
 propensity = pi_X,
 level = 0.90
)

Description

Compute a contrast using a BART model by making two sets of outcome predictions and taking their difference. This function provides the flexibility to compute any contrast of interest by specifying covariates, leaf basis, and random effects bases / IDs for both sides of a two term contrast. For simplicity, we refer to the subtrahend of the contrast as the "control" or Y0 term and the minuend of the contrast as the Y1 term, though the requested contrast need not match the "control vs treatment" terminology of a classic two-treatment causal inference problem. We mirror the function calls and terminology of the predict.bartmodel function, labeling each prediction data term with a 1 to denote its contribution to the treatment prediction of a contrast and 0 to denote inclusion in the control prediction.

Only valid when there is either a mean forest or a random effects term in the BART model.

Usage

computeContrastBARTModel(
  object,
  X_0,
  X_1,
  leaf_basis_0 = NULL,
  leaf_basis_1 = NULL,
  rfx_group_ids_0 = NULL,
  rfx_group_ids_1 = NULL,
  rfx_basis_0 = NULL,
  rfx_basis_1 = NULL,
  type = "posterior",
  scale = "linear"
)

Arguments

Value

Contrast matrix or vector, depending on whether type = "mean" or "posterior".

Examples

n <- 100
p <- 5
X <- matrix(runif(n*p), ncol = p)
W <- matrix(runif(n*1), ncol = 1)
f_XW <- (
    ((0 <= X[,1]) & (0.25 > X[,1])) * (-7.5*W[,1]) +
    ((0.25 <= X[,1]) & (0.5 > X[,1])) * (-2.5*W[,1]) +
    ((0.5 <= X[,1]) & (0.75 > X[,1])) * (2.5*W[,1]) +
    ((0.75 <= X[,1]) & (1 > X[,1])) * (7.5*W[,1])
)
noise_sd <- 1
y <- f_XW + rnorm(n, 0, noise_sd)
test_set_pct <- 0.2
n_test <- round(test_set_pct*n)
n_train <- n - n_test
test_inds <- sort(sample(1:n, n_test, replace = FALSE))
train_inds <- (1:n)[!((1:n) %in% test_inds)]
X_test <- X[test_inds,]
X_train <- X[train_inds,]
W_test <- W[test_inds,]
W_train <- W[train_inds,]
y_test <- y[test_inds]
y_train <- y[train_inds]
bart_model <- bart(X_train = X_train, leaf_basis_train = W_train, y_train = y_train,
                   num_gfr = 10, num_burnin = 0, num_mcmc = 10)
contrast_test <- computeContrastBARTModel(
    bart_model,
    X_0 = X_test,
    X_1 = X_test,
    leaf_basis_0 = matrix(0, nrow = n_test, ncol = 1),
    leaf_basis_1 = matrix(1, nrow = n_test, ncol = 1),
    type = "posterior",
    scale = "linear"
)

Description

Compute a contrast using a BCF model by making two sets of outcome predictions and taking their difference. For simple BCF models with binary treatment, this will yield the same prediction as requesting terms = "cate" in the predict.bcfmodel function. For more general models, such as models with continuous / multivariate treatments or an additive random effects term with a coefficient on the treatment, this function provides the flexibility to compute a any contrast of interest by specifying covariates, treatment, and random effects bases and IDs for both sides of a two term contrast. For simplicity, we refer to the subtrahend of the contrast as the "control" or Y0 term and the minuend of the contrast as the Y1 term, though the requested contrast need not match the "control vs treatment" terminology of a classic two-arm experiment. We mirror the function calls and terminology of the predict.bcfmodel function, labeling each prediction data term with a 1 to denote its contribution to the treatment prediction of a contrast and 0 to denote inclusion in the control prediction.

Usage

computeContrastBCFModel(
  object,
  X_0,
  X_1,
  Z_0,
  Z_1,
  propensity_0 = NULL,
  propensity_1 = NULL,
  rfx_group_ids_0 = NULL,
  rfx_group_ids_1 = NULL,
  rfx_basis_0 = NULL,
  rfx_basis_1 = NULL,
  type = "posterior",
  scale = "linear"
)

Arguments

Value

List of prediction matrices or single prediction matrix / vector, depending on the terms requested.

Examples

n <- 500
p <- 5
X <- matrix(runif(n*p), ncol = p)
mu_x <- (
    ((0 <= X[,1]) & (0.25 > X[,1])) * (-7.5) +
    ((0.25 <= X[,1]) & (0.5 > X[,1])) * (-2.5) +
    ((0.5 <= X[,1]) & (0.75 > X[,1])) * (2.5) +
    ((0.75 <= X[,1]) & (1 > X[,1])) * (7.5)
)
pi_x <- (
    ((0 <= X[,1]) & (0.25 > X[,1])) * (0.2) +
    ((0.25 <= X[,1]) & (0.5 > X[,1])) * (0.4) +
    ((0.5 <= X[,1]) & (0.75 > X[,1])) * (0.6) +
    ((0.75 <= X[,1]) & (1 > X[,1])) * (0.8)
)
tau_x <- (
    ((0 <= X[,2]) & (0.25 > X[,2])) * (0.5) +
    ((0.25 <= X[,2]) & (0.5 > X[,2])) * (1.0) +
    ((0.5 <= X[,2]) & (0.75 > X[,2])) * (1.5) +
    ((0.75 <= X[,2]) & (1 > X[,2])) * (2.0)
)
Z <- rbinom(n, 1, pi_x)
noise_sd <- 1
y <- mu_x + tau_x*Z + rnorm(n, 0, noise_sd)
test_set_pct <- 0.2
n_test <- round(test_set_pct*n)
n_train <- n - n_test
test_inds <- sort(sample(1:n, n_test, replace = FALSE))
train_inds <- (1:n)[!((1:n) %in% test_inds)]
X_test <- X[test_inds,]
X_train <- X[train_inds,]
pi_test <- pi_x[test_inds]
pi_train <- pi_x[train_inds]
Z_test <- Z[test_inds]
Z_train <- Z[train_inds]
y_test <- y[test_inds]
y_train <- y[train_inds]
mu_test <- mu_x[test_inds]
mu_train <- mu_x[train_inds]
tau_test <- tau_x[test_inds]
tau_train <- tau_x[train_inds]
bcf_model <- bcf(X_train = X_train, Z_train = Z_train, y_train = y_train,
                 propensity_train = pi_train, num_gfr = 10,
                 num_burnin = 0, num_mcmc = 10)
tau_hat_test <- computeContrastBCFModel(
    bcf_model, X_0=X_test, X_1=X_test, Z_0=rep(0, n_test), Z_1=rep(1, n_test),
    propensity_0 = pi_test, propensity_1 = pi_test
)

Description

Wrapper around a C++ nlohmann::json object

This class is intended for advanced use cases in which users require detailed control of sampling algorithms and data structures. Minimal input validation and error checks are performed – users are responsible for providing the correct inputs. For tutorials on the "proper" usage of the stochtree's advanced workflow, we provide several vignettes at https://stochtree.ai/

Public fields

Methods

Public methods

• CppJson$new()

• CppJson$add_forest()

• CppJson$add_random_effects()

• CppJson$add_scalar()

• CppJson$add_integer()

• CppJson$add_boolean()

• CppJson$add_string()

• CppJson$add_vector()

• CppJson$add_integer_vector()

• CppJson$add_string_vector()

• CppJson$add_list()

• CppJson$add_string_list()

• CppJson$get_scalar()

• CppJson$get_integer()

• CppJson$get_boolean()

• CppJson$get_string()

• CppJson$get_vector()

• CppJson$get_integer_vector()

• CppJson$get_string_vector()

• CppJson$get_numeric_list()

• CppJson$get_string_list()

• CppJson$return_json_string()

• CppJson$save_file()

• CppJson$load_from_file()

• CppJson$load_from_string()

Method new()

Create a new CppJson object.

Usage

CppJson$new()

Returns

A new CppJson object.

Method add_forest()

Convert a forest container to json and add to the current CppJson object

Usage

CppJson$add_forest(forest_samples)

Arguments

Returns

None

Method add_random_effects()

Convert a random effects container to json and add to the current CppJson object

Usage

CppJson$add_random_effects(rfx_samples)

Arguments

Returns

None

Method add_scalar()

Add a scalar to the json object under the name "field_name" (with optional subfolder "subfolder_name")

Usage

CppJson$add_scalar(field_name, field_value, subfolder_name = NULL)

Arguments

Returns

None

Method add_integer()

Add a scalar to the json object under the name "field_name" (with optional subfolder "subfolder_name")

Usage

CppJson$add_integer(field_name, field_value, subfolder_name = NULL)

Arguments

Returns

None

Method add_boolean()

Add a boolean value to the json object under the name "field_name" (with optional subfolder "subfolder_name")

Usage

CppJson$add_boolean(field_name, field_value, subfolder_name = NULL)

Arguments

Returns

None

Method add_string()

Add a string value to the json object under the name "field_name" (with optional subfolder "subfolder_name")

Usage

CppJson$add_string(field_name, field_value, subfolder_name = NULL)

Arguments

Returns

None

Method add_vector()

Add a vector to the json object under the name "field_name" (with optional subfolder "subfolder_name")

Usage

CppJson$add_vector(field_name, field_vector, subfolder_name = NULL)

Arguments

Returns

None

Method add_integer_vector()

Add an integer vector to the json object under the name "field_name" (with optional subfolder "subfolder_name")

Usage

CppJson$add_integer_vector(field_name, field_vector, subfolder_name = NULL)

Arguments

Returns

None

Method add_string_vector()

Add an array to the json object under the name "field_name" (with optional subfolder "subfolder_name")

Usage

CppJson$add_string_vector(field_name, field_vector, subfolder_name = NULL)

Arguments

Returns

None

Method add_list()

Add a list of vectors (as an object map of arrays) to the json object under the name "field_name"

Usage

CppJson$add_list(field_name, field_list)

Arguments

Returns

None

Method add_string_list()

Add a list of vectors (as an object map of arrays) to the json object under the name "field_name"

Usage

CppJson$add_string_list(field_name, field_list)

Arguments

Returns

None

Method get_scalar()

Retrieve a scalar value from the json object under the name "field_name" (with optional subfolder "subfolder_name")

Usage

CppJson$get_scalar(field_name, subfolder_name = NULL)

Arguments

Returns

None

Method get_integer()

Retrieve a integer value from the json object under the name "field_name" (with optional subfolder "subfolder_name")

Usage

CppJson$get_integer(field_name, subfolder_name = NULL)

Arguments

Returns

None

Method get_boolean()

Retrieve a boolean value from the json object under the name "field_name" (with optional subfolder "subfolder_name")

Usage

CppJson$get_boolean(field_name, subfolder_name = NULL)

Arguments

Returns

None

Method get_string()

Retrieve a string value from the json object under the name "field_name" (with optional subfolder "subfolder_name")

Usage

CppJson$get_string(field_name, subfolder_name = NULL)

Arguments

Returns

None

Method get_vector()

Retrieve a vector from the json object under the name "field_name" (with optional subfolder "subfolder_name")

Usage

CppJson$get_vector(field_name, subfolder_name = NULL)

Arguments

Returns

None

Method get_integer_vector()

Retrieve an integer vector from the json object under the name "field_name" (with optional subfolder "subfolder_name")

Usage

CppJson$get_integer_vector(field_name, subfolder_name = NULL)

Arguments

Returns

None

Method get_string_vector()

Retrieve a character vector from the json object under the name "field_name" (with optional subfolder "subfolder_name")

Usage

CppJson$get_string_vector(field_name, subfolder_name = NULL)

Arguments

Returns

None

Method get_numeric_list()

Reconstruct a list of numeric vectors from the json object stored under "field_name"

Usage

CppJson$get_numeric_list(field_name, key_names)

Arguments

Returns

None

Method get_string_list()

Reconstruct a list of string vectors from the json object stored under "field_name"

Usage

CppJson$get_string_list(field_name, key_names)

Arguments

Returns

None

Method return_json_string()

Convert a JSON object to in-memory string

Usage

CppJson$return_json_string()

Returns

JSON string

Method save_file()

Save a json object to file

Usage

CppJson$save_file(filename)

Arguments

Returns

None

Method load_from_file()

Load a json object from file

Usage

CppJson$load_from_file(filename)

Arguments

Returns

None

Method load_from_string()

Load a json object from string

Usage

CppJson$load_from_string(json_string)

Arguments

Returns

None

Description

Wrapper around a C++ random number generator object (for reproducibility). The class persists a C++ random number generator throughout an R session to ensure a given seed generates the same outputs (on the same OS). If no seed is provided, the C++ random number generator is initialized using std::random_device.

Public fields

Methods

Public methods

• CppRNG$new()

Method new()

Create a new CppRNG object.

Usage

CppRNG$new(random_seed = -1)

Arguments

Returns

A new CppRNG object.

Description

Create an R class that wraps a C++ random number generator

Usage

createCppRNG(random_seed = -1)

Arguments

Value

CppRng object

Examples

rng <- createCppRNG(1234)
rng <- createCppRNG()

Description

Create a forest

This function is intended for advanced use cases in which users require detailed control of sampling algorithms and data structures. Minimal input validation and error checks are performed – users are responsible for providing the correct inputs. For tutorials on the "proper" usage of the stochtree's advanced workflow, we provide several vignettes at https://stochtree.ai/

Usage

createForest(
  num_trees,
  leaf_dimension = 1,
  is_leaf_constant = FALSE,
  is_exponentiated = FALSE
)

Arguments

Value

Forest object

Examples

num_trees <- 100
leaf_dimension <- 2
is_leaf_constant <- FALSE
is_exponentiated <- FALSE
forest <- createForest(num_trees, leaf_dimension, is_leaf_constant, is_exponentiated)

Description

Create a forest dataset object

This function is intended for advanced use cases in which users require detailed control of sampling algorithms and data structures. Minimal input validation and error checks are performed – users are responsible for providing the correct inputs. For tutorials on the "proper" usage of the stochtree's advanced workflow, we provide several vignettes at https://stochtree.ai/

Usage

createForestDataset(covariates, basis = NULL, variance_weights = NULL)

Arguments

Value

ForestDataset object

Examples

covariate_matrix <- matrix(runif(10*100), ncol = 10)
basis_matrix <- matrix(rnorm(3*100), ncol = 3)
weight_vector <- rnorm(100)
forest_dataset <- createForestDataset(covariate_matrix)
forest_dataset <- createForestDataset(covariate_matrix, basis_matrix)
forest_dataset <- createForestDataset(covariate_matrix, basis_matrix, weight_vector)

Description

Create a forest model object

This function is intended for advanced use cases in which users require detailed control of sampling algorithms and data structures. Minimal input validation and error checks are performed – users are responsible for providing the correct inputs. For tutorials on the "proper" usage of the stochtree's advanced workflow, we provide several vignettes at https://stochtree.ai/

Usage

createForestModel(forest_dataset, forest_model_config, global_model_config)

Arguments

Value

ForestModel object

Examples

num_trees <- 100
n <- 100
p <- 10
alpha <- 0.95
beta <- 2.0
min_samples_leaf <- 2
max_depth <- 10
feature_types <- as.integer(rep(0, p))
X <- matrix(runif(n*p), ncol = p)
forest_dataset <- createForestDataset(X)
forest_model_config <- createForestModelConfig(feature_types=feature_types,
                                               num_trees=num_trees, num_features=p,
                                               num_observations=n, alpha=alpha, beta=beta,
                                               min_samples_leaf=min_samples_leaf,
                                               max_depth=max_depth, leaf_model_type=1)
global_model_config <- createGlobalModelConfig(global_error_variance=1.0)
forest_model <- createForestModel(forest_dataset, forest_model_config, global_model_config)

Description

Create a forest model config object

This function is intended for advanced use cases in which users require detailed control of sampling algorithms and data structures. Minimal input validation and error checks are performed – users are responsible for providing the correct inputs. For tutorials on the "proper" usage of the stochtree's advanced workflow, we provide several vignettes at https://stochtree.ai/

Usage

createForestModelConfig(
  feature_types = NULL,
  sweep_update_indices = NULL,
  num_trees = NULL,
  num_features = NULL,
  num_observations = NULL,
  variable_weights = NULL,
  leaf_dimension = 1,
  alpha = 0.95,
  beta = 2,
  min_samples_leaf = 5,
  max_depth = -1,
  leaf_model_type = 1,
  leaf_model_scale = NULL,
  variance_forest_shape = 1,
  variance_forest_scale = 1,
  cloglog_forest_shape = 2,
  cloglog_forest_rate = 2,
  cutpoint_grid_size = 100,
  num_features_subsample = NULL
)

Arguments

Value

ForestModelConfig object

Examples

config <- createForestModelConfig(num_trees = 10, num_features = 5, num_observations = 100)

Description

Create a container of forest samples

This function is intended for advanced use cases in which users require detailed control of sampling algorithms and data structures. Minimal input validation and error checks are performed – users are responsible for providing the correct inputs. For tutorials on the "proper" usage of the stochtree's advanced workflow, we provide several vignettes at https://stochtree.ai/

Usage

createForestSamples(
  num_trees,
  leaf_dimension = 1,
  is_leaf_constant = FALSE,
  is_exponentiated = FALSE
)

Arguments

Value

ForestSamples object

Examples

num_trees <- 100
leaf_dimension <- 2
is_leaf_constant <- FALSE
is_exponentiated <- FALSE
forest_samples <- createForestSamples(num_trees, leaf_dimension, is_leaf_constant, is_exponentiated)

Description

Create a global model config object

This function is intended for advanced use cases in which users require detailed control of sampling algorithms and data structures. Minimal input validation and error checks are performed – users are responsible for providing the correct inputs. For tutorials on the "proper" usage of the stochtree's advanced workflow, we provide several vignettes at https://stochtree.ai/

Usage

createGlobalModelConfig(global_error_variance = 1)

Arguments

Value

GlobalModelConfig object

Examples

config <- createGlobalModelConfig(global_error_variance = 100)

Description

Create an outcome object

This function is intended for advanced use cases in which users require detailed control of sampling algorithms and data structures. Minimal input validation and error checks are performed – users are responsible for providing the correct inputs. For tutorials on the "proper" usage of the stochtree's advanced workflow, we provide several vignettes at https://stochtree.ai/

Usage

createOutcome(outcome)

Arguments

Value

Outcome object

Examples

X <- matrix(runif(10*100), ncol = 10)
y <- -5 + 10*(X[,1] > 0.5) + rnorm(100)
outcome <- createOutcome(y)

Description

Create a RandomEffectSamples object

This function is intended for advanced use cases in which users require detailed control of sampling algorithms and data structures. Minimal input validation and error checks are performed – users are responsible for providing the correct inputs. For tutorials on the "proper" usage of the stochtree's advanced workflow, we provide several vignettes at https://stochtree.ai/

Usage

createRandomEffectSamples(num_components, num_groups, random_effects_tracker)

Arguments

Value

RandomEffectSamples object

Examples

n <- 100
rfx_group_ids <- sample(1:2, size = n, replace = TRUE)
rfx_basis <- matrix(rep(1.0, n), ncol=1)
num_groups <- length(unique(rfx_group_ids))
num_components <- ncol(rfx_basis)
rfx_tracker <- createRandomEffectsTracker(rfx_group_ids)
rfx_samples <- createRandomEffectSamples(num_components, num_groups, rfx_tracker)

Description

Create a random effects dataset object

This function is intended for advanced use cases in which users require detailed control of sampling algorithms and data structures. Minimal input validation and error checks are performed – users are responsible for providing the correct inputs. For tutorials on the "proper" usage of the stochtree's advanced workflow, we provide several vignettes at https://stochtree.ai/

Usage

createRandomEffectsDataset(group_labels, basis, variance_weights = NULL)

Arguments

Value

RandomEffectsDataset object

Examples

rfx_group_ids <- sample(1:2, size = 100, replace = TRUE)
rfx_basis <- matrix(rnorm(3*100), ncol = 3)
weight_vector <- rnorm(100)
rfx_dataset <- createRandomEffectsDataset(rfx_group_ids, rfx_basis)
rfx_dataset <- createRandomEffectsDataset(rfx_group_ids, rfx_basis, weight_vector)

Description

Create a RandomEffectsModel object

This function is intended for advanced use cases in which users require detailed control of sampling algorithms and data structures. Minimal input validation and error checks are performed – users are responsible for providing the correct inputs. For tutorials on the "proper" usage of the stochtree's advanced workflow, we provide several vignettes at https://stochtree.ai/

Usage

createRandomEffectsModel(num_components, num_groups)

Arguments

Value

RandomEffectsModel object

Examples

n <- 100
rfx_group_ids <- sample(1:2, size = n, replace = TRUE)
rfx_basis <- matrix(rep(1.0, n), ncol=1)
num_groups <- length(unique(rfx_group_ids))
num_components <- ncol(rfx_basis)
rfx_model <- createRandomEffectsModel(num_components, num_groups)

Description

Create a RandomEffectsTracker object

This function is intended for advanced use cases in which users require detailed control of sampling algorithms and data structures. Minimal input validation and error checks are performed – users are responsible for providing the correct inputs. For tutorials on the "proper" usage of the stochtree's advanced workflow, we provide several vignettes at https://stochtree.ai/

Usage

createRandomEffectsTracker(rfx_group_indices)

Arguments

Value

RandomEffectsTracker object

Examples

n <- 100
rfx_group_ids <- sample(1:2, size = n, replace = TRUE)
rfx_basis <- matrix(rep(1.0, n), ncol=1)
num_groups <- length(unique(rfx_group_ids))
num_components <- ncol(rfx_basis)
rfx_tracker <- createRandomEffectsTracker(rfx_group_ids)

Description

The functions in this group are designed to handle data preprocessing for stochastic forest models. For example, factor-valued columns in data frames are either one-hot encoded or converted to integer indices before the dataframe is converted to a standard matrix format for sampling. This preprocessing routine defines a set of "steps" that must be repeated on out-of-sample datasets before predictions can be obtained from a sampling model.

preprocessTrainData preprocesses covariates for the forest sampler routines, depending on the input type. DataFrames will be preprocessed based on their column types (numeric columns are not modified, ordered factors are integer coded, and unordered factors are one-hot encoded). Matrices are unmodified (assuming all columns are numeric). This function also records and returns a "metadata" list with preprocessing details to ensure that other datasets can be preprocessing identically.

preprocessPredictionData preprocesses covariates for the forest sampler routines, based on the steps outlined in a metadata list produced by preprocessTrainData.

These procedures are handled internally in the bart() and bcf() functions, but they are provided in stochtree as convenience functions for users writing custom samplers. Furthermore, while R lists can be serialized to RDS format, we offer a number of JSON serialization routines for the metadata list produced by preprocessTrainData for consistency with the broader serialization approach of stochtree (see BARTSerialization and BCFSerialization).

Following the API for serializing bartmodel and bcfmodel objects, we can convert metadata to JSON or JSON strings via savePreprocessorToJson and savePreprocessorToJsonString. Similarly, we can reload a metadata list from JSON or JSON strings via createPreprocessorFromJson and createPreprocessorFromJsonString.

Usage

preprocessTrainData(input_data)

preprocessPredictionData(input_data, metadata)

savePreprocessorToJson(object)

savePreprocessorToJsonString(object)

createPreprocessorFromJson(json_object)

createPreprocessorFromJsonString(json_string)

Arguments

Value

preprocessTrainData returns a list with transformed matrix data and a "metadata" list with details on the preprocessing procedures applied. preprocessPredictionData returns a matrix reflecting the data transformations specified in the provided metadata list.

savePreprocessorToJson return an object of type CppJson. savePreprocessorToJsonString returns a string dump of the preprocessor's JSON representation.

createPreprocessorFromJson and createPreprocessorFromJsonString both return metadata lists.

Examples

# Check that running the same data through `preprocessTrainData`
# and `preprocessPredictionData` yields the same result
n <- 100
x1 <- rnorm(n)
x2 <- factor(sample(1:3, n, replace = TRUE), ordered = TRUE)
x3 <- factor(sample(1:3, n, replace = TRUE), ordered = FALSE)
df1 <- data.frame(x1 = x1, x2 = x2, x3 = x3)
df2 <- data.frame(x1 = x1, x2 = x2, x3 = x3)
preprocess_train_list <- preprocessTrainData(df1)
df1_process <- preprocess_train_list$data
df1_metadata <- preprocess_train_list$metadata
df2_process <- preprocessPredictionData(df2, df1_metadata)
all.equal(df1_process, df2_process)

# Save to in-memory JSON
metadata_json <- savePreprocessorToJson(df1_metadata)
# Save to JSON string
metadata_json_string <- savePreprocessorToJsonString(df1_metadata)

# Reload metadata list from in-memory JSON object
metadata_roundtrip <- createPreprocessorFromJson(metadata_json)
# Reload metadata list from JSON string
metadata_roundtrip <- createPreprocessorFromJsonString(metadata_json_string)

cov_df <- data.frame(x1 = 1:5, x2 = 5:1, x3 = 6:10)
metadata <- list(num_ordered_cat_vars = 0, num_unordered_cat_vars = 0,
                 num_numeric_vars = 3, numeric_vars = c("x1", "x2", "x3"))
X_preprocessed <- preprocessPredictionData(cov_df, metadata)

Description

Generic function for extracting parameter samples from a model object (BCF, BART, etc...)

Usage

extractParameter(object, term)

Arguments

Value

Parameter sample array

Examples

n <- 100
p <- 10
X <- matrix(runif(n*p), ncol = p)
rfx_group_ids <- sample(1:2, size = n, replace = TRUE)
rfx_basis <- rep(1.0, n)
y <- (-5 + 10*(X[,1] > 0.5)) + (-2*(rfx_group_ids==1)+2*(rfx_group_ids==2)) + rnorm(n)
bart_model <- bart(X_train=X, y_train=y, rfx_group_ids_train=rfx_group_ids,
                   rfx_basis_train = rfx_basis, num_gfr=0, num_mcmc=10)
sigma2_samples <- extractParameter(bart_model, "sigma2")

Description

Extract a vector, matrix or array of parameter samples from a BART model by name. Random effects are handled by a separate getRandomEffectSamples function due to the complexity of the random effects parameters. If the requested model term is not found, an error is thrown. The following conventions are used for parameter names:

• Global error variance: "sigma2", "global_error_scale", "sigma2_global"

• Leaf scale: "sigma2_leaf", "leaf_scale"

• In-sample mean function predictions: "y_hat_train"

• Test set mean function predictions: "y_hat_test"

• In-sample variance forest predictions: "sigma2_x_train", "var_x_train"

• Test set variance forest predictions: "sigma2_x_test", "var_x_test"

• Ordinal model cutpoints (valid only for ordinal cloglog models): "cloglog_cutpoints", "cutpoints"

Usage

## S3 method for class 'bartmodel'
extractParameter(object, term)

Arguments

Value

Array of parameter samples. If the underlying parameter is a scalar, this will be a vector of length num_samples. If the underlying parameter is vector-valued, this will be (parameter_dimension x num_samples) matrix, and if the underlying parameter is multidimensional, this will be an array of dimension (parameter_dimension_1 x parameter_dimension_2 x ... x num_samples).

Examples

n <- 100
p <- 5
X <- matrix(runif(n*p), ncol = p)
f_XW <- (
    ((0 <= X[,1]) & (0.25 > X[,1])) * (-7.5) +
    ((0.25 <= X[,1]) & (0.5 > X[,1])) * (-2.5) +
    ((0.5 <= X[,1]) & (0.75 > X[,1])) * (2.5) +
    ((0.75 <= X[,1]) & (1 > X[,1])) * (7.5)
)
snr <- 3
group_ids <- rep(c(1,2), n %/% 2)
rfx_coefs <- matrix(c(-1, -1, 1, 1), nrow=2, byrow=TRUE)
rfx_basis <- cbind(1, runif(n, -1, 1))
rfx_term <- rowSums(rfx_coefs[group_ids,] * rfx_basis)
E_y <- f_XW + rfx_term
y <- E_y + rnorm(n, 0, 1)*(sd(E_y)/snr)
test_set_pct <- 0.2
n_test <- round(test_set_pct*n)
n_train <- n - n_test
test_inds <- sort(sample(1:n, n_test, replace = FALSE))
train_inds <- (1:n)[!((1:n) %in% test_inds)]
X_test <- X[test_inds,]
X_train <- X[train_inds,]
y_test <- y[test_inds]
y_train <- y[train_inds]
rfx_group_ids_test <- group_ids[test_inds]
rfx_group_ids_train <- group_ids[train_inds]
rfx_basis_test <- rfx_basis[test_inds,]
rfx_basis_train <- rfx_basis[train_inds,]
rfx_term_test <- rfx_term[test_inds]
rfx_term_train <- rfx_term[train_inds]
bart_model <- bart(X_train = X_train, y_train = y_train, X_test = X_test,
                   rfx_group_ids_train = rfx_group_ids_train,
                   rfx_group_ids_test = rfx_group_ids_test,
                   rfx_basis_train = rfx_basis_train,
                   rfx_basis_test = rfx_basis_test,
                   num_gfr = 10, num_burnin = 0, num_mcmc = 10)
sigma2_samples <- extractParameter(bart_model, "sigma2")

Description

Extract a vector, matrix or array of parameter samples from a BCF model by name. Random effects are handled by a separate getRandomEffectSamples function due to the complexity of the random effects parameters. If the requested model term is not found, an error is thrown. The following conventions are used for parameter names:

• Global error variance: "sigma2", "global_error_scale", "sigma2_global"

• Prognostic forest leaf scale: "sigma2_leaf_mu", "leaf_scale_mu", "mu_leaf_scale"

• Treatment effect forest leaf scale: "sigma2_leaf_tau", "leaf_scale_tau", "tau_leaf_scale"

• Adaptive coding parameters: "adaptive_coding" (returns both the control and treated parameters jointly, with control in the first row and treated in the second row)

• In-sample mean function predictions: "y_hat_train"

• Test set mean function predictions: "y_hat_test"

• In-sample treatment effect forest predictions: "tau_hat_train"

• Test set treatment effect forest predictions: "tau_hat_test"

• Treatment effect intercept: "tau_0", "treatment_intercept", "tau_intercept"

• In-sample variance forest predictions: "sigma2_x_train", "var_x_train"

• Test set variance forest predictions: "sigma2_x_test", "var_x_test"

Usage

## S3 method for class 'bcfmodel'
extractParameter(object, term)

Arguments

Value

Array of parameter samples. If the underlying parameter is a scalar, this will be a vector of length num_samples. If the underlying parameter is vector-valued, this will be (parameter_dimension x num_samples) matrix, and if the underlying parameter is multidimensional, this will be an array of dimension (parameter_dimension_1 x parameter_dimension_2 x ... x num_samples).

Examples

n <- 500
p <- 5
X <- matrix(runif(n*p), ncol = p)
mu_x <- (
    ((0 <= X[,1]) & (0.25 > X[,1])) * (-7.5) +
    ((0.25 <= X[,1]) & (0.5 > X[,1])) * (-2.5) +
    ((0.5 <= X[,1]) & (0.75 > X[,1])) * (2.5) +
    ((0.75 <= X[,1]) & (1 > X[,1])) * (7.5)
)
pi_x <- (
    ((0 <= X[,1]) & (0.25 > X[,1])) * (0.2) +
    ((0.25 <= X[,1]) & (0.5 > X[,1])) * (0.4) +
    ((0.5 <= X[,1]) & (0.75 > X[,1])) * (0.6) +
    ((0.75 <= X[,1]) & (1 > X[,1])) * (0.8)
)
tau_x <- (
    ((0 <= X[,2]) & (0.25 > X[,2])) * (0.5) +
    ((0.25 <= X[,2]) & (0.5 > X[,2])) * (1.0) +
    ((0.5 <= X[,2]) & (0.75 > X[,2])) * (1.5) +
    ((0.75 <= X[,2]) & (1 > X[,2])) * (2.0)
)
Z <- rbinom(n, 1, pi_x)
E_XZ <- mu_x + Z*tau_x
snr <- 3
rfx_group_ids <- rep(c(1,2), n %/% 2)
rfx_coefs <- matrix(c(-1, -1, 1, 1), nrow=2, byrow=TRUE)
rfx_basis <- cbind(1, runif(n, -1, 1))
rfx_term <- rowSums(rfx_coefs[rfx_group_ids,] * rfx_basis)
y <- E_XZ + rfx_term + rnorm(n, 0, 1)*(sd(E_XZ)/snr)
bcf_model <- bcf(X_train = X, y_train = y, Z_train = Z,
                 rfx_group_ids_train = rfx_group_ids,
                 rfx_basis_train = rfx_basis,
                 num_gfr = 10, num_burnin = 0, num_mcmc = 10)
sigma2_samples <- extractParameter(bcf_model, "sigma2")

Description

Wrapper around a C++ class that stores a single ensemble of decision trees (often treated as the "active forest" / current state of a forest term in a sampling loop in R)

This class is intended for advanced use cases in which users require detailed control of sampling algorithms and data structures. Minimal input validation and error checks are performed – users are responsible for providing the correct inputs. For tutorials on the "proper" usage of the stochtree's advanced workflow, we provide several vignettes at https://stochtree.ai/

Public fields

Methods

Public methods

• Forest$new()

• Forest$merge_forest()

• Forest$add_constant()

• Forest$multiply_constant()

• Forest$predict()

• Forest$predict_raw()

• Forest$set_root_leaves()

• Forest$prepare_for_sampler()

• Forest$adjust_residual()

• Forest$num_trees()

• Forest$leaf_dimension()

• Forest$is_leaf_constant()

• Forest$is_exponentiated()

• Forest$add_numeric_split_tree()

• Forest$get_tree_leaves()

• Forest$get_tree_split_counts()

• Forest$get_forest_split_counts()

• Forest$tree_max_depth()

• Forest$average_max_depth()

• Forest$is_empty()

Method new()

Create a new Forest object.

Usage

Forest$new(
  num_trees,
  leaf_dimension = 1,
  is_leaf_constant = FALSE,
  is_exponentiated = FALSE
)

Arguments

Returns

A new Forest object.

Method merge_forest()

Create a larger forest by merging the trees of this forest with those of another forest

Usage

Forest$merge_forest(forest)

Arguments

Method add_constant()

Add a constant value to every leaf of every tree in an ensemble. If leaves are multi-dimensional, constant_value will be added to every dimension of the leaves.

Usage

Forest$add_constant(constant_value)

Arguments

Method multiply_constant()

Multiply every leaf of every tree by a constant value. If leaves are multi-dimensional, constant_multiple will be multiplied through every dimension of the leaves.

Usage

Forest$multiply_constant(constant_multiple)

Arguments

Method predict()

Predict forest on every sample in forest_dataset

Usage

Forest$predict(forest_dataset)

Arguments

Returns

vector of predictions with as many rows as in forest_dataset

Method predict_raw()

Predict "raw" leaf values (without being multiplied by basis) for every sample in forest_dataset

Usage

Forest$predict_raw(forest_dataset)

Arguments

Returns

Array of predictions for each observation in forest_dataset and each sample in the ForestSamples class with each prediction having the dimensionality of the forests' leaf model. In the case of a constant leaf model or univariate leaf regression, this array is a vector (length is the number of observations). In the case of a multivariate leaf regression, this array is a matrix (number of observations by leaf model dimension, number of samples).

Method set_root_leaves()

Set a constant predicted value for every tree in the ensemble. Stops program if any tree is more than a root node.

Usage

Forest$set_root_leaves(leaf_value)

Arguments

Method prepare_for_sampler()

Set a constant predicted value for every tree in the ensemble. Stops program if any tree is more than a root node.

Usage

Forest$prepare_for_sampler(
  dataset,
  outcome,
  forest_model,
  leaf_model_int,
  leaf_value
)

Arguments

Method adjust_residual()

Adjusts residual based on the predictions of a forest

This is typically run just once at the beginning of a forest sampling algorithm. After trees are initialized with constant root node predictions, their root predictions are subtracted out of the residual.

Usage

Forest$adjust_residual(dataset, outcome, forest_model, requires_basis, add)

Arguments

Method num_trees()

Return number of trees in each ensemble of a Forest object

Usage

Forest$num_trees()

Returns

Tree count

Method leaf_dimension()

Return output dimension of trees in a Forest object

Usage

Forest$leaf_dimension()

Returns

Leaf node parameter size

Method is_leaf_constant()

Return constant leaf status of trees in a Forest object

Usage

Forest$is_leaf_constant()

Returns

TRUE if leaves are constant, FALSE otherwise

Method is_exponentiated()

Return exponentiation status of trees in a Forest object

Usage

Forest$is_exponentiated()

Returns

TRUE if leaf predictions must be exponentiated, FALSE otherwise

Method add_numeric_split_tree()

Add a numeric (i.e. X[,i] <= c) split to a given tree in the ensemble

Usage

Forest$add_numeric_split_tree(
  tree_num,
  leaf_num,
  feature_num,
  split_threshold,
  left_leaf_value,
  right_leaf_value
)

Arguments

Method get_tree_leaves()

Retrieve a vector of indices of leaf nodes for a given tree in a given forest

Usage

Forest$get_tree_leaves(tree_num)

Arguments

Method get_tree_split_counts()

Retrieve a vector of split counts for every training set variable in a given tree in the forest

Usage

Forest$get_tree_split_counts(tree_num, num_features)

Arguments

Method get_forest_split_counts()

Retrieve a vector of split counts for every training set variable in the forest

Usage

Forest$get_forest_split_counts(num_features)

Arguments

Method tree_max_depth()

Maximum depth of a specific tree in the forest

Usage

Forest$tree_max_depth(tree_num)

Arguments

Returns

Maximum leaf depth

Method average_max_depth()

Average the maximum depth of each tree in the forest

Usage

Forest$average_max_depth()

Returns

Average maximum depth

Method is_empty()

When a forest object is created, it is "empty" in the sense that none of its component trees have leaves with values. There are two ways to "initialize" a Forest object. First, the set_root_leaves() method simply initializes every tree in the forest to a single node carrying the same (user-specified) leaf value. Second, the prepare_for_sampler() method initializes every tree in the forest to a single node with the same value and also propagates this information through to a ForestModel object, which must be synchronized with a Forest during a forest sampler loop.

Usage

Forest$is_empty()

Returns

TRUE if a Forest has not yet been initialized with a constant root value, FALSE otherwise if the forest has already been initialized / grown.

Description

Wrapper around a C++ dataset class used to sample a forest. A dataset consists of three matrices / vectors: covariates, bases, and variance weights. Both the basis vector and variance weights are optional.

This class is intended for advanced use cases in which users require detailed control of sampling algorithms and data structures. Minimal input validation and error checks are performed – users are responsible for providing the correct inputs. For tutorials on the "proper" usage of the stochtree's advanced workflow, we provide several vignettes at https://stochtree.ai/

Public fields

Methods

Public methods

• ForestDataset$new()

• ForestDataset$update_basis()

• ForestDataset$update_variance_weights()

• ForestDataset$num_observations()

• ForestDataset$num_covariates()

• ForestDataset$num_basis()

• ForestDataset$get_covariates()

• ForestDataset$get_basis()

• ForestDataset$get_variance_weights()

• ForestDataset$has_basis()

• ForestDataset$has_variance_weights()

• ForestDataset$has_auxiliary_dimension()

• ForestDataset$add_auxiliary_dimension()

• ForestDataset$get_auxiliary_data_value()

• ForestDataset$set_auxiliary_data_value()

• ForestDataset$get_auxiliary_data_vector()

Method new()

Create a new ForestDataset object.

Usage

ForestDataset$new(covariates, basis = NULL, variance_weights = NULL)

Arguments

Returns

A new ForestDataset object.

Method update_basis()

Update basis matrix in a dataset

Usage

ForestDataset$update_basis(basis)

Arguments

Method update_variance_weights()

Update variance_weights in a dataset

Usage

ForestDataset$update_variance_weights(variance_weights, exponentiate = F)

Arguments

Method num_observations()

Return number of observations in a ForestDataset object

Usage

ForestDataset$num_observations()

Returns

Observation count

Method num_covariates()

Return number of covariates in a ForestDataset object

Usage

ForestDataset$num_covariates()

Returns

Covariate count

Method num_basis()

Return number of bases in a ForestDataset object

Usage

ForestDataset$num_basis()

Returns

Basis count

Method get_covariates()

Return covariates as an R matrix

Usage

ForestDataset$get_covariates()

Returns

Covariate data

Method get_basis()

Return bases as an R matrix

Usage

ForestDataset$get_basis()

Returns

Basis data

Method get_variance_weights()

Return variance weights as an R vector

Usage

ForestDataset$get_variance_weights()

Returns

Variance weight data

Method has_basis()

Whether or not a dataset has a basis matrix

Usage

ForestDataset$has_basis()

Returns

True if basis matrix is loaded, false otherwise

Method has_variance_weights()

Whether or not a dataset has variance weights

Usage

ForestDataset$has_variance_weights()

Returns

True if variance weights are loaded, false otherwise

Method has_auxiliary_dimension()

Whether or not a dataset has auxiliary data stored at the dimension indicated

Usage

ForestDataset$has_auxiliary_dimension(dim_idx)

Arguments

Returns

True if auxiliary data has been allocated for dim_idx False otherwise

Method add_auxiliary_dimension()

Initialize a new dimension / lane of auxiliary data and allocate data in its place

Usage

ForestDataset$add_auxiliary_dimension(dim_size)

Arguments

Returns

None

Method get_auxiliary_data_value()

Retrieve auxiliary data value

Usage

ForestDataset$get_auxiliary_data_value(dim_idx, element_idx)

Arguments

Returns

Floating point value stored in the requested auxiliary data space

Method set_auxiliary_data_value()

Set auxiliary data value

Usage

ForestDataset$set_auxiliary_data_value(dim_idx, element_idx, value)

Arguments

Returns

None

Method get_auxiliary_data_vector()

Retrieve entire auxiliary data vector

Usage

ForestDataset$get_auxiliary_data_vector(dim_idx)

Arguments

Returns

Vector of all of the auxiliary data stored at dimension dim_idx

Description

Decision tree ensembles can be represented in part by a "kernel" function whose distance metric is based on the extent to which two observations are mapped to the same leaf nodes. This function group offers utilities for evaluating this kernel.

computeForestLeafIndices computes and return a vector representation of a forest's leaf predictions for every observation in a dataset. The resulting vector has a "tree-major" format that can be easily re-represented as as a CSR sparse matrix: elements are organized so that the first n elements correspond to leaf predictions for all n observations in a dataset for the first tree in an ensemble, the next n elements correspond to predictions for the second tree and so on. The "data" for each element corresponds to a uniquely mapped column index that corresponds to a single leaf of a single tree (i.e. if tree 1 has 3 leaves, its column indices range from 0 to 2, and then tree 2's leaf indices begin at 3, etc...).

computeForestLeafVariances returns each forest's leaf node scale parameters. If leaf scale is not sampled for the forest in question, the function throws an error that the leaf model does not have a stochastic scale parameter.

computeForestMaxLeafIndex computes and returns the largest possible leaf index computable by computeForestLeafIndices for the forests in a designated forest sample container.

These functions are intended for advanced use cases in which users require detailed control of sampling algorithms and data structures. Minimal input validation and error checks are performed – users are responsible for providing the correct inputs. For tutorials on the "proper" usage of the stochtree's advanced workflow, we provide several vignettes at https://stochtree.ai/

Usage

computeForestLeafIndices(
  model_object,
  covariates,
  forest_type = NULL,
  propensity = NULL,
  forest_inds = NULL
)

computeForestLeafVariances(model_object, forest_type, forest_inds = NULL)

computeForestMaxLeafIndex(model_object, forest_type = NULL, forest_inds = NULL)

Arguments

Value

computeForestLeafIndices returns a vector of size num_obs * num_trees, where num_obs = nrow(covariates) and num_trees is the number of trees in the relevant forest of model_object.

computeForestLeafVariances returns a vector of size length(forest_inds) with the leaf scale parameter for each requested forest.

computeForestMaxLeafIndex returns a vector containing the largest possible leaf index computable by computeForestLeafIndices for the forests in a designated forest sample container.

Examples

X <- matrix(runif(10*100), ncol = 10)
y <- -5 + 10*(X[,1] > 0.5) + rnorm(100)
bart_model <- bart(X, y, num_gfr=0, num_mcmc=10)
leaf_indices <- computeForestLeafIndices(bart_model, X, "mean")
leaf_indices <- computeForestLeafIndices(bart_model, X, "mean", 0)
leaf_indices <- computeForestLeafIndices(bart_model, X, "mean", c(1,3,9))
leaf_variances <- computeForestLeafVariances(bart_model, "mean")
leaf_variances <- computeForestLeafVariances(bart_model, "mean", 0)
leaf_variances <- computeForestLeafVariances(bart_model, "mean", c(1,3,5))
max_leaf_index <- computeForestMaxLeafIndex(bart_model, "mean")
max_leaf_index <- computeForestMaxLeafIndex(bart_model, "mean", 0)
max_leaf_index <- computeForestMaxLeafIndex(bart_model, "mean", c(1,3,9))

Description

Wraps the C++ data structures needed to sample an ensemble of decision trees and exposes functionality to run a forest sampler (using either MCMC or the grow-from-root algorithm).

This class is intended for advanced use cases in which users require detailed control of sampling algorithms and data structures. Minimal input validation and error checks are performed – users are responsible for providing the correct inputs. For tutorials on the "proper" usage of the stochtree's advanced workflow, we provide several vignettes at https://stochtree.ai/

Public fields

Methods

Public methods

• ForestModel$new()

• ForestModel$sample_one_iteration()

• ForestModel$get_cached_forest_predictions()

• ForestModel$propagate_basis_update()

• ForestModel$propagate_residual_update()

• ForestModel$update_alpha()

• ForestModel$update_beta()

• ForestModel$update_min_samples_leaf()

• ForestModel$update_max_depth()

• ForestModel$get_alpha()

• ForestModel$get_beta()

• ForestModel$get_min_samples_leaf()

• ForestModel$get_max_depth()

Method new()

Create a new ForestModel object.

Usage

ForestModel$new(
  forest_dataset,
  feature_types,
  num_trees,
  n,
  alpha,
  beta,
  min_samples_leaf,
  max_depth = -1
)

Arguments

Returns

A new ForestModel object.

Method sample_one_iteration()

Run a single iteration of the forest sampling algorithm (MCMC or GFR)

Usage

ForestModel$sample_one_iteration(
  forest_dataset,
  residual,
  forest_samples,
  active_forest,
  rng,
  forest_model_config,
  global_model_config,
  num_threads = -1,
  keep_forest = TRUE,
  gfr = TRUE
)

Arguments

Method get_cached_forest_predictions()

Extract an internally-cached prediction of a forest on the training dataset in a sampler.

Usage

ForestModel$get_cached_forest_predictions()

Returns

Vector with as many elements as observations in the training dataset

Method propagate_basis_update()

Propagates basis update through to the (full/partial) residual by iteratively (a) adding back in the previous prediction of each tree, (b) recomputing predictions for each tree (caching on the C++ side), (c) subtracting the new predictions from the residual.

This is useful in cases where a basis (for e.g. leaf regression) is updated outside of a tree sampler (as with e.g. adaptive coding for binary treatment BCF). Once a basis has been updated, the overall "function" represented by a tree model has changed and this should be reflected through to the residual before the next sampling loop is run.

Usage

ForestModel$propagate_basis_update(dataset, outcome, active_forest)

Arguments

Method propagate_residual_update()

Update the current state of the outcome (i.e. partial residual) data by subtracting the current predictions of each tree. This function is run after the Outcome class's update_data method, which overwrites the partial residual with an entirely new stream of outcome data.

Usage

ForestModel$propagate_residual_update(residual)

Arguments

Returns

None

Method update_alpha()

Update alpha in the tree prior

Usage

ForestModel$update_alpha(alpha)

Arguments

Returns

None

Method update_beta()

Update beta in the tree prior

Usage

ForestModel$update_beta(beta)

Arguments

Returns

None

Method update_min_samples_leaf()

Update min_samples_leaf in the tree prior

Usage

ForestModel$update_min_samples_leaf(min_samples_leaf)

Arguments

Returns

None

Method update_max_depth()

Update max_depth in the tree prior

Usage

ForestModel$update_max_depth(max_depth)

Arguments

Returns

None

Method get_alpha()

Update alpha in the tree prior

Usage

ForestModel$get_alpha()

Returns

Value of alpha in the tree prior

Method get_beta()

Update beta in the tree prior

Usage

ForestModel$get_beta()

Returns

Value of beta in the tree prior

Method get_min_samples_leaf()

Query min_samples_leaf in the tree prior

Usage

ForestModel$get_min_samples_leaf()

Returns

Value of min_samples_leaf in the tree prior

Method get_max_depth()

Query max_depth in the tree prior

Usage

ForestModel$get_max_depth()

Returns

Value of max_depth in the tree prior

Description

Object used to get / set parameters and other model configuration options for a forest model in the "low-level" stochtree interface. The "low-level" stochtree interface enables a high degreee of sampler customization, in which users employ R wrappers around C++ objects like ForestDataset, Outcome, CppRng, and ForestModel to run the Gibbs sampler of a BART model with custom modifications. ForestModelConfig allows users to specify / query the parameters of a forest model they wish to run.

This class is intended for advanced use cases in which users require detailed control of sampling algorithms and data structures. Minimal input validation and error checks are performed – users are responsible for providing the correct inputs. For tutorials on the "proper" usage of the stochtree's advanced workflow, we provide several vignettes at https://stochtree.ai/

Value

Vector of integer-coded feature types (integers where 0 = numeric, 1 = ordered categorical, 2 = unordered categorical)

Vector of (0-indexed) indices of trees to update in a sweep

Vector specifying sampling probability for all p covariates in ForestDataset

Number of trees in a forest

Number of features in a forest model training set

Number of observations in a forest model training set

Root node split probability in tree prior

Depth prior penalty in tree prior

Minimum number of samples in a tree leaf

Maximum depth of any tree in the ensemble in the model

Integer coded leaf model type

Scale parameter used in Gaussian leaf models

Shape parameter for IG leaf models

Scale parameter for IG leaf models

Shape parameter for conditional gamma component of cloglog leaf models

Rate parameter for conditional gamma component of cloglog leaf models

Number of unique cutpoints to consider

Number of features to subsample for the GFR algorithm

Public fields