Nothing
# October 26, 2018
#
#' Backwards Outcome Weighted Learning.
#'
#' Function performs a single step
#' of the bowl method. Multiple decision points can be analyzed by
#' repeated calls, as is done for qLearn() and optimalClass().
#'
#' @references Yingqi Zhao, Donglin Zeng, Eric B. Laber, Michael R. Kosorok
#' (2015)
#' New statistical learning methods for estimating optimal dynamic
#' treatment regimes.
#' Journal of the American Statistical Association,
#' 110:510, 583--598.
#'
#' @param ... Used primarily to require named input. However, inputs for
#' the optimization methods can be sent through the ellipsis. If surrogate
#' is hinge, the optimization method is dfoptim::hjk(). For all other
#' surrogates, stats::optim() is used.
#' @param moPropen An object of class modelObj or modelObjSubset, which
#' defines the model and
#' R methods to be used to obtain parameter estimates and
#' predictions for the propensity for tx.
#' See ?moPropen for details.
#' @param data A data frame of the covariates and tx histories.
#' @param reward The response vector.
#' @param txName A character object.
#' The column header of \emph{data} that corresponds to the tx covariate
#' @param regime A formula object or a list of formula objects.
#' The covariates to be included in the decision function/kernel.
#' If a list is provided,
#' this specifies that there is an underlying subset structure -- fSet must
#' then be defined. For subsets, the name of each element of the list must
#' correspond to the name of a subset. If a regime is to be estimated using
#' multiple subsets combined, each subset must be included in the name and
#' separated by a comma (no spaces).
#' @param response A numeric vector.
#' The same as reward above. Allows for naming convention followed in most
#' DynTxRegime methods.
#' @param BOWLObj NULL or \code{\link{BOWL-class}} object returned from
#' previous call to bowl(). If NULL, indicates that the function call is
#' for the first STEP of the BOWL algorithm (i.e., the final decision point).
#' If a \code{\link{BOWL-class}} object, assumed that the object was
#' returned by the preceding step of the BOWL algorithm.
#' @param lambdas A numeric object or a numeric vector object giving the
#' penalty tuning parameter(s). If more than 1 is provided,
#' the set of tuning parameter values to be considered in the
#' cross-validation algorithm (note that cvFolds must be positive in this case).
#' @param cvFolds If cross-validation is to be used to select the tuning
#' parameters and/or kernel parameters, the number of folds.
#' @param kernel A character object.
#' Must be one of \{'linear', 'poly', 'radial'\}
#' @param kparam A numeric object. \cr
#' If kernel = linear, kparam is ignored. \cr
#' If kernel = poly, kparam is the degree of the polynomial. \cr
#' If kernel = radial, kparam is the inverse bandwidth of the
#' kernel. If a vector of bandwidth parameters is given,
#' cross-validation will be used to select the parameter
#' (note that cvFolds must be positive in this case).
#' @param fSet A function or NULL defining subset structure. See ?fSet for
#' details.
#' @param surrogate The surrogate 0-1 loss function. Must be one of
#' \{'logit', 'exp', 'hinge', 'sqhinge', 'huber'\}.
#' @param verbose An integer or logical. If 0, no screen prints are generated. If 1,
#' screen prints are generated with the exception of optimization results
#' obtained in iterative algorithm. If 2, all screen prints are generated.
#'
#' @return a \code{\link{BOWL-class}} object
#'
#' @include R_class_BOWL.R
#'
#' @family statistical methods
#' @family weighted learning methods
#' @family multiple decision point methods
#'
#' @export
#'
#' @examples
#'
#' # Load and process data set
#' data(bmiData)
#'
#' # define the negative 12 month change in BMI from baseline
#' y12 <- -100*(bmiData[,6L] - bmiData[,4L])/bmiData[,4L]
#'
#' # define the negative 4 month change in BMI from baseline
#' y4 <- -100*(bmiData[,5L] - bmiData[,4L])/bmiData[,4L]
#'
#' # reward for second stage
#' rewardSS <- y12 - y4
#'
#' #### Second-stage regression
#'
#' # Constant propensity model
#' moPropen <- buildModelObj(model = ~1,
#' solver.method = 'glm',
#' solver.args = list('family'='binomial'),
#' predict.method = 'predict.glm',
#' predict.args = list(type='response'))
#'
#' fitSS <- bowl(moPropen = moPropen,
#' data = bmiData, reward = rewardSS, txName = 'A2',
#' regime = ~ parentBMI + month4BMI)
#'
#' ##Available methods
#'
#' # Coefficients of the propensity score regression
#' coef(fitSS)
#'
#' # Description of method used to obtain object
#' DTRstep(fitSS)
#'
#' # Estimated value of the optimal treatment regime for training set
#' estimator(fitSS)
#'
#' # Value object returned by propensity score regression method
#' fitObject(fitSS)
#'
#' # Summary of optimization routine
#' optimObj(fitSS)
#'
#' # Estimated optimal treatment for training data
#' optTx(fitSS)
#'
#' # Estimated optimal treatment for new data
#' optTx(fitSS, bmiData)
#'
#' # Plots if defined by propensity regression method
#' dev.new()
#' par(mfrow = c(2,4))
#'
#' plot(fitSS)
#' plot(fitSS, suppress = TRUE)
#'
#' # Value object returned by propensity score regression method
#' propen(fitSS)
#'
#' # Parameter estimates for decision function
#' regimeCoef(fitSS)
#'
#' # Show main results of method
#' show(fitSS)
#'
#' # Show summary results of method
#' summary(fitSS)
#'
#' #### First-stage regression
#'
#' # Constant propensity model
#' fitFS <- bowl(moPropen = moPropen,
#' data = bmiData, reward = y4, txName = 'A1',
#' regime = ~ gender + parentBMI,
#' BOWLObj = fitSS, lambdas = c(0.5, 1.0), cvFolds = 4L)
#'
#' ##Available methods for fitFS are as shown above for fitSS
#'
#' # Results of the cross-validation
#' cvInfo(fitFS)
#'
bowl <- function(...,
moPropen,
data,
reward,
txName,
regime,
response,
BOWLObj = NULL,
lambdas = 2.0,
cvFolds = 0L,
kernel = "linear",
kparam = NULL,
fSet = NULL,
surrogate = 'hinge',
verbose = 2L) {
# verify moPropen provided and is modelObj or modelObjSubset
if (missing(x = moPropen)) moPropen <- NULL
if (is.null(x = moPropen)) stop("moPropen must be provided")
moPropen <- .checkModelObjOrListModelObjSubset(object = moPropen,
nm = 'moPropen')
# if subset structure specified in moPropen, ensure fSet is a function
if (is(object = moPropen, class2 = "ModelObj_SubsetList")) {
if (!is.function(x = fSet)) {
stop("if subset structure in moPropen, fSet must be provided")
}
}
# data must be provided as a data.frame object.
data <- .verifyDataFrame(data = data)
# verify that a reward vector is provided.
if (missing(x = response)) response <- NULL
if (missing(x = reward)) reward <- NULL
response <- .verifyResponseReward(response = response, reward = reward)
# response must be a vector
response <- .verifyVectorResponse(response = response) + 1e-6
# verify treatment is appropriately coded.
data <- .checkTxData(txName = txName, data = data)
# regime must be formula or a list of formula
regime <- .verifyRegime(regime = regime, fSet = fSet)
# if BOWLObj not provided default to NULL
if (missing(x = BOWLObj)) BOWLObj <- NULL
# BOWLObj must be BOWL or NULL
if (!is(object = BOWLObj, class2 = "BOWL") && !is.null(x = BOWLObj)) {
stop(paste("BOWLObj must be a BOWL object from a previous",
"call or NULL if analysis is for the final",
"decision point"))
}
if (is.list(x = kernel)) {
if (!is.function(x = fSet)) {
stop("fSet must be a function when using multiple kernels")
}
if (any(is.null(x = names(x = kernel)))) {
stop("if multiple kernels, kernel must be a named list")
}
if (any(nchar(x = names(x = kernel)) == 0L)) {
stop("if multiple kernels, kernel must be a named list")
}
if (is.list(x = kparam)) {
if (!all(names(x = kparam) %in% names(x = kernel)) |
!all(names(x = kernel) %in% names(x = kparam))) {
stop("names of kernel and kparam list elements do not match")
}
}
if (!is.list(x = regime)) {
stop("a list of regimes is required for multiple kernels")
}
if (!all(names(x = regime) %in% names(x = kernel)) |
!all(names(x = kernel) %in% names(x = regime))) {
stop("names of kernel and regime list elements do not match")
}
kernelObj <- list()
cvHold <- NULL
for (i in 1L:length(x = kernel)) {
kname <- names(x = kernel)[i]
# verify cross-validation quantities
cvVerified <- .verifyCV(lambdas = lambdas,
cvFolds = cvFolds,
kparam = kparam[[ kname ]])
if (!is.null(x = cvVerified$cvFolds)) {
cvHold <- cvVerified$cvFolds
}
# define kernel
kernelObj[[ kname ]] <- .newKernelObj(kernel = tolower(kernel[[ i ]]),
data = data,
model = regime[[ kname ]],
kparam = cvVerified$kparam)@kernel
}
kernelObj <- new("SubsetList", kernelObj)
cvVerified$cvFolds <- cvVerified$cvFolds
} else {
# verify cross-validation quantities
cvVerified <- .verifyCV(lambdas = lambdas,
cvFolds = cvFolds,
kparam = kparam)
# define kernel
kernelObj <- .newKernelObj(kernel = tolower(kernel),
data = data,
model = regime,
kparam = cvVerified$kparam)@kernel
}
# fSet must be NULL or a function.
# fSet can be a function with non-list regime and modelObj moPropen when
# single tx options exist for some patients.
if (!is.function(x = fSet) && !is.null(x = fSet)) {
stop("fSet must be NULL or a function")
}
# define surrogate
surrogate <- .verifySurrogate(surrogate = surrogate)
# verbose must 0, 1, 2
if (is.logical(x = verbose)) verbose <- verbose*1L
verbose <- as.integer(x = round(x = verbose))
if (verbose > 2L) verbose <- 2L
if (verbose < 0L) verbose <- 0L
obj <- .newBOWL(BOWLObj = BOWLObj,
moPropen = moPropen,
fSet = fSet,
data = data,
response = response,
txName = txName,
lambdas = cvVerified$lambdas,
cvFolds = cvVerified$cvFolds,
kernel = kernelObj,
surrogate = surrogate,
suppress = verbose,
guess = NULL, ...)
obj@analysis@call <- match.call()
return( obj )
}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.