Support for flexsurv models (#781)

* Inital work to support flexsurv models

* Add model to sanity_model.R

* Add models to CSV file

* Allow for multiple predictions per covariate pattern

* Add class prediction types and whitelist valid dot args

* flexsurv: finalize

* bump

---------

Co-authored-by: Vincent Arel-Bundock <[email protected]>

DESCRIPTION

@@ -77,6 +77,7 @@ Suggests:
     equivalence,
     estimatr,
     fixest,
+    flexsurv,
     fmeffects,
     fontquiver,
     future,
@@ -231,6 +232,7 @@ Collate:
     'methods_dataframe.R'
     'methods_dbarts.R'
     'methods_fixest.R'
+    'methods_flexsurv.R'
     'methods_gamlss.R'
     'methods_glmmTMB.R'
     'methods_glmx.R'

NAMESPACE

@@ -57,6 +57,7 @@ S3method(get_predict,coxph)
 S3method(get_predict,crch)
 S3method(get_predict,default)
 S3method(get_predict,fixest)
+S3method(get_predict,flexsurvreg)
 S3method(get_predict,gamlss)
 S3method(get_predict,glimML)
 S3method(get_predict,glm)
@@ -127,6 +128,7 @@ S3method(set_coef,betareg)
 S3method(set_coef,crch)
 S3method(set_coef,data.frame)
 S3method(set_coef,default)
+S3method(set_coef,flexsurvreg)
 S3method(set_coef,gamlss)
 S3method(set_coef,glimML)
 S3method(set_coef,glm)

NEWS.md

@@ -7,6 +7,10 @@ Breaking changes:
 * `datagrid()` no longer includes the response variable by default when it is not explicitly specified by the user. Use the new `response` argument to include it.
 * `datagrid(grid_type="balanced")` returns factors ordered by level rather than by order of appearance in the original data.
 
+New modeling packages supported:
+
+* `flexsurv`: Thanks to @mattwarkentin for code contributions in PR #781. https://cran.r-project.org/web/packages/flexsurv/index.html
+
 New:
 
 * `wts=TRUE` tries to retrieves weights used in a weighted fit such as `lm()` with the `weights` argument or a model fitted using the `survey` package. Thanks to @ngreifer for feature request 

R/methods_flexsurv.R

@@ -0,0 +1,46 @@
+#' @rdname set_coef
+#' @export
+set_coef.flexsurvreg <- function(model, coefs, ...) {
+  out <- model
+  out$res[, 1] <- coefs
+  out$coefficients <- coefs
+  return(out)
+}
+
+
+#' @rdname get_predict
+#' @export
+get_predict.flexsurvreg <- function(model, newdata, type, ...) {
+  preds <- stats::predict(
+    object = model,
+    newdata = newdata,
+    type = type,
+    ...
+  )
+
+  if (ncol(preds) == 1L) {
+    if (names(preds) == '.pred') {
+      gp <- unlist(lapply(preds$.pred, function(x) {x[, 1, drop = TRUE]}))
+      val <- unlist(lapply(preds$.pred, function(x) {x[, 2, drop = TRUE]}))
+
+      out <- data.frame(
+        rowid = seq_len(nrow(preds)),
+        group = as.vector(gp),
+        estimate = as.vector(val)
+      )
+      return(out)
+    }
+    out <- data.frame(
+      rowid = seq_len(nrow(preds)),
+      estimate = as.vector(preds[, 1, drop = TRUE])
+    )
+    return(out)
+  }
+
+  out <- data.frame(
+    rowid = seq_len(nrow(preds)),
+    group = as.vector(preds[, 1, drop = TRUE]),
+    estimate = as.vector(preds[, 2, drop = TRUE])
+  )
+  out
+}

R/sanity_dots.R

@@ -40,6 +40,9 @@ sanity_dots <- function(model, calling_function = NULL, ...) {
     valid[["bife"]] <- c("alpha_new", "corrected") # nlme::lme
     valid[["process_error"]] <-  # mvgam::mvgam
 
+    # flexsurv
+    valid[["flexsurvreg"]] <- c("times", "p", "start")
+
     white_list <- c(
         "conf.int", "modeldata", "internal_call", "df",
         "transform", "comparison", "side", "delta", "null", "equivalence", "draw",

R/sanity_model.R

@@ -43,6 +43,7 @@ sanity_model_supported_class <- function(model) {
         c("clmm2", "clm2"),
         "coxph",
         "crch",
+        "flexsurvreg", # package: flexsurv
         "fixest",
         "flic",
         "flac",

R/type_dictionary.R

@@ -51,6 +51,15 @@ hxlr,scale
 hxlr,density
 ivpml,pr
 ivpml,xb
+flexsurvreg,response,
+flexsurvreg,mean,
+flexsurvreg,link,
+flexsurvreg,lp,
+flexsurvreg,linear,
+flexsurvreg,rmst,
+flexsurvreg,survival,
+flexsurvreg,hazard,
+flexsurvreg,cumhaz,
 fixest,invlink(link)
 fixest,response
 fixest,link

data-raw/supported_models.csv

@@ -24,6 +24,8 @@ DCchoice,oohbchoice,TRUE,TRUE,,,,,,
 estimatr,lm_lin,TRUE,TRUE,,,,,,
 estimatr,lm_robust,TRUE,TRUE,TRUE,TRUE,TRUE,U,TRUE,TRUE
 estimatr,iv_robust,TRUE,TRUE,TRUE,TRUE,U,U,U,U
+flexsurv,flexsurvreg,TRUE,TRUE,,,,,,
+flexsurv,flexsurvspline,TRUE,TRUE,,,,,,
 fixest,feols,TRUE,TRUE,TRUE,TRUE,U,U,U,U
 fixest,feglm,TRUE,TRUE,,,U,U,U,U
 fixest,fenegbin,TRUE,TRUE,,,U,U,U,U
@@ -46,8 +48,8 @@ MASS,glmmPQL,TRUE,TRUE,,,U,U,TRUE,TRUE
 MASS,glm.nb,TRUE,TRUE,TRUE,TRUE,TRUE,TRUE,TRUE,TRUE
 MASS,polr,TRUE,TRUE,TRUE,TRUE,FALSE,FALSE,TRUE,TRUE
 MASS,rlm,TRUE,TRUE,,,TRUE,TRUE,TRUE,TRUE
-mclogit,mblogit,TRUE,TRUE,,,U,U,U,U 
-mclogit,mclogit,TRUE,TRUE,,,U,U,U,U 
+mclogit,mblogit,TRUE,TRUE,,,U,U,U,U
+mclogit,mclogit,TRUE,TRUE,,,U,U,U,U
 MCMCglmm,MCMCglmm,TRUE,TRUE,U,U,U,U,U,U
 mgcv,gam,TRUE,TRUE,,,U,U,TRUE,TRUE
 mgcv,bam,TRUE,TRUE,,,U,U,TRUE,FALSE

inst/tinytest/test-pkg-flexsurv.R

@@ -0,0 +1,12 @@
+source("helpers.R")
+if (!requiet("flexsurv")) exit_file("flexsurv")
+
+
+mod <- flexsurvreg(formula = Surv(futime, fustat) ~ age + ecog.ps, data = ovarian, dist= "gengamma")
+x <- avg_slopes(mod)
+expect_inherits(x, "slopes")
+x <- predictions(mod)
+expect_inherits(x, "predictions")
+x <- comparisons(mod)
+expect_inherits(x, "comparisons")
+