sensbayescomp: Verifying Optimality of Bayesian Compound DP-optimal Designs

In ICAOD: Optimal Designs for Nonlinear Statistical Models by Imperialist Competitive Algorithm (ICA)

Description

This function plot the sensitivity (derivative) function given an approximate (continuous) design and calculate the efficiency lower bound (ELB) for Bayesian DP-optimal designs. Let x belongs to χ that denotes the design space. Based on the general equivalence theorem, generally, a design ξ* is optimal if and only if the value of its sensitivity (derivative) function be non-positive for all x in χ and it only reaches zero when x belong to the support of ξ* (be equal to one of the design point). Therefore, the user can look at the sensitivity plot and the ELB and decide whether the design is optimal or close enough to the true optimal design (ELB tells us that without knowing the latter).

Usage

sensbayescomp(
  formula,
  predvars,
  parvars,
  family = gaussian(),
  x,
  w,
  lx,
  ux,
  fimfunc = NULL,
  prior = list(),
  prob,
  alpha,
  sens.control = list(),
  sens.bayes.control = list(),
  crt.bayes.control = list(),
  plot_3d = c("lattice", "rgl"),
  plot_sens = TRUE,
  npar = NULL,
  calculate_criterion = TRUE,
  silent = FALSE
)

Arguments

formula	A linear or nonlinear model formula. A symbolic description of the model consists of predictors and the unknown model parameters. Will be coerced to a formula if necessary.
predvars	A vector of characters. Denotes the predictors in the formula.
parvars	A vector of characters. Denotes the unknown parameters in the formula.
family	A description of the response distribution and the link function to be used in the model. This can be a family function, a call to a family function or a character string naming the family. Every family function has a link argument allowing to specify the link function to be applied on the response variable. If not specified, default links are used. For details see family. By default, a linear gaussian model gaussian() is applied.
x	A vector of candidate design (support) points. When is not set to NULL (default), the algorithm only finds the optimal weights for the candidate points in x. Should be set when the user has a finite number of candidate design points and the purpose is to find the optimal weight for each of them (when zero, they will be excluded from the design). For design points with more than one dimension, see 'Details' of sensminimax.
w	Vector of the corresponding design weights for x.
lx	Vector of lower bounds for the predictors. Should be in the same order as predvars.
ux	Vector of upper bounds for the predictors. Should be in the same order as predvars.
fimfunc	A function. Returns the FIM as a matrix. Required when formula is missing. See 'Details' of minimax.
prior	An object of class cprior. User can also use one of the functions uniform, normal, skewnormal or student to create the prior. See 'Details' of bayes.
prob	Either formula or a function. When function, its argument are x and param, and they are the same as the arguments in fimfunc. prob as a function takes the design points and vector of parameters and returns the probability of success at each design points. See 'Examples'.
alpha	A value between 0 and 1. Compound or combined DP-criterion is the product of the efficiencies of a design with respect to D- and average P- optimality, weighted by alpha.
sens.control	Control Parameters for Calculating the ELB. For details, see sens.control.
sens.bayes.control	A list. Control parameters to verify the general equivalence theorem. For details, see sens.bayes.control.
crt.bayes.control	A list. Control parameters to approximate the integral in the Bayesian criterion at a given design over the parameter space. For details, see crt.bayes.control.
plot_3d	Which package should be used to plot the sensitivity (derivative) function for two-dimensional design space. Defaults to "lattice".
plot_sens	Plot the sensitivity (derivative) function? Defaults to TRUE.
npar	Number of model parameters. Used when fimfunc is given instead of formula to specify the number of model parameters. If not specified correctly, the sensitivity (derivative) plot may be shifted below the y-axis. When NULL (default), it will be set to length(parvars) or prior$npar when missing(formula).
calculate_criterion	Calculate the optimality criterion? See 'Details' of sensminimax.
silent	Do not print anything? Defaults to FALSE.

Details

Depending on the complexity of the problem at hand, sometimes, the CPU time can be considerably reduced by choosing a set of less conservative values for the tuning parameters tol and maxEval in the function sens.bayes.control when its method component is equal to "cubature". Similarly, this applies when sens.bayes.control$method = "quadrature". In general, if the CPU time matters, the user should find an appropriate speed-accuracy trade-off for her/his own problem. See 'Examples' for more details.

Note

The default values of the tuning parameters in sens.bayes.control are set in a way that having accurate plots for the sensitivity (derivative) function and calculating the ELB to a high precision to be the primary goals, although the process may take too long. The user should choose a set of less conservative values via the argument sens.bayes.control when the CPU-time is too long or matters.

See Also

bayescomp

Examples

##########################################
# Verifing the DP-optimality of a design
# The logistic model with two predictors
##########################################

# The design points and corresponding weights are as follows:
# Point1     Point2     Point3     Point4     Point5     Point6     Point7
# 0.07410  -0.31953   -1.00000     1.00000  -1.00000     1.00000    0.30193
# -1.00000  1.00000   -1.00000     1.00000   -0.08251   -1.00000    1.00000
# Weight1   Weight2    Weight3    Weight4    Weight5    Weight6    Weight7
# 0.020      0.275      0.224      0.131      0.092      0.156      0.103

# It should be given to the function as two seperate vectors:
x1 <- c(0.07409639, -0.3195265, -1, 1, -1, 1, 0.3019317, -1, 1, -1, 1, -0.08251169, -1, 1)
w1 <- c(0.01992863, 0.2745394, 0.2236575, 0.1312331, 0.09161503, 0.1561454, 0.1028811)


p <- c(1, -2, 1, -1)

## Not run: 
sensbayescomp(formula = ~exp(b0+b1*x1+b2*x2+b3*x1*x2)/(1+exp(b0+b1*x1+b2*x2+b3*x1*x2)),
              predvars = c("x1", "x2"),
              parvars = c("b0", "b1", "b2", "b3"),
              family = binomial(),
              x = x1, w = w1,
              lx = c(-1, -1), ux = c(1, 1),
              prior = uniform(p -1.5, p + 1.5),
              prob = ~1-1/(1+exp(b0 + b1 * x1 + b2 * x2 + b3 * x1 * x2)),
              alpha = .5, plot_3d = "rgl",
              sens.bayes.control = list(cubature = list(tol = 1e-3, maxEval = 1000)))

## End(Not run)

ICAOD documentation built on Oct. 23, 2020, 6:40 p.m.