How to find the root process of a custom plot function inside a package

Question

In my understanding, some customised object classes have a customised plot function, which allows to use plot(object_name) and get the plot which you are (most probably) interested to.

For example, I can get from the help page of pmcalibration::pmcalibration a code like the one below

library(pmcalibration)
library(survival)
  
data('transplant', package="survival")
transplant <- na.omit(transplant)
transplant = subset(transplant, futime > 0)
transplant$ltx <- as.numeric(transplant$event == "ltx")
  
cph <- coxph(Surv(futime, ltx) ~ age + sex + abo + year, data = transplant)
  
time <- 100
newd <- transplant; newd$futime <- time; newd$ltx <- 1
p <- 1 - exp(-predict(cph, type = "expected", newdata=newd))
y <- with(transplant, Surv(futime, ltx))
  
cal <- pmcalibration(y = y, p = p, smooth = "rcs", nk=5, ci = "pw", time = time)

class(cal)
  [1] "pmcalibration"

plot(cal)

As shown, cal has class "pmcalibration", and it is enough to provide it as the argument for the plot function to get the calibration plot you are almost certainly interested to after using the pmcalibration function

My question is: how can I retrieve what has the plot function done here? How can I get which values and variables inside cal have been used to create the plot? What is happening underneath?

For this specific example, I can see that a hint might be given by the existence of the section cal$plot; however, I still cannot understand which variables have been used and which functions have been applied to the values stored there to finally be able to create the plot.

Answer 1

You could use sloop::s3_dispatch to see which method is used,

> sloop::s3_dispatch(plot(cal))
=> plot.pmcalibration
 * plot.default

and with a good guess of the package, inspect (likely non-exported) source code.

> pmcalibration:::plot.pmcalibration
function (x, conf_level = 0.95, ...) 
{
    dots <- list(...)
    pdat <- get_cc(x, conf_level = conf_level)
    if ("xlim" %in% names(dots)) 
        xlim <- dots[["xlim"]]
    else xlim <- c(0, 1)
    if ("ylim" %in% names(dots)) 
        ylim <- dots[["ylim"]]
    else ylim <- c(0, 1)
    if ("xlab" %in% names(dots)) 
        xlab <- dots[["xlab"]]
    else xlab <- "Predicted Probability"
    if ("ylab" %in% names(dots)) 
        ylab <- dots[["ylab"]]
    else ylab <- "Estimated Probability"
    plot(x = pdat$p, y = pdat$p_c, type = "l", xlim = xlim, ylim = ylim, 
        xlab = xlab, ylab = ylab)
    abline(0, 1, lty = 2, col = "grey")
    if ("lower" %in% colnames(pdat)) {
        lines(x = pdat$p, y = pdat$lower, lty = 2)
        lines(x = pdat$p, y = pdat$upper, lty = 2)
    }
}
<bytecode: 0x5559b1f0d810>
<environment: namespace:pmcalibration>
>