Permutation test for mob() tree in partykit

Question

I have data of $N\approx 1200$ whereby treated and control individuals have been matched (via full-matching) as a pre-processing step. This matching step induces correlations between treated and control individuals. The next step in the analysis is to use mob() to identify any subgroups in my data based on some partitioning variables.

However, the variable selection step for mob() that is currently implemented in the partykit package only uses M-fluctuation tests based on normal approximation. For my data, I think this approximation isn't appropriate due to the correlation from the pre-processing step. Is there a way to do permutation tests via mob_control()? I only see it in the ctree_control function.

Answer 1

At the moment, it is not possible to combine mob() in partykit with permutation tests a la ctree(). We are working on a new implementation that is more flexible but we are not quite there, yet.

However, you can also combine ctree() with models fitted in every node. There are differences between the algorithms but often the results are similar. We've written up some comments regarding the differences in Schlosser et al. (2019, arXiv:1906.10179). In that paper we use only the asymptotic tests, though. Below I include a short worked example.

Most importantly, though, I'm not sure whether permutation tests are the best solution for your problem. They can be much better for small sample sizes but this does not appear to be an issue in your setup. While the permutation tests do not require independence, they still require exachangability and it is not clear to me whether that is fulfilled with your correlation pattern. So before abandoning mob() it might be worth to consider the cluster argument or supply some sandwich covariance matrix. Yet another option might be to combine it with mixed effects models as implemented in glmertree.

Finally, the worked example. One of the illustration for mob() combined with a binary logit glm() from the partykit examples is for the PimaIndiansDiabetes data from mlbench:

library("partykit")
data("PimaIndiansDiabetes", package = "mlbench")
mb <- glmtree(diabetes ~ glucose | mass + age,
  data = PimaIndiansDiabetes, family = binomial)
mb
## Generalized linear model tree (family: binomial)
## 
## Model formula:
## diabetes ~ glucose | pregnant + mass + pedigree + age
## 
## Fitted party:
## [1] root
## |   [2] mass <= 26.3: n = 167
## |       (Intercept)     glucose 
## |       -9.95150963  0.05870786 
## |   [3] mass > 26.3
## |   |   [4] age <= 30: n = 304
## |   |       (Intercept)     glucose 
## |   |       -6.70558554  0.04683748 
## |   |   [5] age > 30: n = 297
## |   |       (Intercept)     glucose 
## |   |       -2.77095386  0.02353582 
## 
## Number of inner nodes:    2
## Number of terminal nodes: 3
## Number of parameters per node: 2
## Objective function (negative log-likelihood): 355.4578

For fitting a similar tree with ctree() you need to supply a transformation function that fits the model and extracts the corresponding scores. Also, you need to pass both the actual dependent variable and the regressor variables through the dependent variable of the ctree():

logitscore <- function(y, x, start = NULL, weights, offset, estfun = TRUE, object = FALSE, ...) {
    if (is.null(weights)) weights <- rep(1, NROW(y))
    s <- sandwich::estfun(glm(diabetespos ~ glucose, data = as.data.frame(y),
      family = binomial, weights = weights, subset = weights > 0))
    r <- matrix(0, nrow = length(weights), ncol = NCOL(s))
    r[weights > 0, ] <- s
    list(estfun = r, converged = TRUE)
}
ct <- ctree(diabetes + glucose ~ mass + age,
  data = PimaIndiansDiabetes, ytrafo = logitscore)
ct
## Model formula:
## ~diabetes + glucose + (mass + age)
## 
## Fitted party:
## [1] root
## |   [2] mass <= 27.3: *
## |   [3] mass > 27.3
## |   |   [4] age <= 30: *
## |   |   [5] age > 30: *
## 
## Number of inner nodes:    2
## Number of terminal nodes: 3

So you get the desired tree structure but you would still have to refit the model in the nodes yourself if you want to have the coefficients etc.