mhmkin: Easy specification of ill-defined parms

The argument 'no_random_effect' now accepts an illparms.mhmkin object

1 files changed, 20 insertions, 18 deletions

diff --git a/R/mhmkin.R b/R/mhmkin.R
index 1f29dc40..6a61e8bb 100644
--- a/R/mhmkin.R
+++ b/R/mhmkin.R
@@ -12,13 +12,13 @@
 #' degradation models to the same data
 #' @param backend The backend to be used for fitting. Currently, only saemix is
 #' supported
-#' @param no_random_effect Default is NULL and will be passed to [saem]. If
-#' you specify "auto", random effects are only included if the number
-#' of datasets in which the parameter passed the t-test is at least 'auto_ranef_threshold'.
-#' Beware that while this may make for convenient model reduction or even
-#' numerical stability of the algorithm, it will likely lead to
-#' underparameterised models.
-#' @param auto_ranef_threshold See 'no_random_effect.
+#' @param no_random_effect Default is NULL and will be passed to [saem]. If a
+#' character vector is supplied, it will be passed to all calls to [saem],
+#' regardless if the corresponding parameter is in the model. Alternatively,
+#' an object of class [illparms.mhmkin] can be specified. This has to have
+#' the same dimensions as the return object of the current call. In this way,
+#' ill-defined parameters found in corresponding simpler versions of the
+#' degradation model can be specified.
 #' @param algorithm The algorithm to be used for fitting (currently not used)
 #' @param \dots Further arguments that will be passed to the nonlinear mixed-effects
 #' model fitting function.
@@ -51,7 +51,7 @@ mhmkin.mmkin <- function(objects, ...) {
 #' @export
 #' @rdname mhmkin
 mhmkin.list <- function(objects, backend = "saemix", algorithm = "saem",
-  no_random_effect = NULL, auto_ranef_threshold = 3,
+  no_random_effect = NULL,
   ...,
   cores = if (Sys.info()["sysname"] == "Windows") 1 else parallel::detectCores(), cluster = NULL)
 {
@@ -102,22 +102,24 @@ mhmkin.list <- function(objects, backend = "saemix", algorithm = "saem",
     deg_index <- w[1]
     error_index <- w[2]
     mmkin_row <- objects[[error_index]][deg_index, ]
-    if (identical(no_random_effect, "auto")) {
-      ip <- illparms(mmkin_row)
-      ipt <- table(unlist(lapply(as.list(ip), strsplit, ", ")))
-      no_ranef <- names(ipt)[which((length(ds) - ipt) <= auto_ranef_threshold)]
-      transparms <- rownames(mmkin_row[[1]]$start_transformed)
-      bparms <- rownames(mmkin_row[[1]]$start)
-      names(transparms) <- bparms
-      no_ranef_trans <- transparms[no_ranef]
+    if (is(no_random_effect, "illparms.mhmkin")) {
+      if (identical(dim(no_random_effect), dim(fit_indices))) {
+        no_ranef_split <- strsplit(no_random_effect[[fit_index]], ", ")
+        no_ranef <- sapply(no_ranef_split, function(x) {
+          gsub("sd\\((.*)\\)", "\\1", x)
+        })
+      } else {
+        stop("Dimensions of illparms.mhmkin object given as 'no_random_effect' are not suitable")
+      }
     } else {
-      no_ranef_trans <- no_random_effect
+      no_ranef <- no_random_effect
     }
     res <- try(do.call(backend_function,
-        args = c(list(mmkin_row), dot_args, list(no_random_effect = no_ranef_trans))))
+        args = c(list(mmkin_row), dot_args, list(no_random_effect = no_ranef))))
     return(res)
   }
 
+
   fit_time <- system.time({
     if (is.null(cluster)) {
       results <- parallel::mclapply(as.list(1:n.fits), fit_function,