est_quant_TwoStg_ipwe: Estimate the marginal quantile response of a specific dynamic...
In QTOCen: Quantile-Optimal Treatment Regimes with Censored Data

Description Usage Arguments Examples

Assume we have binary treatment options for two sequential stages with a fixed time duration between them. This means for each subject in the target population if the censored survival time or the time-to-event is beyond the timepoint of the second treatment.

This function evaluates a given dynamic treatment regime and returns the estimated marginal quantile response.

We assume the space of two-stage treatment regimes is a cartesian product of two single-stage linear treatment regime space.

The user facing function that applies this function is IPWE_Qopt_DTR_IndCen.

est_quant_TwoStg_ipwe(n, beta, sign_beta1.stg1, sign_beta1.stg2, txVec1,
  txVec2_na_omit, s_Diff_Time, nvars.stg1, nvars.stg2, p.data1, p.data2,
  censor_y, delta, ELG, w_di_vec, tau, check_complete = TRUE,
  Penalty.level = 0)

`n`	the sample size
`beta`	the vector of coefficients indexing a two-stage treatment regime
`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 \|β_1\| = 1 scale normalizaion.
`sign_beta1.stg2`	Default is NULL. Similar to `sign_beta1.stg1`.
`txVec1`	the vector of treatment received at the first stage
`txVec2_na_omit`	the vector of second stage treatment for patients who indeed second stage treatment
`s_Diff_Time`	the length of time between the first stage treatment and the second stage 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
`censor_y`	Numeric vector. The censored survival times from all observed data, i.e. `censor_y = min(Y, C)`
`delta`	Numeric vector. The censoring indicators from all observed data. We use 1 for uncensored, 0 for censored.
`ELG`	the boolean vector of whether patients get the second stage treatment
`w_di_vec`	the inverse probability weight for two stage experiments
`tau`	a value between 0 and 1. This is the quantile of interest.
`check_complete`	logical. Since this value estimation method is purely nonparametric, we need at least one unit in collected data such that the observed treatment assignment is the same what the regime parameter suggests. If `check_complete` is `TRUE`. It will check if any observation satisfies this criterion. When none observation satisfies, a message is printed to console to raise users awareness that the input regime parameter `beta` does not agree with any observed treatment level assignment. Then a sufficiently small number is returned from this function, to keep the genetic algorithm running smoothly.
`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.

##########################################################################
# Note: the preprocessing steps prior to calling est_quant_TwoStg_ipwe() #
# are wrapped up in IPWE_Qopt_DTR_IndCen().                              #
# w_di_vec is the inverse probability weight for two stage experiments   #
# We recommend users to use function IPWE_Qopt_DTR_IndCen() directly.    #
# Below is a simple customized calculation of the weight that only works #
# for this example                                                       #
##########################################################################

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"]
# get none-na second stage treatment levels in data
txVec2 <- D[, "a1"]
txVec2_na_omit <- txVec2[which(!is.na(txVec2))]

# 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
 }
}


qhat <- est_quant_TwoStg_ipwe(n=n, beta=c(2.5,2.8), 
             sign_beta1.stg1 = FALSE, sign_beta1.stg2=FALSE,
             txVec1=txVec1, txVec2_na_omit=txVec2_na_omit, s_Diff_Time=1, 
             nvars.stg1=2, nvars.stg2=2, 
             p.data1=p.data1, 
             p.data2=p.data2, 
             censor_y=D$censor_y, 
             delta=D$delta, 
             ELG=ELG, w_di_vec=D$w_di_vec, 
             tau=0.3)