BayesRegDTR (version 1.1.2)

BayesLinRegDTR.model.fit: Main function for fitting a Bayesian likelihood-based linear regression model

Description

Fits the Bayesian likelihood-based linear model to obtain an estimated posterior distribution of the optimal treatment option at a user-specified prediction stage. Uses backward induction and dynamic programming theory for computing expected values.

Usage

BayesLinRegDTR.model.fit(
  Dat.train,
  Dat.pred,
  n.train,
  n.pred,
  num_stages,
  num_treats,
  p_list,
  t,
  R = 30,
  tau = 0.01,
  B = 10000,
  nu0 = 3,
  V0 = mapply(diag, p_list, SIMPLIFY = FALSE),
  alph = 1,
  gam = 1,
  showBar = TRUE
)

Value

GCV_results

An array of dimension \(n.pred \times num\_treats[t] \times B\), indicating the expected value under each treatment option at stage t.

post.prob

An \(n.pred \times num\_treats[t]\) matrix of the posterior probability that each treatment type at stage t is optimal

MC_draws.train

A list of Monte Carlo draws containing:

• sigmat_B_list - A list of length num_stages with each element a vector of size \(B \times p\_list[t]\)

• Wt_B_list - A list of length num_stages with each element a matrix of size \(B \times p\_list[t]\)

• beta_B - A list of length B

• sigmay_2B - A list of length B

Arguments

Dat.train

Training data in format returned by generate_dataset: organised as a list of \(\{y, X_1, X_2..., X_{num\_stages}, A\}\) where y is a vector of the final outcomes, \(X_1, X_2..., X_{num\_stages}\) is a list of matrices of the intermediate covariates and A is an \(n.train \times num\_stages\) matrix of the assigned treatments, where num_stages is the total number of stages

Dat.pred

Prediction data in format returned by generate_dataset: organised as a list of \(\{X_1, X_2..., X_t, A\}\) where \(X_1, X_2..., X_t\) is a list of matrices of the intermediate covariates and A is an \(n.pred \times (t-1)\) matrix of the assigned treatments, where t is the prediction stage

n.train

Number of samples/individuals in the training data

n.pred

Number of samples/individuals in the prediction data

num_stages

Total number of stages

num_treats

Vector of number of treatment options at each stage

p_list

Vector of intermediate covariate dimensions for each stage

t

Prediction stage t, where t \(\leq\) num_stages

R

Draw size from distribution of intermediate covariates. default: 30

tau

Normal prior scale parameter for regression coefficients. Should be specified with a small value. default: 0.01

B

Number of MC draws from posterior of regression parameters. default 10000

nu0

Inverse-Wishart prior degrees of freedom for regression error Vcov matrix. Ignored if using a univariate dataset. default: 3

V0

List of Inverse-Wishart prior scale matrix for regression error Vcov matrix. Ignored if using a univariate dataset. default: list of identity matrices

alph

Inverse-Gamma prior shape parameter for regression error variance of y. default: 1

gam

Inverse-Gamma prior rate parameter for regression error variance of y. default: 1

showBar

Whether to show a progress bar. Uses API from progressr and future for parallel integration deafult: TRUE

Details

Utilises a future framework, so to enable parallel processing and register a parallel backend, plan and registerDoFuture must be called first.

Additionally, progress bars use progressr API, and a non-default progress bar (e.g. cli) is recommended. See below or registerDoFuture and handlers for examples.

Note that to have a progress bar for the parallel sections, future must be used. To turn off the immediate warnings, use options(BRDTR_warn_imm = FALSE).

Examples

# Code does not run within 10 seconds, so don't run
# \donttest{
# -----------------------------
# Set Up Parallelism & Progress Bar
# -----------------------------
progressr::handlers("cli")          # Set handler to something with title/text
numCores <- parallel::detectCores() # Detect number of cores, use max
future::plan(future::multisession,  # Or plan(multicore, workers) on Unix
            workers = numCores)     # Set number of cores to use
doFuture::registerDoFuture()        # Or doParallel::registerDoParallel()
                                    # if no progress bar is needed and future
                                    # is unwanted

## UVT
# -----------------------------
# Initialise Inputs
# -----------------------------
num_stages  <- 5
t           <- 3
p_list      <- rep(1, num_stages)
num_treats  <- rep(2, num_stages)
n.train     <- 5000
n.pred      <- 10

# -----------------------------
# Generate Dataset
# -----------------------------
Dat.train  <- generate_dataset(n.train,  num_stages, p_list, num_treats)
Dat.pred  <- generate_dataset(n.pred,  num_stages, p_list, num_treats)
Dat.pred  <- Dat.pred[-1]
Dat.pred[[num_stages+1]]  <- Dat.pred[[num_stages+1]][1:n.pred, 1:(t-1), drop = FALSE]

# -----------------------------
# Main
# -----------------------------
gcv_uvt <- BayesLinRegDTR.model.fit(Dat.train, Dat.pred, n.train, n.pred,
                                    num_stages, num_treats,
                                    p_list, t, R = 30,
                                    tau = 0.01, B = 500, nu0 = NULL,
                                    V0 = NULL, alph = 3, gam = 4)

## MVT
# -----------------------------
# Initialise Inputs
# -----------------------------
num_stages  <- 3
t           <- 2
p_list      <- rep(2, num_stages)
num_treats  <- rep(2, num_stages)
n.train     <- 5000
n.pred      <- 10

# -----------------------------
# Generate Dataset
# -----------------------------
Dat.train <- generate_dataset(n.train, num_stages, p_list, num_treats)
Dat.pred  <- generate_dataset(n.pred,  num_stages, p_list, num_treats)
Dat.pred  <- Dat.pred[-1]
Dat.pred[[num_stages+1]]  <- Dat.pred[[num_stages+1]][1:n.pred, 1:(t-1), drop = FALSE]

# -----------------------------
# Main
# -----------------------------
gcv_res <- BayesLinRegDTR.model.fit(Dat.train, Dat.pred, n.train, n.pred,
                                    num_stages, num_treats,
                                    p_list, t, R = 30,
                                    tau = 0.01, B = 500, nu0 = 3,
                                    V0 = mapply(diag, p_list, SIMPLIFY = FALSE),
                                    alph = 3, gam = 4)
# }