How can I write a function in R to calculate and plot a Kaplan-Meier curve depending on the independent variable?

I have a dataset my_data for survival analysis with

time until a certain event of interest,
status coded as 0 and 1, indicating whether the event has occurred or not at time, and
several factors (factor1, factor2, ...) that might influence the occurrence of the event.

I analyze the data in R, using the survival and survminer package, à la

surv_object1 <- survfit(Surv(time, status) ~ factor1, data = my_data)
ggsurvplot(surv_object1)

surv_object2 <- survfit(Surv(time, status) ~ factor2, data = my_data)
ggsurvplot(surv_object2)

etc.

Instead of writing multiple lines of code, I would like to write a function that expects a factor as a variable and returns a survival plot.

So far, I have tried writing

survival_plot <- function(factor) {
    
    survival_object <- survfit(Surv(time, status) ~ factor, data = my_data)
    surv_plot <- ggsurvplot(survival_object)
    
    return(surv_plot)
}

However, when running

survival_plot(factor1)

it returns

Error in eval(predvars, data, env) : object 'factor1' not found

I have already double checked the spelling of my factors, made sure, that the columns that I name do exist in my dataframe, and tried !!sym{}, but without any luck.

I also tried the solution pointed out here, but again without any success.

Even though I can write the lines as pointed out above, I am still wondering what I am doing wrong. Any help would be much appreciated!

A minimal working example should be

library(survival)
library(survminer)

survival_plot <- function(factor) {
    
    s_obj <- survfit(Surv(rfstime, status) ~ factor, data = survival::gbsg)
    surv_plot <- ggsurvplot(s_obj)
    
    return(surv_plot)
}

survival_plot(meno)

And for the solution of this thread

survival_function <- function(strata_x) {
  
  gbsg$s <- Surv(gbsg$rfstime, gbsg$status)
  survFormula <- as.formula(paste("s ~", strata_x))
  my_survfit <- survfit(survFormula, data = gbsg)
  my_survfit$call$formula <- survFormula
  ggsurvplot(my_survfit, data = gbsg)

}

Solution

ggsurvplot will look in s_obj$call$formula to find what model was fit.

You can e.g. make your own call, and use eval(parse()), or overwrite the wrong call. These both work:

survival_plot <- function(factor) {
  call <- glue::glue('s_obj <- survfit(Surv(rfstime, status) ~ {factor}, data = survival::gbsg)')
  eval(parse(text = call))
  surv_plot <- ggsurvplot(s_obj)
  
  return(surv_plot)
}


survival_plot2 <- function(factor) {
  form <- as.formula(glue::glue('Surv(rfstime, status) ~ {factor}'))
  
  s_obj <- survfit(form, data = survival::gbsg)
  s_obj$call$formula <- form
  surv_plot <- ggsurvplot(s_obj)
  
  return(surv_plot)
}

Use as

survival_plot('meno')
survival_plot2('meno')