Support covariates in endpoints()

1 files changed, 27 insertions, 6 deletions

diff --git a/R/endpoints.R b/R/endpoints.R
index 9ea0e598..0aafd728 100644
--- a/R/endpoints.R
+++ b/R/endpoints.R
@@ -14,10 +14,11 @@
 #' another object that has list components mkinmod containing an [mkinmod]
 #' degradation model, and two numeric vectors, bparms.optim and bparms.fixed,
 #' that contain parameter values for that model.
-#' @param \dots Possibility to specify values for the covariates in the model.
-#' In case more than one vector is given, they either have to be of the same
-#' length, or of length one, in which case the respective covariate values are
-#' recycled.
+#' @param covariates Numeric vector with covariate values for all variables in
+#' any covariate models in the object. If given, it overrides 'covariate_quantile'.
+#' @param covariate_quantile This argument only has an effect if the fitted
+#' object has covariate models. If so, the default is to show endpoints
+#' for the median of the covariate values (50th percentile).
 #' @importFrom stats optimize
 #' @return A list with a matrix of dissipation times named distimes, and, if
 #' applicable, a vector of formation fractions named ff and, if the SFORB model
@@ -38,14 +39,34 @@
 #'   }
 #'
 #' @export
-endpoints <- function(fit, ..., covariates = mean) {
+endpoints <- function(fit, covariates = NULL, covariate_quantile = 0.5) {
   mkinmod <- fit$mkinmod
   obs_vars <- names(mkinmod$spec)
 
-  degparms <- c(fit$bparms.optim, fit$bparms.fixed)
+  if (!is.null(fit$covariate_models)) {
+    if (is.null(covariates)) {
+      covariates = as.data.frame(
+        apply(fit$covariates, 2, quantile,
+          covariate_quantile, simplify = FALSE))
+    } else {
+      covariates <- data.frame("User" = covariates)
+    }
+    degparms_trans <- parms(fit, covariates = covariates)[, 1]
+    if (inherits(fit, "saem.mmkin") & (fit$transformations == "saemix")) {
+      degparms <- degparms_trans
+    } else {
+      degparms <- backtransform_odeparms(degparms_trans,
+        fit$mkinmod,
+        transform_rates = fit$transform_rates,
+        transform_fractions = fit$transform_fractions)
+    }
+  } else {
+    degparms <- c(fit$bparms.optim, fit$bparms.fixed)
+  }
 
   # Set up object to return
   ep <- list()
+  ep$covariates <- covariates
   ep$ff <- vector()
   ep$SFORB <- vector()
   ep$distimes <- data.frame(