power.glm.random.a0: Power/type I error calculation for generalized linear models...

Description Usage Arguments Details Value References See Also Examples

View source: R/main_func.R

Description

Power/type I error calculation using normalized power priors for generalized linear models with random a_0

Usage

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
power.glm.random.a0(
  data.type,
  data.link,
  data.size,
  n = 1,
  historical,
  samp.prior.beta,
  samp.prior.var,
  prior.a0.shape1 = 1,
  prior.a0.shape2 = 1,
  a0.coefficients,
  lower.limits = rep(-100, 50),
  upper.limits = rep(100, 50),
  slice.widths = rep(0.1, 50),
  delta = 0,
  gamma = 0.95,
  nMC = 10000,
  nBI = 250,
  N = 10000
)

Arguments

data.type

Character string specifying the type of response. The options are "Normal", "Bernoulli", "Binomial", "Poisson" and "Exponential".

data.link

Character string specifying the link function. The options are "Logistic", "Probit", "Log", "Identity-Positive", "Identity-Probability" and "Complementary Log-Log". Does not apply if data.type is "Normal".

data.size

Sample size of the simulated datasets.

n

(For binomial data only) vector of integers specifying the number of subjects who have a particular value of the covariate vector. If the data is binary and all covariates are discrete, collapsing Bernoulli data into a binomial structure can make the slice sampler much faster.

historical

List of historical dataset(s). East historical dataset is stored in a list which contains two named elements: y0 and x0.

  • y0 is a vector of responses.

  • x0 is a matrix of covariates. x0 should NOT have the treatment indicator. Apart from missing the treatent/control indicator, x0 should have the same set of covariates in the same order as x.

For binomial data, an additional element n0 is required.

  • n0 is vector of integers specifying the number of subjects who have a particular value of the covariate vector.

samp.prior.beta

Matrix of possible values of β to sample (with replacement) from. Each row is a possible β vector (a realization from the sampling prior for β), where the first element is the coefficient for the intercept and the second element is the coefficient for the treatment indicator. The length of the vector should be equal to the total number of parameters, i.e. P+2 where P is the number of columns of x0 in historical.

samp.prior.var

Vector of possible values of σ^2 to sample (with replacement) from. Only applies if data.type is "Normal". The vector contains realizations from the sampling prior (e.g. inverse-gamma distribution) for σ^2.

prior.a0.shape1

First shape parameter of the beta prior for a_0. The default is 1.

prior.a0.shape2

Second shape parameter of the beta prior for a_0. The default is 1.

a0.coefficients

Vector of coefficients for a_0 returned by the function normalizing.constant. This is necessary for estimating the normalizing constant for the normalized power prior. Does not apply if data.type is "Normal".

lower.limits

Vector of lower limits for parameters to be used by the slice sampler. If data.type is "Normal", slice sampling is used for a_0, and the length of the vector should be equal to the number of historical datasets. For all other data types, slice sampling is used for β and a_0. The first P+1 elements apply to the sampling of β and the rest apply to the sampling of a_0. The length of the vector should be equal to the sum of the total number of parameters (i.e. P+1 where P is the number of covariates) and the number of historical datasets. The default is -100 for all parameters (may not be appropriate for all situations).

upper.limits

Vector of upper limits for parameters to be used by the slice sampler. If data.type is "Normal", slice sampling is used for a_0, and the length of the vector should be equal to the number of historical datasets. For all other data types, slice sampling is used for β and a_0. The first P+1 elements apply to the sampling of β and the rest apply to the sampling of a_0. The length of the vector should be equal to the sum of the total number of parameters (i.e. P+1 where P is the number of covariates) and the number of historical datasets. The default is 100 for all parameters (may not be appropriate for all situations).

slice.widths

Vector of initial slice widths used by the slice sampler. If data.type is "Normal", slice sampling is used for a_0, and the length of the vector should be equal to the number of historical datasets. For all other data types, slice sampling is used for β and a_0. The first P+1 elements apply to the sampling of β and the rest apply to the sampling of a_0. The length of the vector should be equal to the sum of the total number of parameters (i.e. P+1 where P is the number of covariates) and the number of historical datasets. The default is 0.1 for all parameter (may not be appropriate for all situations).

