Update formula

The estimate and predict functions of learner class objects are by design immutable. Both functions can be inspected as part of the return values of the summary method

lr_xgboost$summary()$estimate
#> function (x, y, ...) 
#> {
#>     d <- xgboost::xgb.DMatrix(x, label = y)
#>     res <- do.call(xgboost::xgboost, c(list(data = d), list(...)), 
#>         )
#>     return(res)
#> }
#> <bytecode: 0x9c00e3508>
#> <environment: 0x9c0804938>

Rare situations may arise where one wants to update the formula argument. This supported and implemented via the update method

lr_xgboost$update(y ~ age + sex)

The design matrix that results from by the specified formula can be inspected via

head(lr_xgboost$design(pbc)$x)
#>        age sexf
#> 1 58.76523    1
#> 2 56.44627    1
#> 3 70.07255    0
#> 4 54.74059    1
#> 5 38.10541    1
#> 6 66.25873    1

and the response vector via

head(lr_xgboost$response(pbc))
#> 1 2 3 4 5 6 
#> 1 0 0 0 0 0

See ?learner for the differences between learner$design and learner$response and ?design for more details about the construction of the design matrix and response vector.

Cross-validation

The {targeted} package provides a generic implementation for repeated k-fold cross-validation with learner class objects. Parallelization is supported via the {future} and {parallel} packages (see ?cv for more details).

# future::plan("multicore")
lrs <- list(
  glm = learner_glm(y ~ age + age, family = "binomial"),
  gam = learner_gam(y ~ s(age) + s(bili), family = "binomial")
)
 # 2 times repeated 5-fold cross-validation
cv(lrs, data = pbc, rep = 2, nfolds = 5)
#> Call: cv.default(object = lrs, data = pbc, nfolds = 5, rep = 2)
#> 
#> 5-fold cross-validation with 2 repetitions
#> 
#>            mse       mae
#> glm 0.10736745 0.2127521
#> gam 0.09894063 0.1888139

Fit function with formula and data arguments

The first general case covers fit functions which expect a formula and data argument. Both arguments are used then by the fitting routine to construct the design matrix and response vector. Statistical models are usually implement with such an interface, with examples including stats::glm, mgcv::gam and earth::earth. The learner R6 class supports these fitting routines by checking if the provided estimate function has a formula and data argument. If it does, then the formula and data argument are passed on to the estimate function without further modifications. This is exemplified in the following for stats::glm, where we define a new constructor to return a learner class object that fits a generalized model

new_glm <- function(formula, ...) {
  learner$new(
    formula = formula,
    estimate = stats::glm,
    predict = stats::predict,
    predict.args = list(type = "response"),
    estimate.args = list(...),
    info = "new glm learner" # optional
  )
}
lr <- new_glm(y ~ age, family = "binomial")
lr
#> ────────── learner object ──────────
#> new glm learner 
#> 
#> Estimate arguments: family=binomial 
#> Predict arguments: type=response 
#> Formula: y ~ age

It can be seen that the optional arguments of new_glm define the estimate.args of a constructed learner object. When estimating the model with

lr$estimate(pbc)

the parameters defined in estimate.args are passed on with the formula and data to stats::glm (i.e. the function specified via the estimate argument). Thus, the above is equivalent to

fit <- glm(y ~ age, family = "binomial", data = pbc)
all(coef(fit) == coef(lr$fit))
#> [1] TRUE

The code further instructs to construct a learner object that uses stats::predict to make predictions. The learner object similarly passes on the predict.args to stats::predict for

lr$predict(head(pbc))
#>          1          2          3          4          5          6 
#> 0.16171479 0.14624529 0.25614862 0.13566571 0.06272415 0.22071201

which is equivalent to

predict(fit, newdata = head(pbc), type = "response")
#>          1          2          3          4          5          6 
#> 0.16171479 0.14624529 0.25614862 0.13566571 0.06272415 0.22071201

Indeed, the documentation of ?learner reveals that the predict method always requires an object and newdata argument. An important implementation detail is that the learner R6 class allows to overrule the defined estimate.args and predict.args in the method calls.

lr$estimate(pbc, family = "poisson")
lr$predict(head(pbc), type = "link")
#>         1         2         3         4         5         6 
#> -1.837279 -1.939001 -1.341276 -2.013822 -2.743535 -1.508572

Fit function with x and y arguments

Most fitting routines for machine learning models expect a design matrix and response vector as inputs. The R6 learner class supports these fitting functions by internally processing the formula argument via the targeted::design function. Take the example of grf::probability_forest, which expects a X (design matrix) and Y (response vector) as inputs.

new_grf <- function(formula, ...) {
  learner$new(
    formula = formula,
    estimate = function(x, y, ...) grf::probability_forest(X = x, Y = y, ...),
    predict = function(object, newdata) {
      predict(object, newdata = newdata)$predictions
    },
    estimate.args = list(...),
    info = "grf::probability_forest"
  )
}
lr <- new_grf(as.factor(y) ~ age + bili, num.trees = 100)
lr$estimate(pbc)

