Gene_Quantile_CenIPWE_DTR {QTOCen}

A low-level function for the generic optimization step in estimating dynamic Quanilte-optimal treatment regime for censored data

Description

This function supports wrapper functions for two stage Quantile-optimal treatment regime estimation, namely IPWE_Qopt_DTR_IndCen.

Usage

Gene_Quantile_CenIPWE_DTR(data, max, tau, regimeClass.stg1,
  regimeClass.stg2, s_Diff_Time, txVec1, txVec2, nvars.stg1, nvars.stg2,
  p.data1, p.data2, sign_beta1.stg1, sign_beta1.stg2, p_level, cluster,
  s.tol, it.num, pop.size, Domains1 = NULL, Domains2 = NULL,
  Penalty.level = 0)

Arguments

data

raw data.frame

max

Maximization (TRUE) or Minimizing (FALSE). Determines if genoud minimizes or maximizes the objective function.

tau

a quantile level of interest

regimeClass.stg1

the class of treatment regimes for stage one

regimeClass.stg2

the class of treatment regimes for stage two

s_Diff_Time

the length of time between the first stage treatment and the second stage treatment

txVec1

the vector of treatment received at the first stage

txVec2

the vector of treatment received at the second stage, it expects entries to be NA for patients who did not receive the second treatment

nvars.stg1

number of coeffients for the decision rule of the first stage

nvars.stg2

number of coeffients for the decision rule of the second stage

p.data1

the design matrix to be used for decision in stage one

p.data2

the design matrix to be used for decision in stage two

sign_beta1.stg1

Is sign of the coefficient for the first non-intercept variable for the first stage known? Default is NULL, meaning user does not have contraint on the sign; FALSE if the coefficient for the first continuous variable is fixed to be -1; TRUE if 1. We can make the search space discrete because we employ |\beta_1| = 1 scale normalizaion.

sign_beta1.stg2

Default is NULL. Similar to sign_beta1.stg1.

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).

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).

Domains1

This is optional. If not NULL, please provide the two-column matrix for the searching range of coeffients in stage one. The coefficient taking value of positive/negative one should not be included.

Domains2

This is optional. If not NULL, please provide the two-column matrix for the searching range of coeffients in stage two. The coefficient taking value of positive/negative one should not be included.

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.

Examples

library(survival)
# Simulate data
n=200
s_Diff_Time = 1
D <- simJLSDdata(n, case="a")

# give regime classes
regimeClass.stg1 <- as.formula(a0~x0)
regimeClass.stg2 <- as.formula(a1~x1)

# extract columns that matches each stage's treatment regime formula
p.data1 <- model.matrix(regimeClass.stg1, D)

# p.data2 would only contain observations with non-null value.
p.data2 <- model.matrix(regimeClass.stg2, D)

txVec1 <- D[, "a0"]
txVec2 <- D[, "a1"]

# Eligibility flag
ELG <- (D$censor_y  >  s_Diff_Time)

# Build weights
D$deltaC <- 1 - D$delta
survfit_all <- survfit(Surv(censor_y, event = deltaC)~1, data=D)
survest <- stepfun(survfit_all$time, c(1, survfit_all$surv))
D$ghat <- survest(D$censor_y)
g_s_Diff_Time <- survest(s_Diff_Time)
D$w_di_vec <- rep(-999, n)
for(i in 1:n){
  if (!ELG[i]) {
      D$w_di_vec[i] <- 0.5 * D$ghat[i]} else {
         D$w_di_vec[i] <- 0.5* D$ghat[i] * 0.5
 }
}


 
fit1  <- Gene_Quantile_CenIPWE_DTR(data=D, max=TRUE,
  tau=0.3,
  regimeClass.stg1 = regimeClass.stg1,
  regimeClass.stg2 = regimeClass.stg2,
  s_Diff_Time = s_Diff_Time,
  txVec1 = txVec1,
  txVec2 = txVec2,
  nvars.stg1=2,
  nvars.stg2=2,
  p.data1=p.data1,
  p.data2=p.data2,
  sign_beta1.stg1=FALSE,
  sign_beta1.stg2=NULL,
  p_level=1,
  cluster=FALSE,
  s.tol=1e-6,
  it.num=5,
  pop.size=6000,
  Domains1 = NULL,
  Domains2 = NULL,
  Penalty.level = 0
  )