# Partial dependence plot for glm in r -- why linear?

I'd like to understand why my partial dependence plots for a logistic regression model simply show up as straight lines -- even when I'd expect basically a threshold effect from a covariate. I know partial dependence plots are typical of machine learning, but the (excellent) description by the authors of the pdp] package suggest glms are fair game. So why does the relationship between outcome and effort (below) appear to be linear?

Here's a dummy dataset. Note that I forced higher values of effort for outcomes corresponding to 1 (a "win"). Also note that sometimes the algorithm won't converge -- if that's the case, just generate new data.

library(pdp)
library(randomForest)

# Sample game data
outcome <- as.vector(cbind(rep(0,25), rep(1,25)))
effort <- as.vector(cbind(rnorm(25, 25, 5), rnorm(25, 50, 10)))
skill <- rnorm(50, 50, 20)
game <- cbind(outcome, effort, skill) %>% as.data.frame()

# Simple glm
mod <- glm(outcome ~ effort + skill, data = game, family = binomial(link = "logit"))
summary(mod)
partial(mod, pred.var = c("effort"), plot = TRUE)

Call:
glm(formula = outcome ~ effort + skill, family = binomial(link = "logit"),
data = game)

Deviance Residuals:
Min        1Q    Median        3Q       Max
-1.26979  -0.13985  -0.00751   0.01736   2.34734

Coefficients:
Estimate Std. Error z value Pr(>|z|)
(Intercept) -15.12758    5.73393  -2.638  0.00833 **
effort        0.50174    0.19218   2.611  0.00903 **
skill        -0.05414    0.05142  -1.053  0.29231


Clearly, effort is going to be a strong predictor -- with way more wins (1s) associated with higher effort (given my data assignments). However, the partial dependence plot looks like this:

partial(mod, pred.var = c("effort"), plot = TRUE)


If I use a random forest instead, that threshold effect shows up. (Yes, I know it throws a warning about using <5 unique response values in regression. It also shows up if you force outcome to be a factor.)

rf <- randomForest(outcome ~ effort + skill, data = game)
partial(rf, pred.var = c("effort"), plot = TRUE)


My primary question here is not about which model is a better fit, but why the partial dependence is apparently linear with the logistic regression? Why doesn't that 30-40 range pop out as a threshold in the glm plot? Is that truly representing the relationship between game and effort in the model?

Thanks for any insights!

A partial dependence plot for a logistic-type model is constructed by setting all but one feature to fixed, static values, varying the remaining feature throughout a range, and plotting:

$$t \mapsto \log \left( \frac{p}{1-p} \right)$$

Where $$p$$ is the (probability) prediction for your model when the varied feature is set to the value $$t$$. Note that, in particular, the $$y$$-axis of a partial dependency plot is measured on the log-odds scale, not the probability scale.

For a standard logistic regression, the functional form of your model is:

$$\log \left( \frac{p}{1-p} \right) = \beta_0 + \beta_1 x_1 + \cdots + \beta_k x_k$$

So the form of the partial dependence plot is:

$$t \mapsto \beta_j t + \text{constant}$$

where $$j$$ is the index of the feature you are constructing the partial dependence plot of. This is why you get a line, the slope of that line is the parameter estimate $$\hat \beta_j$$ in the regression.

In a random forest the functional form of your model is:

$$p = \text{average} \left( T_0(x), T_1(x), \ldots, T_{\text{n_trees}}(x) \right)$$

where the $$T(x)$$'s are the probability predictions from your individual classification trees. So the partial dependence plot is the unwieldy:

$$t \mapsto \log \left( \frac{p}{1-p} \right) = \frac{\text{average} \left( T_0(t), T_1(t), \ldots, T_{\text{n_trees}}(t) \right)}{1 - \text{average} \left( T_0(t), T_1(t), \ldots, T_{\text{n_trees}}(t) \right)}$$

This can be a very complicated, non-linear function of any individual feature, resulting in a vast multitude of possible shapes for the partial dependence plots. The fact that you are seeing a soft threshold shape is due to the particulars of the problem you are solving, not something structural about partial dependence plots.

• Thanks @Matthew-Drury. That makes really good sense. Glad to be reminded that the y-axis is on log-odds scale, too. Feb 27 '19 at 20:10
• It should be, though looking at that random forest plot, it looks like R is being inconsistent and the logistic one is on the log-odds scale, but the forest is on the probability scale... : ( I just noticed that... Feb 27 '19 at 20:13
• Perhaps because rf is forcing a regression (vs. classification) but didn't figure out the underlying distribution? Feb 27 '19 at 20:15

The pdp package allows converting back to probability scale by setting "prob = TRUE". You may also add the rugs to display the distribution of predictors.

partial(mod, pred.var = c("effort"), plot = TRUE, prob = TRUE, rug = TRUE)