Gene_Quantile_CenIPWE: A low-level function for the generic optimization step in...
In QTOCen: Quantile-Optimal Treatment Regimes with Censored Data

This function supports several user facing functions for Quantile-optimal treatment regime estimation, namely

IPWE_Qopt_IndCen(), IPWE_Qopt_DTR_IndCen(), IPWE_Qopt_DepCen_trt(), and IPWE_Qopt_DepCen_general().

It implements the genetic algorithm based policy-search method with inverse probability weighting for censored data, such that the estimator is cube root consistent under the assumption that the propensity score model and the model for the survival distriution of the censoring time variable are both correct.

1
2
3

Gene_Quantile_CenIPWE(data_aug, tau, p_level, regimeClass,
  cluster = FALSE, s.tol = 1e-04, it.num = 8, pop.size = 5000,
  Domains = NULL, sign_beta1 = NULL, Penalty.level = 0)

`data_aug`	a data.frame of the observed data after preprocessing. It should include be augmented with a new column: `epsi` for the composite weights.
`tau`	a value between 0 and 1. This is the quantile of interest.
`p_level`	choose between 0,1,2,3 to indicate different levels of output from the genetic function. Specifically, 0 (minimal printing), 1 (normal), 2 (detailed), and 3 (debug).
`regimeClass`	a formula specifying the class of treatment regimes to search, e.g. if `regimeClass = a~x1+x2`, and then this function will search the class of treatment regimes of the form d(x)=I(β_0 +β_1 x1 + β_2 * x2 > 0).* Polynomial arguments are also supported.
`cluster`	default is FALSE, meaning do not use parallel computing for the genetic algorithm(GA).
`s.tol`	tolerance level for the GA algorithm. This is input for parameter `solution.tolerance` in function `rgenoud::genoud`.
`it.num`	the maximum GA iteration number
`pop.size`	an integer with the default set to be 3000. This is roughly the number individuals for the first generation in the genetic algorithm (`rgenoud::genoud`).
`Domains`	default is NULL. Otherwise, the object should be a `nvars *2` matrix used as the space of parameters, which will be supplied to `rgenoud::genoud`.
`sign_beta1`	logical. FALSE if the coefficient for the first continuous variable is fixed to be negative one; TRUE if positive one.
`Penalty.level`	the level that determines which objective function to use. `Penalty.level = 0` indicates no regularization; `Penalty.level = 1` indicates the value function estimation minus the means absolute average coefficient is the output, which is useful trick to achieve uniqueness of estimated optimal TR when resolution of input response is low.

GenerateData <- function(n)
{
  x1 <- runif(n, min=-0.5,max=0.5)
  x2 <- runif(n, min=-0.5,max=0.5)
  error <- rnorm(n, sd= 1)
  ph <- rep(0.5,n)
  a <- rbinom(n = n, size = 1, prob=ph)
  c <- 1.5 +  + runif(n = n, min=0, max=2)
  cmplt_y <-  pmin(2+x1+x2 +  a*(1 - x1 - x2) +  (0.2 + a*(1+x1+x2)) * error, 4.4)
  censor_y <- pmin(cmplt_y, c)
  delta <- as.numeric(c > cmplt_y)
  return(data.frame(x1=x1,x2=x2,a=a, censor_y = censor_y, delta=delta))
}
n <- 100
data <- GenerateData(n)

# preprocessing
data_aug <- data
data_aug$ph <- rep(mean(data$a), n)
data_aug$deltaC <- 1 - data_aug$delta
library(survival)
survfit_all <- survfit(Surv(censor_y, event = deltaC)~1, data=data_aug)
survest <- stepfun(survfit_all$time, c(1, survfit_all$surv))
data_aug$ghat <- survest(data_aug$censor_y)
data_aug$epsi <- (data_aug$ph * data_aug$a + (1 - data_aug$ph) * (1 - data_aug$a)) * data_aug$ghat

# estimate the median-optimal treatment regime

quantopt_fit <- Gene_Quantile_CenIPWE(data_aug=data_aug,tau=0.5,
                                      p_level=1, regimeClass=a~x1+x2^2, 
                                      sign_beta1=FALSE)