Compared to the previous case, the pbc data object is not directly passed on to the fit function. Instead,targeted::design constructs the design matrix and response vector from the defined formula argument. As shown previously,

dsgn <- lr$design(pbc)

can be used to inspect the design object. The x and y attributes of the returned design object are then passed on to the fit function.

Fit function with single data argument

To support ensemble/meta learners, it is also possible to construct a learner without providing a formula argument.

new_sl <- function(learners, ...) {
  learner$new(
    info = "new superlearner",
    estimate = superlearner,
    predict = targeted:::predict.superlearner,
    estimate.args = c(list(learners = learners), list(...))
  )
}
lrs <- list(
  glm = learner_glm(y ~ age, family = "binomial"),
  gam = learner_gam(y ~ s(age), family = "binomial")
)
lr <- new_sl(lrs, nfolds = 2)
lr$estimate(pbc)
lr
#> ────────── learner object ──────────
#> new superlearner 
#> 
#> Estimate arguments: learners=<list>, nfolds=2 
#> Predict arguments:   
#> Formula: NULL 
#> ─────────────────────────────────────
#>         score    weight
#> glm 0.1081636 0.6177615
#> gam 0.1082369 0.3822385

In this case, all arguments provided to lr$estimate are joined together with the specified estimate.args and passed on to the defined estimate function (i.e. superlearner).

Specials in formula

Certain learner functions allow for specifying special terms in the formula. For example, let’s consider an aggregated dataset that includes only the response variable, treatment, and sex:

library("data.table")
dd <- data.table(pbc)[!is.na(trt), .(.N), by=.(y,trt,sex)]
print(dd)
#>        y   trt    sex     N
#>    <int> <int> <fctr> <int>
#> 1:     1     1      f    13
#> 2:     0     1      f   124
#> 3:     0     1      m    19
#> 4:     0     2      f   123
#> 5:     1     2      f    16
#> 6:     0     2      m    12
#> 7:     1     2      m     3
#> 8:     1     1      m     2

Next, we fit a Naive Bayes classifier using the weights special term to specify the frequency weights for the estimation

lr <- learner_naivebayes(y ~ trt + sex + weights(N))
lr$estimate(dd)
lr$predict(dd)
#> [1] 0.09138473 0.09138473 0.12139346 0.11913520 0.11913520 0.15668515 0.15668515
#> [8] 0.12139346

Here weights.numeric is simply the identity function

targeted:::weights.numeric
#> function (object, ...) 
#> object
#> <bytecode: 0x9c1e5d380>
#> <environment: namespace:targeted>

To illustrate how to define a custom learner that utilizes special terms, consider the following example where we introduce a strata term. In this case, we introduce an estimation method that fits a linear regression model for each value of the strata variable, as well as a corresponding prediction method

est <- function(formula, data, strata, ...)
  lapply(levels(strata), \(s) lm(formula, data[which(strata==s),]))

pred <- function(object, newdata, strata, ...) {
  res <- numeric(length(strata))
  for (i in seq_along(levels(strata))) {
    idx <- which(strata == levels(strata)[i])
    res[idx] <- predict(object[[i]], newdata[idx, ], ...)
  }
  return(res)
}

The new learner is now defined by including the argument `specials = “strata”, which ensures that the strata variable is correctly passed to both the estimate and predict functions

lr <- learner$new(y ~ sex + strata(trt),
                  estimate=est, predict=pred, specials = "strata")
des <- lr$design(head(pbc))
des
#> ────────── design object ──────────
#> 
#> response (length: 6)     
#> 1   1
#> 2   0
#> ---  
#> 5   0
#> 6   0
#> 
#> specials
#>  - strata [factor]
#> 
#> design matrix (dim: 6, 1)
#>     sexf
#> 1   1   
#> 2   1   
#> ---     
#> 5   1   
#> 6   1
des$strata
#>     1     2     3     4     5     6 
#> trt=1 trt=1 trt=1 trt=1 trt=2 trt=2 
#> Levels: trt=1 trt=2

lr$estimate(pbc)
lr
#> ────────── learner object ──────────
#> Estimate arguments:   
#> Predict arguments:   
#> Formula: y ~ sex + strata(trt) 
#> ─────────────────────────────────────
#> [[1]]
#> 
#> Call:
#> lm(formula = formula, data = data[which(strata == s), ])
#> 
#> Coefficients:
#> (Intercept)         sexf  
#>   0.0952381   -0.0003476  
#> 
#> 
#> [[2]]
#> 
#> Call:
#> lm(formula = formula, data = data[which(strata == s), ])
#> 
#> Coefficients:
#> (Intercept)         sexf  
#>     0.20000     -0.08489
lr$predict(head(pbc))
#> [1] 0.09489051 0.09489051 0.09523810 0.09489051 0.11510791 0.11510791