augX {quantoptr}

Generate Pseudo-Responses Based on Conditional Quantile Regression Models

Description

This function supports the DR_Qopt function. For every observation, we generate pseudo-observations corresponding to treatment 0 and 1 respectively based on working conditional quantile models.

Usage

augX(raw.data, length.out = 200, txVec, moCondQuant_0, moCondQuant_1,
  nlCondQuant_0 = FALSE, nlCondQuant_1 = FALSE, start_0 = NULL,
  start_1 = NULL, clnodes)

Arguments

raw.data

A data frame, must contain all the variables that appear in moCondQuant_0 and moCondQuant_1.

length.out

an integer greater than 1. If one of the conditional quantile model is set to be nonlinear, this argument will be triggered and we will fit length.out models across quantiles equally spaced between 0.001 and 0.999. The larger this value, the more refined the performance of this method. Default is 200.

txVec

a numeric vector of observed treatment levels coded 0L and 1L.

moCondQuant_0

A formula, used to specify the formula for the conditional quantile function when treatment = 0.

moCondQuant_1

A formula, used to specify the formula for the conditional quantile function when treatment = 1.

nlCondQuant_0

logical. When nlCondQuant_0 = TRUE, it is indicated that moCondQuant_0 is nonlinear. The default value of this variable is FALSE.

nlCondQuant_1

logical. When nlCondQuant_1 = TRUE, it is indicated that moCondQuant_1 is nonlinear. The default value of this variable is FALSE.

start_0

either a list object, providing the starting value in estimating the parameters in the nonlinear conditional quantile model, given that treatment=0. Default is NULL, corresponding to the case when nlCondQuant_0=FALSE.

start_1

either a list object, providing the starting value in estimating the parameters in the nonlinear conditional quantile model, given that treatment=0. Default is NULL, corresponding to the case when nlCondQuant_1=FALSE.

clnodes

Either a cluster object to enable parallel computation or NULL. If NULL, no parallel computation will be used.

Details

This function implements the algorithm to generate individual level pseudo responses for two treatment levels respectively.

For each observation, two independent random variables from {unif}[0,1] are generated. Denote them by u_0 and u_1. Approximately, this function then estimates the u_0th quantile of this observation were treatment level 0 is applied via the conditional u_0th quantile regression. This estimated quantile will be the pseudo-response for treatment 0. Similarly, this function the pseudo-response for treatment 1 will be estimated and returned.

See the reference paper for a more formal explanation.

Value

It returns a list object, consisting of the following elements:

1. y.a.0, the vector of estimated individual level pseudo outcomes, given the treatment is 0;

2. y.a.1, the vector of estimated individual level pseudo outcomes, given the treatment is 1;

3. nlCondQuant_0, logical, indicating whether the y.a.0 is generated based on a nonlinear conditional quantile model.

4. nlCondQuant_1, logical, indicating whether the y.a.1 is generated based on a nonlinear conditional quantile model.

References

Wang L, Zhou Y, Song R and Sherwood B (2017). “Quantile-Optimal Treatment Regimes.” Journal of the American Statistical Association.

Examples

ilogit <- function(x) exp(x)/(1 + exp(x))
GenerateData.DR <- function(n)
{
  x1 <- runif(n,min=-1.5,max=1.5)
  x2 <- runif(n,min=-1.5,max=1.5)
  tp <- ilogit( 1 - 1*x1^2 - 1* x2^2)
  a <-rbinom(n,1,tp)
  y <- a * exp(0.11 - x1- x2) + x1^2 + x2^2 +  a*rgamma(n, shape=2*x1+3, scale = 1) +
       (1-a)*rnorm(n, mean = 2*x1 + 3, sd = 0.5)
  return(data.frame(x1=x1,x2=x2,a=a,y=y))
}
regimeClass = as.formula(a ~ x1+x2)
moCondQuant_0 = as.formula(y ~ x1+x2+I(x1^2)+I(x2^2))
moCondQuant_1 = as.formula(y ~ exp( 0.11 - x1 - x2)+ x1^2 + p0 + p1*x1
+ p2*x1^2 + p3*x1^3 +p4*x1^4 )
start_1 = list(p0=0, p1=1.5, p2=1, p3 =0,p4=0)

## Not run: 
n<-200
testdata <- GenerateData.DR(n)
fit1 <- augX(raw.data=testdata, txVec = testdata$a,
             moCondQuant_0=moCondQuant_0, moCondQuant_1=moCondQuant_1,
             nlCondQuant_0=FALSE,   nlCondQuant_1=TRUE,
             start_1=start_1, 
             clnodes=NULL)  
 
# How to use parallel computing in AugX(): ##
 
# on Mac OSX/linux
 clnodes <- parallel::makeForkCluster(nnodes =getOption("mc.cores",2))
 fit2 <- augX(raw.data=testdata, length.out = 5, txVec = testdata$a,
             moCondQuant_0=moCondQuant_0, moCondQuant_1=moCondQuant_1,
             nlCondQuant_0=FALSE,   nlCondQuant_1=TRUE,
             start_1=start_1, 
             clnodes=clnodes)  
  
# on Windows
 clnodes <- parallel::makeCluster(2, type="PSOCK")
 fit3 <- augX(raw.data=testdata, length.out = 5, txVec = testdata$a,
             moCondQuant_0=moCondQuant_0, moCondQuant_1=moCondQuant_1,
             nlCondQuant_0=FALSE,   nlCondQuant_1=TRUE,
             start_1=start_1, 
             clnodes=clnodes)  
 
## End(Not run)