How to change the line type in the calibration.plot of the GBM package?

61 Views Asked by At

I have a model that I calibrated using a logistic regression. I plotted the data using the calibration.plot from the gbm package. Now that I am submitting the paper for publication, the publisher requested the figures to be black and white, and therefore I need to change the line types. However, passing lty to the line.par argument doesn't do anything. Is there any other suggestion?

Example code:

library(gbm)
set.seed(123)
data <- data.frame(y = sample(c(0,1), size = 100, replace = TRUE),
                   pred = runif(100, min = 0, max = 1))
calibrate.plot(data$y, data$pred, 
               distribution = 'bernouli', 
               replace = TRUE, 
               shade.col = NA, 
               line.par = list(col ="black", lty= 2))
2

There are 2 best solutions below

1
Allan Cameron On

Unfortunately the red line is hard-coded in the function.

If you type calibrate.plot into your console, you will see that the last line is

    abline(0, 1, col = "red")

So there are no parameters that can be passed into the function to make it draw anything other than a red line.

However, the other thing that is obvious from looking at the code is what it does. It fits a binomial glm using natural cubic splines, gets the predictions along the range of the predictors and plots these, along with a rug plot and an abline.

The following code is therefore completely equivalent:

library(splines)

newdata <- data.frame(pred = seq(min(data$y), max(data$y), len = 100))

newdata$y <- glm(y ~ ns(pred, df = 6), data = data, family = binomial) |>
  predict(newdata = newdata, type = "response")

plot(newdata$pred, newdata$y, type = "l", ylab = "Observed average",
     xlab = "Predicted value", xlim = range(newdata$pred), ylim = range(newdata$pred) ) 
quantile.rug(data$pred, side = 1)
abline(0, 1, col = "gray")

enter image description here

We can see this gives us a black-and-white equivalent of calibrate.plot

calibrate.plot(data$y, data$pred, 
               distribution = 'bernouli', 
               replace = TRUE, 
               shade.col = NA, 
               col = "black",
               line.par = list(col ="black", lty= 2))

enter image description here

0
Sayf Said On

It seems the line function within the gbm::calibrate.plot only calls for a color argument to be passed from line.par argument. I made a clone function that takes line.par$lty argument to be passed to the line function within the gbm::calibration.plot.

clone_calibrate.plot <- function (y, p, distribution = "bernoulli", replace = TRUE, line.par = list(col = "black"), 
                                shade.col = "lightyellow", shade.density = NULL, rug.par = list(side = 1), 
                                xlab = "Predicted value", ylab = "Observed average", xlim = NULL, 
                                ylim = NULL, knots = NULL, df = 6, ...) 
{
  if (!requireNamespace("splines", quietly = TRUE)) {
    stop("The splines package is needed for this function to work. Please ", 
         "install it.", call. = FALSE)
  }
  data <- data.frame(y = y, p = p)
  if (is.null(knots) && is.null(df)) {
    stop("Either knots or df must be specified")
  }
  if ((df != round(df)) || (df < 1)) {
    stop("df must be a positive integer")
  }
  if (distribution == "bernoulli") {
    family1 <- binomial
  }
  else if (distribution == "poisson") {
    family1 <- poisson
  }
  else {
    family1 <- gaussian
  }
  gam1 <- glm(y ~ splines::ns(p, df = df, knots = knots), data = data, 
              family = family1)
  x <- seq(min(p), max(p), length = 200)
  yy <- predict(gam1, newdata = data.frame(p = x), se.fit = TRUE, 
                type = "response")
  x <- x[!is.na(yy$fit)]
  yy$se.fit <- yy$se.fit[!is.na(yy$fit)]
  yy$fit <- yy$fit[!is.na(yy$fit)]
  if (!is.na(shade.col)) {
    se.lower <- yy$fit - 2 * yy$se.fit
    se.upper <- yy$fit + 2 * yy$se.fit
    if (distribution == "bernoulli") {
      se.lower[se.lower < 0] <- 0
      se.upper[se.upper > 1] <- 1
    }
    if (distribution == "poisson") {
      se.lower[se.lower < 0] <- 0
    }
    if (is.null(xlim)) {
      xlim <- range(se.lower, se.upper, x)
    }
    if (is.null(ylim)) {
      ylim <- range(se.lower, se.upper, x)
    }
  }
  else {
    if (is.null(xlim)) {
      xlim <- range(yy$fit, x)
    }
    if (is.null(ylim)) {
      ylim <- range(yy$fit, x)
    }
  }
  if (replace) {
    plot(0, 0, type = "n", xlab = xlab, ylab = ylab, xlim = xlim, 
         ylim = ylim, ...)
  }
  if (!is.na(shade.col)) {
    polygon(c(x, rev(x), x[1L]), c(se.lower, rev(se.upper), 
                                   se.lower[1L]), col = shade.col, border = NA, density = shade.density)
  }
  lines(x, yy$fit, col = line.par$col, lty= line.par$lty) ## Added here
  quantile.rug(p, side = rug.par$side)
  abline(0, 1, col = "black") # The diagonal line can be changed too
}

Now passing the parameters in the question's stem to the clone function should change the line type.

clone_calibrate.plot(
  data$y,
  data$pred, 
  distribution = 'bernouli', 
  replace = TRUE, 
  shade.col = NA, 
  line.par = list(col ="black", lty= 2)
)