Chapter 4 The TMLE Framework
Jeremy Coyle
Based on the tmle3 R package.
Updated: 2019-05-22
4.1 Learning Objectives
- Use
tmle3to estimate an Average Treatment Effect (ATE) - Understand
tmle3“Specs” - Fit
tmle3for a custom set of parameters - Use the delta method to estimate transformations of parameters
4.2 Easy-Bake Example: tmle3 for ATE
We’ll illustrate the most basic use of TMLE using the WASH Benefits data introduced earlier and estimating an Average Treatment Effect (ATE).
As a reminder, the ATE is identified with the following statistical parameter (under assumptions): \(ATE = \mathbb{E}_0(Y(1)-Y(0)) = \mathbb{E}_0\left(\mathbb{E}_0[Y \mid A=1,W]-\mathbb{E}_0[Y \mid A=0,W] \right)\)
This Easy-Bake implementation consists of the following steps:
- Load the necessary libraries and data
- Define the variable roles
- Create a “Spec” object
- Define the super learners
- Fit the TMLE
- Evaluate the TMLE estimates
0. Load the Data
We’ll use the same WASH Benefits data as the earlier chapters:
library(data.table)
library(tmle3)
library(sl3)
washb_data <- fread("https://raw.githubusercontent.com/tlverse/tlverse-data/master/wash-benefits/washb_data.csv", stringsAsFactors = TRUE)1. Define the variable roles
We’ll use the common \(W\) (covariates), \(A\) (treatment/intervention), \(Y\)
(outcome) data structure. tmle3 needs to know what variables in the dataset
correspond to each of these roles. We use a list of character vectors to tell
it. We call this a “Node List” as it corresponds to the nodes in a Directed
Acyclic Graph (DAG), a way of displaying causal relationships between variables.
node_list <- list(
W = c(
"month", "aged", "sex", "momage", "momedu",
"momheight", "hfiacat", "Nlt18", "Ncomp", "watmin",
"elec", "floor", "walls", "roof", "asset_wardrobe",
"asset_table", "asset_chair", "asset_khat",
"asset_chouki", "asset_tv", "asset_refrig",
"asset_bike", "asset_moto", "asset_sewmach",
"asset_mobile"
),
A = "tr",
Y = "whz"
)Handling Missingness
Currently, missingness in tmle3 is handled in a fairly simple way:
- Missing covariates are median (for continuous) or mode (for discrete) imputed, and additional covariates indicating imputation are generated
- Observations missing either treatment or outcome variables are excluded.
We plan to implement IPCW-TMLE to more efficiently handle missingness in the treatment and outcome variables.
These steps are implemented in the process_missing function in tmle3:
processed <- process_missing(washb_data, node_list)
washb_data <- processed$data
node_list <- processed$node_list2. Create a “Spec” Object
tmle3 is general, and allows most components of the TMLE procedure to be
specified in a modular way. However, most end-users will not be interested in
manually specifying all of these components. Therefore, tmle3 implements a
tmle3_Spec object that bundles a set of components into a specification
that, with minimal additional detail, can be run by an end-user.
We’ll start with using one of the specs, and then work our way down into the
internals of tmle3.
ate_spec <- tmle_ATE(
treatment_level = "Nutrition + WSH",
control_level = "Control"
)3. Define the Relevant Super Learners
Currently, the only other thing a user must define are the sl3 learners used
to estimate the relevant factors of the likelihood: Q and g.
This takes the form of a list of sl3 learners, one for each likelihood factor
to be estimated with sl3:
# choose base learners
lrnr_mean <- make_learner(Lrnr_mean)
lrnr_xgboost <- make_learner(Lrnr_xgboost)
# define metalearners appropriate to data types
ls_metalearner <- make_learner(Lrnr_nnls)
mn_metalearner <- make_learner(
Lrnr_solnp, metalearner_linear_multinomial,
loss_loglik_multinomial
)
sl_Y <- Lrnr_sl$new(
learners = list(lrnr_mean, lrnr_xgboost),
metalearner = ls_metalearner
)
sl_A <- Lrnr_sl$new(
learners = list(lrnr_mean, lrnr_xgboost),
metalearner = mn_metalearner
)
learner_list <- list(A = sl_A, Y = sl_Y)Here, we use a Super Learner as defined in the previous sl3 section. In the
future, we plan to include reasonable default learners.
4. Fit the TMLE
We now have everything we need to fit the tmle using tmle3:
tmle_fit <- tmle3(ate_spec, washb_data, node_list, learner_list)5. Evaluate the Estimates
We can see the summary results by printing the fit object. Alternatively, we can extra results from the summary by indexing into it:
print(tmle_fit)A tmle3_Fit that took 1 step(s)
type param init_est tmle_est
1: ATE ATE[Y_{A=Nutrition + WSH}-Y_{A=Control}] 0.002335129 0.00714325
se lower upper psi_transformed lower_transformed
1: 0.05044591 -0.09172891 0.1060154 0.00714325 -0.09172891
upper_transformed
1: 0.1060154estimates <- tmle_fit$summary$psi_transformed
print(estimates)[1] 0.007143254.3 tmle3 Components
Now that we’ve successfully used a spec to obtain a TML estimate, let’s look under the hood at the components. The spec has a number of functions that generate the objects necessary to define and fit a TMLE.
4.3.1 tmle3_task
First is, a tmle3_Task, analogous to an sl3_Task, containing the data we’re
fitting the TMLE to, as well as an NPSEM generated from the node_list defined
above, describing the variables and their relationships.
tmle_task <- ate_spec$make_tmle_task(washb_data, node_list)tmle_task$npsem$W
tmle3_Node: W
Variables: month, aged, sex, momedu, hfiacat, Nlt18, Ncomp, watmin, elec, floor, walls, roof, asset_wardrobe, asset_table, asset_chair, asset_khat, asset_chouki, asset_tv, asset_refrig, asset_bike, asset_moto, asset_sewmach, asset_mobile, momage, momheight, delta_momage, delta_momheight
Parents:
$A
tmle3_Node: A
Variables: tr
Parents: W
$Y
tmle3_Node: Y
Variables: whz
Parents: A, W4.3.2 Initial Likelihood
Next, is an object representing the likelihood, factorized according to the NPSEM described above:
initial_likelihood <- ate_spec$make_initial_likelihood(
tmle_task,
learner_list
)
print(initial_likelihood)W: Lf_emp
A: LF_fit
Y: LF_fitThese components of the likelihood indicate how the factors were estimated: the
marginal distribution of \(W\) was estimated using NP-MLE, and the conditional
distributions of \(A\) and \(Y\) were estimated using sl3 fits (as defined with
the learner_list) above.
We can use this in tandem with the tmle_task object to obtain likelihood
estimates for each observation:
initial_likelihood$get_likelihoods(tmle_task) W A Y
1: 0.0002129925 0.2477793 -0.6650462
2: 0.0002129925 0.2547205 -0.6374359
3: 0.0002129925 0.2592492 -0.6249582
4: 0.0002129925 0.2805867 -0.6042583
5: 0.0002129925 0.2536626 -0.5464700
---
4691: 0.0002129925 0.1350039 -0.4648466
4692: 0.0002129925 0.1261634 -0.4840708
4693: 0.0002129925 0.1264105 -0.5704188
4694: 0.0002129925 0.1758334 -0.8245943
4695: 0.0002129925 0.1299581 -0.54349674.3.3 Targeted Likelihood (updater)
We also need to define a “Targeted Likelihood” object. This is a special type
of likelihood that is able to be updated using an tmle3_Update object. This
object defines the update strategy (e.g. submodel, loss function, CV-TMLE or
not, etc).
targeted_likelihood <- Targeted_Likelihood$new(initial_likelihood)When constructing the targeted likelihood, you can specify different update
options. See the documentation for tmle3_Update for details of the different
options. For example, you can disable CV-TMLE (the default in tmle3) as
follows:
targeted_likelihood_no_cv <-
Targeted_Likelihood$new(initial_likelihood,
updater = list(cvtmle = FALSE)
)4.3.4 Parameter Mapping
Finally, we need to define the parameters of interest. Here, the spec defines a single parameter, the ATE. In the next section, we’ll see how to add additional parameters.
tmle_params <- ate_spec$make_params(tmle_task, targeted_likelihood)
print(tmle_params)[[1]]
Param_ATE: ATE[Y_{A=Nutrition + WSH}-Y_{A=Control}]4.3.5 Putting it all together
Having used the spec to manually generate all these components, we can now
manually fit a tmle3:
tmle_fit_manual <- fit_tmle3(
tmle_task, targeted_likelihood, tmle_params,
targeted_likelihood$updater
)
print(tmle_fit_manual)A tmle3_Fit that took 1 step(s)
type param init_est tmle_est
1: ATE ATE[Y_{A=Nutrition + WSH}-Y_{A=Control}] 0.00247001 0.007240568
se lower upper psi_transformed lower_transformed
1: 0.05045259 -0.09164469 0.1061258 0.007240568 -0.09164469
upper_transformed
1: 0.1061258The result is equivalent to fitting using the tmle3 function as above.
4.4 Fitting tmle3 with multiple parameters
Above, we fit a tmle3 with just one parameter. tmle3 also supports fitting
multiple parameters simultaneously. To illustrate this, we’ll use the
tmle_TSM_all spec:
tsm_spec <- tmle_TSM_all()
targeted_likelihood <- Targeted_Likelihood$new(initial_likelihood)
all_tsm_params <- tsm_spec$make_params(tmle_task, targeted_likelihood)
print(all_tsm_params)[[1]]
Param_TSM: E[Y_{A=Control}]
[[2]]
Param_TSM: E[Y_{A=Handwashing}]
[[3]]
Param_TSM: E[Y_{A=Nutrition}]
[[4]]
Param_TSM: E[Y_{A=Nutrition + WSH}]
[[5]]
Param_TSM: E[Y_{A=Sanitation}]
[[6]]
Param_TSM: E[Y_{A=WSH}]
[[7]]
Param_TSM: E[Y_{A=Water}]This spec generates a Treatment Specific Mean (TSM) for each level of the exposure variable. Note that we must first generate a new targeted likelihood, as the old one was targeted to the ATE. However, we can recycle the initial likelihood we fit above, saving us a super learner step.
4.4.1 Delta Method
We can also define parameters based on Delta Method Transformations of other parameters. For instance, we can estimate a ATE using the delta method and two of the above TSM parameters:
ate_param <- define_param(
Param_delta, targeted_likelihood,
delta_param_ATE,
list(all_tsm_params[[1]], all_tsm_params[[4]])
)
print(ate_param)Param_delta: E[Y_{A=Nutrition + WSH}] - E[Y_{A=Control}]This can similarly be used to estimate other derived parameters like Relative Risks, and Population Attributable Risks
4.4.2 Fit
We can now fit a TMLE simultaneously for all TSM parameters, as well as the above defined ATE parameter
all_params <- c(all_tsm_params, ate_param)
tmle_fit_multiparam <- fit_tmle3(
tmle_task, targeted_likelihood, all_params,
targeted_likelihood$updater
)
print(tmle_fit_multiparam)A tmle3_Fit that took 1 step(s)
type param init_est tmle_est
1: TSM E[Y_{A=Control}] -0.59814753 -0.623969131
2: TSM E[Y_{A=Handwashing}] -0.61084934 -0.641321765
3: TSM E[Y_{A=Nutrition}] -0.60631283 -0.617745710
4: TSM E[Y_{A=Nutrition + WSH}] -0.59567752 -0.616611780
5: TSM E[Y_{A=Sanitation}] -0.59145033 -0.590441751
6: TSM E[Y_{A=WSH}] -0.53190979 -0.447828018
7: TSM E[Y_{A=Water}] -0.57962980 -0.536874708
8: ATE E[Y_{A=Nutrition + WSH}] - E[Y_{A=Control}] 0.00247001 0.007357351
se lower upper psi_transformed lower_transformed
1: 0.02970541 -0.68219066 -0.5657476 -0.623969131 -0.68219066
2: 0.04220838 -0.72404866 -0.5585949 -0.641321765 -0.72404866
3: 0.04240669 -0.70086129 -0.5346301 -0.617745710 -0.70086129
4: 0.04093170 -0.69683644 -0.5363871 -0.616611780 -0.69683644
5: 0.04251325 -0.67376618 -0.5071173 -0.590441751 -0.67376618
6: 0.04511973 -0.53626106 -0.3593950 -0.447828018 -0.53626106
7: 0.03933466 -0.61396923 -0.4597802 -0.536874708 -0.61396923
8: 0.05044396 -0.09151099 0.1062257 0.007357351 -0.09151099
upper_transformed
1: -0.5657476
2: -0.5585949
3: -0.5346301
4: -0.5363871
5: -0.5071173
6: -0.3593950
7: -0.4597802
8: 0.10622574.5 Exercise
Follow the steps below to estimate an average treatment effect using data from
the Collaborative Perinatal Project (CPP), available in the sl3 package. To simplify this example, we define a binary intervention variable, parity01
– an indicator of having one or more children before the current child and a
binary outcome, haz01 – an indicator of having an above average height for
age.
Work with a buddy/team. You have 20 minutes.
In the etherpad, submit your group’s answers to the following:
- Interpret the
tmle3fit both causally and statistically. - Did your group face any challenges?
- Any additional comments/questions about this
tmle3section of the workshop?
# load the data set
data(cpp)
cpp <- cpp[!is.na(cpp[, "haz"]), ]
cpp$parity01 <- as.numeric(cpp$parity > 0)
cpp[is.na(cpp)] <- 0
cpp$haz01 <- as.numeric(cpp$haz > 0)- Define the variable roles \((W,A,Y)\) by creating a list of these nodes.
Include the following baseline covariates in \(W\):
apgar1,apgar5,gagebrth,mage,meducyrs,sexn. Both \(A\) and \(Y\) are specified above. - Define a
tmle3_Specobject for the ATE,tmle_ATE(). - Using the same base learning libraries defined above, specify
sl3base learners for estimation of \(Q = E(Y|A,Y)\) and \(g=P(A|W)\). - Define the metalearner like below
metalearner <- make_learner(Lrnr_solnp,
loss_function = loss_loglik_binomial,
learner_function = metalearner_logistic_binomial)- Define one super learner for estimating \(Q\) and another for estimating \(g\). Use the metalearner above for both \(Q\) and \(g\) super learners.
- Create a list of the two super learners defined in Step 5 and call this
object
learner_list. The list names should beA(defining the super learner for estimating \(g\)) andY(defining the super learner for estimating \(Q\)). - Fit the tmle with the
tmle3function by specifying (1) thetmle3_Spec, which we defined in Step 2; (2) the data; (3) the list of nodes, which we specified in Step 1; and (4) the list of super learners for estimating \(g\) and \(Q\), which we defined in Step 6. Note: Like before, you will need to make a data copy to deal withdata.tableweirdness (cpp2 <- data.table::copy(cpp)) and usecpp2as the data.
4.6 Summary
tmle3 is a general purpose framework for generating TML estimates. The
easiest way to use it is to use a predefined spec, allowing you to just fill in
the blanks for the data, variable roles, and sl3 learners. However, digging
under the hood allows users to specify a wide range of TMLEs. In the next
sections, we’ll see how this framework can be used to estimate advanced
parameters such as optimal treatments and shift interventions.