bayes_sim_unknownvar {bayesassurance}

Bayesian Simulation with Composite Sampling

Description

Approximates the Bayesian assurance of a one-sided hypothesis test through Monte Carlo sampling with the added assumption that the variance is unknown.

Usage

bayes_sim_unknownvar(
  n,
  p = NULL,
  u,
  C,
  R,
  Xn = NULL,
  Vn,
  Vbeta_d,
  Vbeta_a_inv,
  mu_beta_d,
  mu_beta_a,
  a_sig_a,
  b_sig_a,
  a_sig_d,
  b_sig_d,
  alt = "greater",
  alpha,
  mc_iter
)

Arguments

n

sample size (either vector or scalar)

p

column dimension of design matrix Xn. If Xn = NULL, p must be specified to denote the column dimension of the default design matrix generated by the function.

u

a scalar or vector to evaluate

u'\beta > C,

where \beta is an unknown parameter that is to be estimated.

C

constant value to be compared to when evaluating u'\beta > C

R

number of iterations we want to pass through to check for satisfaction of the analysis stage objective. The proportion of those iterations meeting the analysis objective corresponds to the approximated Bayesian assurance.

Xn

design matrix that characterizes where the data is to be generated from. This is specifically given by the normal linear regression model

yn = Xn\beta + \epsilon,

\epsilon ~ N(0, \sigma^2 Vn),

where \sigma^2 is unknown in this setting. Note that Xn must have column dimension p.

Vn

an n by n correlation matrix for the marginal distribution of the sample data yn. Takes on an identity matrix when set to NULL.

Vbeta_d

correlation matrix that helps describe the prior information on \beta in the design stage

Vbeta_a_inv

inverse-correlation matrix that helps describe the prior information on \beta in the analysis stage

mu_beta_d

design stage mean

mu_beta_a

analysis stage mean

a_sig_a, b_sig_a

shape and scale parameters of the inverse gamma distribution where variance \sigma^2 is sampled from in the analysis stage

a_sig_d, b_sig_d

shape and scale parameters of the inverse gamma distribution where variance \sigma^2 is sampled from in the design stage

alt

specifies alternative test case, where alt = "greater" tests if u'\beta > C, alt = "less" tests if u'\beta < C, and alt = "two.sided" performs a two-sided test. By default, alt = "greater".

alpha

significance level

mc_iter

number of MC samples evaluated under the analysis objective

Value

a list of objects corresponding to the assurance approximations

• assurance_table: table of sample size and corresponding assurance values

• assur_plot: plot of assurance values

• mc_samples: number of Monte Carlo samples that were generated and evaluated

See Also

bayes_sim for the Bayesian assurance function for known variance.

Examples

## O'Hagan and Stevens (2001) include a series of examples with 
## pre-specified parameters that we will be using to replicate their 
## results through our Bayesian assurance simulation.
## The inputs are as follows:

n <- 285
p <- 4 ## includes two parameter measures (cost and efficacy) for each of 
       ## the two treatments, for a total of p = 4 parameters. 
K <- 20000
C <- 0
u <- as.matrix(c(-K, 1, K, -1))

## Set up correlation matrices
Vbeta_a_inv <- matrix(rep(0, p^2), nrow = p, ncol = p)

sigsq <- 4.04^2
tau1 <- tau2 <- 8700
sig <- sqrt(sigsq)
Vn <- matrix(0, nrow = n*p, ncol = n*p)
Vn[1:n, 1:n] <- diag(n)
Vn[(2*n - (n-1)):(2*n), (2*n - (n-1)):(2*n)] <- (tau1 / sig)^2 * diag(n)
Vn[(3*n - (n-1)):(3*n), (3*n - (n-1)):(3*n)] <- diag(n)
Vn[(4*n - (n-1)):(4*n), (4*n - (n-1)):(4*n)] <- (tau2 / sig)^2 * diag(n)

Vbeta_d <- (1 / sigsq) * 
matrix(c(4, 0, 3, 0, 0, 10^7, 0, 0, 3, 0, 4, 0, 0, 0, 0, 10^7),
nrow = 4, ncol = 4)
mu_beta_d <- as.matrix(c(5, 6000, 6.5, 7200))
mu_beta_a <- as.matrix(rep(0, p))
alpha <- 0.05
epsilon <- 10e-7
a_sig_d <- (sigsq / epsilon) + 2
b_sig_d <- sigsq * (a_sig_d - 1)
a_sig_a <- -p / 2
b_sig_a <- 0

bayesassurance::bayes_sim_unknownvar(n = n, p = 4, 
u = as.matrix(c(-K, 1, K, -1)), C = 0, R = 40,
Xn = NULL, Vn = Vn, Vbeta_d = Vbeta_d, 
Vbeta_a_inv = Vbeta_a_inv, mu_beta_d = mu_beta_d,
mu_beta_a = mu_beta_a, a_sig_a = a_sig_a, b_sig_a = b_sig_a, 
a_sig_d = a_sig_d, b_sig_d = b_sig_d, alt = "two.sided", alpha = 0.05, 
mc_iter = 500)