Learn R Programming

stremr (version 0.4)

getIPWeights: Inverse Probability Weights.

Description

Evaluate the inverse probability weights for up to 3 intervention nodes: CENS, TRT and MONITOR. This is based on the inverse of the propensity score fits for the observed likelihood (g0.C, g0.A, g0.N), multiplied by the indicator of not being censored and the probability of each intervention in intervened_TRT and intervened_MONITOR. Requires column name(s) that specify the counterfactual node values or the counterfactual probabilities of each node being 1 (for stochastic interventions). The output is person-specific data with evaluated weights, wts.DT, only observation-times with non-zero weight are kept Can be one regimen per single run of this block, which are then combined into a list of output datasets with lapply. Alternative is to allow input with several rules/regimens, which are automatically combined into a list of output datasets.

Usage

getIPWeights(OData, intervened_TRT = NULL, intervened_MONITOR = NULL, useonly_t_TRT = NULL, useonly_t_MONITOR = NULL, rule_name = paste0(c(intervened_TRT, intervened_MONITOR), collapse = ""))

Arguments

OData
Input data object created by importData function.
intervened_TRT
Column name in the input data with the probabilities (or indicators) of counterfactual treatment nodes being equal to 1 at each time point. Leave the argument unspecified (NULL) when not intervening on treatment node(s).
intervened_MONITOR
Column name in the input data with probabilities (or indicators) of counterfactual monitoring nodes being equal to 1 at each time point. Leave the argument unspecified (NULL) when not intervening on the monitoring node(s).
useonly_t_TRT
Use for intervening only on some subset of observation and time-specific treatment nodes. Should be a character string with a logical expression that defines the subset of intervention observations. For example, using "TRT == 0" will intervene only at observations with the value of TRT being equal to zero. The expression can contain any variable name that was defined in the input dataset. Leave as NULL when intervening on all observations/time-points.
useonly_t_MONITOR
Same as useonly_t_TRT, but for monitoring nodes.
rule_name
Optional name for the treatment/monitoring regimen.

Value

...

Examples

Run this code
options(stremr.verbose = TRUE)
require("data.table")
set_all_stremr_options(fit.package = "speedglm", fit.algorithm = "glm")

# ----------------------------------------------------------------------
# Simulated Data
# ----------------------------------------------------------------------
data(OdataNoCENS)
OdataDT <- as.data.table(OdataNoCENS, key=c(ID, t))

# define lagged N, first value is always 1 (always monitored at the first time point):
OdataDT[, ("N.tminus1") := shift(get("N"), n = 1L, type = "lag", fill = 1L), by = ID]
OdataDT[, ("TI.tminus1") := shift(get("TI"), n = 1L, type = "lag", fill = 1L), by = ID]

# ----------------------------------------------------------------------
# Define intervention (always treated):
# ----------------------------------------------------------------------
OdataDT[, ("TI.set1") := 1L]
OdataDT[, ("TI.set0") := 0L]

# ----------------------------------------------------------------------
# Import Data
# ----------------------------------------------------------------------
OData <- importData(OdataDT, ID = "ID", t = "t", covars = c("highA1c", "lastNat1", "N.tminus1"),
                    CENS = "C", TRT = "TI", MONITOR = "N", OUTCOME = "Y.tplus1")

# ----------------------------------------------------------------------
# Model the Propensity Scores
# ----------------------------------------------------------------------
gform_CENS <- "C ~ highA1c + lastNat1"
gform_TRT = "TI ~ CVD + highA1c + N.tminus1"
gform_MONITOR <- "N ~ 1"
stratify_CENS <- list(C=c("t < 16", "t == 16"))

# ----------------------------------------------------------------------
# Fit Propensity Scores
# ----------------------------------------------------------------------
OData <- fitPropensity(OData, gform_CENS = gform_CENS,
                        gform_TRT = gform_TRT,
                        gform_MONITOR = gform_MONITOR,
                        stratify_CENS = stratify_CENS)

# ----------------------------------------------------------------------
# IPW Ajusted KM or Saturated MSM
# ----------------------------------------------------------------------
require("magrittr")
AKME.St.1 <- getIPWeights(OData, intervened_TRT = "TI.set1") %>%
             survNPMSM(OData) %$%
             IPW_estimates
AKME.St.1

# ----------------------------------------------------------------------
# Bounded IPW
# ----------------------------------------------------------------------
IPW.St.1 <- getIPWeights(OData, intervened_TRT = "TI.set1") %>%
             survDirectIPW(OData)
IPW.St.1[]

# ----------------------------------------------------------------------
# IPW-MSM for hazard
# ----------------------------------------------------------------------
wts.DT.1 <- getIPWeights(OData = OData, intervened_TRT = "TI.set1", rule_name = "TI1")
wts.DT.0 <- getIPWeights(OData = OData, intervened_TRT = "TI.set0", rule_name = "TI0")
survMSM_res <- survMSM(list(wts.DT.1, wts.DT.0), OData, t_breaks = c(1:8,12,16)-1,)
survMSM_res$St