delta

Prespecified constant that defines the boundary of the null hypothesis. The default is zero.

gamma

Posterior probability threshold for rejecting the null. The null hypothesis is rejected if posterior probability is greater gamma. The default is 0.95.

nMC

Number of iterations (excluding burn-in samples) for the slice sampler or Gibbs sampler. The default is 10,000.

nBI

Number of burn-in samples for the slice sampler or Gibbs sampler. The default is 250.

N

Number of simulated datasets to generate. The default is 10,000.

Details

The user should use the function normalizing.constant to obtain a0.coefficients (does not apply if data.type is "Normal").

samp.prior.beta can be generated using the sampling priors (see example). samp.prior.var is necessary for generating normally distributed data.

If data.type is "Normal", the response y_i is assumed to follow N(x_i'β, τ^{-1}) where x_i is the vector of covariates for subject i. Historical datasets are assumed to have the same precision parameter as the current dataset for computational simplicity. The initial prior for τ is the Jeffery's prior, τ^{-1}. The initial prior for β is the uniform improper prior. Posterior samples for β and τ are obtained through Gibbs sampling. Posterior samples for a_0 are obtained through slice sampling.

For all other data types, posterior samples are obtained through slice sampling. The initial prior for β is the uniform improper prior. The default lower limits for the parameters are -100. The default upper limits for the parameters are 100. The default slice widths for the parameters are 0.1. The defaults may not be appropriate for all situations, and the user can specify the appropriate limits and slice width for each parameter.

If a sampling prior with support in the null space is used, the value returned is a Bayesian type I error rate. If a sampling prior with support in the alternative space is used, the value returned is a Bayesian power.

Value

Power or type I error is returned, depending on the sampling prior used. If data.type is "Normal", average posterior means of β, τ and a_0 are also returned. For all other data types, average posterior means of β and a_0 are also returned.

References

Chen, Ming-Hui, et al. "Bayesian design of noninferiority trials for medical devices using historical data." Biometrics 67.3 (2011): 1163-1170.

Neal, Radford M. Slice sampling. Ann. Statist. 31 (2003), no. 3, 705–767.

See Also

normalizing.constant and glm.random.a0

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
data.type <- "Bernoulli"
data.link <- "Logistic"
data.size <- 100

# Simulate two historical datasets
p <- 3
historical <- list(list(y0=rbinom(data.size,size=1,prob=0.2),
                        x0=matrix(rnorm(p*data.size),ncol=p,nrow=data.size)),
                   list(y0=rbinom(data.size, size=1, prob=0.5),
                        x0=matrix(rnorm(p*data.size),ncol=p,nrow=data.size)))

# Generate sampling priors

# The null hypothesis here is H0: beta_1 >= 0. To calculate power,
# we can provide samples of beta_1 such that the mass of beta_1 < 0.
# To calculate type I error, we can provide samples of beta_1 such that
# the mass of beta_1 >= 0.
samp.prior.beta1 <- rnorm(100, mean=-3, sd=1)
# Here, mass is put on the alternative region, so power is calculated.
samp.prior.beta <- cbind(rnorm(100), samp.prior.beta1, matrix(rnorm(100*p), 100, p))

# Please see function "normalizing.constant" for how to obtain a0.coefficients
# Here, suppose one-degree polynomial regression is chosen by the "normalizing.constant" 
# function. The coefficients are obtained for the intercept, a0_1 and a0_2. 
a0.coefficients <- c(1, 0.5, -1)

nMC <- 100 # nMC should be larger in practice
nBI <- 50
N <- 3 # N should be larger in practice
result <- power.glm.random.a0(data.type=data.type, data.link=data.link,
                              data.size=data.size, historical=historical,
                              samp.prior.beta=samp.prior.beta, a0.coefficients=a0.coefficients,
                              delta=0, nMC=nMC, nBI=nBI, N=N)

BayesPPD documentation built on Sept. 8, 2021, 5:06 p.m.