R/estimateQ1.ltmle.R

In benkeser/drinf: Computes doubly-robust estimators and doubly-robust confidence intervals

#' estimateQ1.ltmle
#' 
#' This function computes the conditional treatment probabilities at time 1
#' 
#' @param L0 A \code{data.frame} featuring covariates measured at baseline.
#' @param L1 A \code{data.frame} featuring time-varying covariates measured at 
#' the first timepoint.
#' @param L2 A \code{vector} outcome of interest
#' @param A0 A \code{vector} treatment delivered at baseline.
#' @param A1 A \code{vector} treatment deliver after \code{L1} is measured.
#' @param folds Vector of cross-validation folds
#' @param validFold Which fold is validation fold 
#' @param abar A \code{vector} of length 2 indicating the treatment assignment 
#' that is of interest. 
#' @param stratify A \code{boolean} indicating whether to pool across treatment
#' nodes or to estimate outcome regression separately in each category.
#' @param SL.Q A \code{vector} or \code{list} specifying the SuperLearner library
#' to be used to estimate the conditional probability of treatment at each time point.  
#' See \code{SuperLearner} package for details.
#' @param glm.Q A \code{character} specifying the right-hand side of the \code{glm} 
#' formula used to estimate the conditional probability of treatment at each time point. 
#' Only used if \code{SL.Q = NULL}.
#' @param SL.Q.options A \code{list} of additional arguments passed to \code{SuperLearner} 
#' for outcome regression fits.
#' @param glm.Q.options A \code{list} of additional arguments passed to \code{glm} for the
#' outcome regression fits. Typically, the \code{family} argument.
#' @param return.models  A \code{boolean} indicating whether the models for g00 should be 
#' returned with the output. 
#' @param ... Other arguments (not currently used).
#' 
#' @return Returns a list with \code{Q2n}, \code{Q1n}, and the estimated model objects if
#' \code{return.models = TRUE}
#' 
#' @importFrom SuperLearner SuperLearner
#' 
#' @export
#' 
#' @examples TO DO : add examples


estimateQ1.ltmle <- function(
    L2, Q2n, L0, A0, A1, folds, validFold, abar, SL.Q, SL.Q.options, glm.Q, 
    glm.Q.options, return.models, verbose, stratify, 
    ...
){
    all1 <- all(folds == 1)
    if(all1){
        train_L0 <- valid_L0 <- L0
        train_A0 <- valid_A0 <- A0
        train_A1 <- valid_A1 <- A1
        train_L2 <- valid_L2 <- L2
        train_Q2n <- valid_Q2n <- Q2n
    }else{
        # training data
        train_L0 <- L0[folds != validFold, , drop = FALSE]
        train_Q2n <- Q2n[folds != validFold]
        train_A0 <- A0[folds != validFold]
        train_A1 <- A1[folds != validFold]
        train_L2 <- L2[folds != validFold]

        # validation data
        valid_L0 <- L0[folds == validFold, , drop = FALSE]
        valid_Q2n <- Q2n[folds == validFold]
        valid_A0 <- A0[folds == validFold]
        valid_A1 <- A1[folds == validFold]
        valid_L2 <- L2[folds == validFold]
    }
    if(is.null(SL.Q) & is.null(glm.Q)){ 
        stop("Specify Super Learner library or GLM formula for g")
    }
    if(!is.null(SL.Q) & !is.null(glm.Q)){
        warning("Super Learner library and GLM formula specified. Proceeding with Super Learner only.")
        glm.Q <- NULL
    }
    # Super Learner
    if(!is.null(SL.Q)){
        # if multiple learners specified, call SuperLearner
        multiAlgos <- (length(SL.Q) > 1 | is.list(SL.Q))
       #-----------
        # Q1n
        #-----------
        Q1mod <- do.call(ifelse(multiAlgos,getFromNamespace("SuperLearner","SuperLearner"),SL.Q),args = c(list(
            Y=eval(parse(text=paste0(
                ifelse(stratify,
                       "train_Q2n[train_A0==abar[1]]",
                       "train_Q2n"
                )
            ))), 
            X=eval(parse(text=paste0(
                ifelse(stratify,
                       "train_L0[train_A0==abar[1],,drop=FALSE]",
                       "cbind(train_L0,train_A0)"
                )
            ))), 
            newX = eval(parse(text=paste0(
                ifelse(stratify,
                       "valid_L0",
                       "cbind(valid_L0,A0=abar[1])"
                )
            ))),
            SL.library=SL.Q,
            obsWeights = rep(1, ifelse(stratify, length(train_Q2n[train_A0==abar[1]]),length(train_Q2n))),
            verbose=verbose), SL.Q.options
        ))
        if(multiAlgos){
            # Super Learner predictions
            Q1n <- Q1mod$SL.predict
        }else{
            Q1n <- Q1mod$pred
        }
        # replace extrapolated predictions with 
        # smallest/largest value
        Q1n[Q1n < min(Q2n)] <- min(Q2n)
        Q1n[Q1n > max(Q2n)] <- max(Q2n)
    } # end Super Learner g call
    
    # if SL.Q == NULL then call glm
    if(!is.null(glm.Q)){
        #-------
        # Q1n
        #-------
        Q1mod <- do.call("glm", c(list(formula = as.formula(paste0(
            ifelse(stratify,
                   "train_Q2n[train_A0==abar[1]] ~",
                   "train_Q2n ~"
            ),glm.Q)),
            data = eval(parse(text=paste0(
                ifelse(stratify,
                       "train_L0[train_A0==abar[1],,drop=FALSE]",
                       "cbind(train_L0,train_A0)"
                ))))), glm.Q.options))
        
        Q1n <- predict(Q1mod, type="response", 
                       newdata = eval(parse(text=paste0(
                           ifelse(stratify,
                                  "valid_L0",
                                  "data.frame(valid_L0, train_A0 = abar[1])"
                           ))))
        )
        
        # replace extrapolated predictions with 
        # smallest/largest value
        Q1n[Q1n < min(Q2n)] <- min(Q2n)
        Q1n[Q1n > max(Q2n)] <- max(Q2n)
    } # end glm call
    
    #-----------
    # returning
    #-----------
    # don't like the output when do.call is used in call
    Q1mod$call <- NULL
    
    out <- list(Q1n = Q1n, Q1mod = NULL)
    if(return.models){
        out$Q1mod <- Q1mod
    }
    return(out)
} # end function