# ----------------------------------------------------------------------
# Sequential G-COMP
# ----------------------------------------------------------------------
t.surv <- c(0:15)
Qforms <- rep.int("Q.kplus1 ~ CVD + highA1c + N + lastNat1 + TI + TI.tminus1", (max(t.surv)+1))
params = list(fit.package = "speedglm", fit.algorithm = "glm")

## Not run: 
# gcomp_est <- fitSeqGcomp(OData, t_periods = t.surv, intervened_TRT = "TI.set1",
#                           Qforms = Qforms, params_Q = params, stratifyQ_by_rule = FALSE)
# gcomp_est[]
# ## End(Not run)
# ----------------------------------------------------------------------
# TMLE
# ----------------------------------------------------------------------
## Not run: 
# tmle_est <- fitTMLE(OData, t_periods = t.surv, intervened_TRT = "TI.set1",
#                     Qforms = Qforms, params_Q = params, stratifyQ_by_rule = TRUE)
# tmle_est[]
# ## End(Not run)

# ----------------------------------------------------------------------
# Running IPW-Adjusted KM with optional user-specified weights:
# ----------------------------------------------------------------------
addedWts_DT <- OdataDT[, c("ID", "t"), with = FALSE]
addedWts_DT[, new.wts := sample.int(10, nrow(OdataDT), replace = TRUE)/10]
survNP_res_addedWts <- survNPMSM(wts.DT.1, OData, weights = addedWts_DT)

# ----------------------------------------------------------------------
# Multivariate Propensity Score Regressions
# ----------------------------------------------------------------------
gform_CENS <- "C + TI + N ~ highA1c + lastNat1"
OData <- fitPropensity(OData, gform_CENS = gform_CENS, gform_TRT = gform_TRT,
                        gform_MONITOR = gform_MONITOR)

# ----------------------------------------------------------------------
# Fitting Propensity scores with Random Forests:
# ----------------------------------------------------------------------
## Not run: 
# set_all_stremr_options(fit.package = "h2o", fit.algorithm = "randomForest")
# require("h2o")
# h2o::h2o.init(nthreads = -1)
# gform_CENS <- "C ~ highA1c + lastNat1"
# OData <- fitPropensity(OData, gform_CENS = gform_CENS,
#                         gform_TRT = gform_TRT,
#                         gform_MONITOR = gform_MONITOR,
#                         stratify_CENS = stratify_CENS)
# 
# # For Gradient Boosting machines:
# set_all_stremr_options(fit.package = "h2o", fit.algorithm = "gbm")
# # Use `H2O-3` distributed implementation of GLM
# set_all_stremr_options(fit.package = "h2o", fit.algorithm = "glm")
# # Use Deep Neural Nets:
# set_all_stremr_options(fit.package = "h2o", fit.algorithm = "deeplearning")
# ## End(Not run)

# ----------------------------------------------------------------------
# Fitting different models with different algorithms
# Fine tuning modeling with optional tuning parameters.
# ----------------------------------------------------------------------
## Not run: 
# params_TRT = list(fit.package = "h2o", fit.algorithm = "gbm", ntrees = 50,
#     learn_rate = 0.05, sample_rate = 0.8, col_sample_rate = 0.8,
#     balance_classes = TRUE)
# params_CENS = list(fit.package = "speedglm", fit.algorithm = "glm")
# params_MONITOR = list(fit.package = "speedglm", fit.algorithm = "glm")
# OData <- fitPropensity(OData,
#             gform_CENS = gform_CENS, stratify_CENS = stratify_CENS, params_CENS = params_CENS,
#             gform_TRT = gform_TRT, params_TRT = params_TRT,
#             gform_MONITOR = gform_MONITOR, params_MONITOR = params_MONITOR)
# ## End(Not run)

# ----------------------------------------------------------------------
# Running TMLE based on the previous fit of the propensity scores.
# Also applying Random Forest to estimate the sequential outcome model
# ----------------------------------------------------------------------
## Not run: 
# t.surv <- c(0:5)
# Qforms <- rep.int("Q.kplus1 ~ CVD + highA1c + N + lastNat1 + TI + TI.tminus1", (max(t.surv)+1))
# params_Q = list(fit.package = "h2o", fit.algorithm = "randomForest",
#                 ntrees = 100, learn_rate = 0.05, sample_rate = 0.8,
#                 col_sample_rate = 0.8, balance_classes = TRUE)
# tmle_est <- fitTMLE(OData, t_periods = t.surv, intervened_TRT = "TI.set1",
#             Qforms = Qforms, params_Q = params_Q,
#             stratifyQ_by_rule = TRUE)
# ## End(Not run)

## Not run: 
# t.surv <- c(0:5)
# Qforms <- rep.int("Q.kplus1 ~ CVD + highA1c + N + lastNat1 + TI + TI.tminus1", (max(t.surv)+1))
# params_Q = list(fit.package = "h2o", fit.algorithm = "randomForest",
#                 ntrees = 100, learn_rate = 0.05, sample_rate = 0.8,
#                 col_sample_rate = 0.8, balance_classes = TRUE)
# tmle_est <- fitTMLE(OData, t_periods = t.surv, intervened_TRT = "TI.set1",
#             Qforms = Qforms, params_Q = params_Q,
#             stratifyQ_by_rule = FALSE)
# ## End(Not run)

Run the code above in your browser using DataLab