From 9eca01988f3507d5b06b26f483634d85908e6554 Mon Sep 17 00:00:00 2001 From: "Matthew L. Fidler" Date: Sat, 2 Dec 2023 21:48:36 -0600 Subject: [PATCH 1/5] Test/add changing variable type with ~ka --- R/piping-ini.R | 66 +++++++++++++++++ R/piping-model.R | 4 ++ R/piping.R | 10 +-- tests/testthat/test-piping-ini.R | 117 +++++++++++++++++++++++++++++++ 4 files changed, 193 insertions(+), 4 deletions(-) diff --git a/R/piping-ini.R b/R/piping-ini.R index 4aa99fb16..797a38624 100644 --- a/R/piping-ini.R +++ b/R/piping-ini.R @@ -408,6 +408,70 @@ assign("iniDf", ret, envir=rxui) invisible() } +#' Handle switching theta to eta and vice versa +#' +#' This is coded as model |> ini(~par) +#' +#' @param expr Expression, this would be the ~par expression +#' @param rxui rxui uncompressed environment +#' @param envir Environment for evaluation (if needed) +#' @return Nothing, called for side effects +#' @noRd +#' @author Matthew L. Fidler +.iniHandleSwitchType <- function(expr, rxui, envir=parent.frame()) { + .var <- as.character(expr[[2]]) + .iniDf <- rxui$iniDf + .w <- which(.iniDf$name == .var) + if (length(.w) != 1L) stop("cannot switch parameter type for '", .var, "'", call.=FALSE) + .theta <- .iniDf[!is.na(.iniDf$ntheta),, drop = FALSE] + .eta <- .iniDf[is.na(.iniDf$ntheta),, drop = FALSE] + if (is.na(.iniDf$ntheta[.w])) { + # switch eta to theta + .neta <- .iniDf$neta1[.w] + .eta <- .eta[.eta$neta1 != .neta,, drop = FALSE] + .eta <- .eta[.eta$neta2 != .neta,, drop = FALSE] + .eta$neta1 <- .eta$neta1 - ifelse(.eta$neta1 < .neta, 0L, 1L) + .eta$neta2 <- .eta$neta2 - ifelse(.eta$neta2 < .neta, 0L, 1L) + .newTheta <- .iniDf[.w, ] + .newTheta$neta1 <- NA_integer_ + .newTheta$neta2 <- NA_integer_ + if (length(.theta$ntheta) == 0L) { + .newTheta$ntheta <- 1L + } else { + .newTheta$ntheta <- max(.theta$ntheta) + 1L + } + .minfo(paste0("convert '", .var, "' from between subject variability to population parameter")) + .theta <- rbind(.theta, .newTheta) + } else { + # switch theta to eta + if (!is.na(.iniDf$err[.w])){ + stop("cannot switch error parameter '", .var, + "' to a different type", call. = FALSE) + } + .ntheta <- .iniDf$ntheta[.w] + .theta <- .theta[.theta$ntheta != .ntheta,, drop = FALSE] + .theta$ntheta <- .theta$ntheta - ifelse(.theta$ntheta < .ntheta, 0L, 1L) + .newEta <- .iniDf[.w, ] + .newEta$ntheta <- NA_integer_ + if (length(.eta$neta1) == 0L) { + .newEta$neta1 <- .newEta$neta2 <- 1L + } else { + .newEta$neta1 <- .newEta$neta2 <- max(.eta$neta1) + 1L + } + .minfo(paste0("convert '", .var, "' from population parameter to between subject variability")) + if (.newEta$est == 0) { + .minfo("old initial estimate is zero, changing to 1") + .newEta$est <- 1 + } else if (.newEta$est < 0) { + .minfo("old initial estimate was negative, changing to positive") + .newEta$est <- -.newEta$est + } + .eta <- rbind(.eta, .newEta) + } + .ini <- rbind(.theta, .eta) + assign("iniDf", .ini, envir=rxui) + invisible() +} #' Update the iniDf of a model #' @@ -464,6 +528,8 @@ expr[[3]] <- eval(as.call(list(quote(`lotri`), as.call(list(quote(`{`), expr)))), envir=envir)[1, 1] .iniHandleFixOrUnfixEqual(expr=expr, rxui=rxui, envir=envir, maxLen=1L) + } else if (.isTildeExpr(expr)) { + .iniHandleSwitchType(expr=expr, rxui=rxui, envir=envir) } else { # Can this error be improved to clarify what is the expression causing the # issue? It needs a single character string representation of something diff --git a/R/piping-model.R b/R/piping-model.R index c666b22e8..531b955fc 100644 --- a/R/piping-model.R +++ b/R/piping-model.R @@ -180,6 +180,10 @@ model.rxModelVars <- model.rxode2 .matchesLangTemplate(expr, str2lang(". = .")) } +.isTildeExpr <- function(expr) { + .matchesLangTemplate(expr, str2lang("~ .")) +} + # get the left hand side of an assignment or endpoint; returns NULL if the input # is not an assignment or endpoint .getLhs <- function(expr) { diff --git a/R/piping.R b/R/piping.R index 7ba76fdd0..93da808b6 100644 --- a/R/piping.R +++ b/R/piping.R @@ -295,10 +295,12 @@ } else if (identical(.quoted[[1]], quote(`as.formula`))) { .quoted <- .quoted[[2]] } else if (identical(.quoted[[1]], quote(`~`))) { - .quoted[[3]] <- .iniSimplifyFixUnfix(.quoted[[3]]) - if (identical(.quoted[[3]], quote(`fix`)) || - identical(.quoted[[3]], quote(`unfix`))) { - .quoted <- as.call(list(quote(`<-`), .quoted[[2]], .quoted[[3]])) + if (length(.quoted) == 3L) { + .quoted[[3]] <- .iniSimplifyFixUnfix(.quoted[[3]]) + if (identical(.quoted[[3]], quote(`fix`)) || + identical(.quoted[[3]], quote(`unfix`))) { + .quoted <- as.call(list(quote(`<-`), .quoted[[2]], .quoted[[3]])) + } } } else if (identical(.quoted[[1]], quote(`$`))) { .tmp <- try(eval(.quoted), silent=TRUE) diff --git a/tests/testthat/test-piping-ini.R b/tests/testthat/test-piping-ini.R index 4003a90d3..6eae63850 100644 --- a/tests/testthat/test-piping-ini.R +++ b/tests/testthat/test-piping-ini.R @@ -569,3 +569,120 @@ test_that("append allows promoting from covariate (#472)", { ) expect_equal(newmod$iniDf$name, c("lka", "ka_dose", "lcl", "lvc", "propSd")) }) + +test_that("change ini type with ~", { + + mod <- function() { + ini({ + lka <- 0.45 + lcl <- 1 + lvc <- 3.45 + propSd <- 0.5 + }) + model({ + ka <- exp(lka) + cl <- exp(lcl) + vc <- exp(lvc) + kel <- cl / vc + d/dt(depot) <- -ka*depot + d/dt(central) <- ka*depot-kel*central + cp <- central / vc + cp ~ prop(propSd) + }) + } + + mod1 <- mod |> ini( ~ lka) + expect_equal(mod1$omega, lotri(lka ~ 0.45)) + + mod2 <- mod1 |> ini( ~ lka) + expect_equal(mod2$omega, NULL) + + expect_error(mod1 |> ini( ~ propSd)) + + expect_error(mod1 |> ini( ~ matt)) + + ## all etas + + mod <- function() { + ini({ + lka ~ 0.45 + lcl ~ 1 + lvc ~ 3.45 + }) + model({ + ka <- exp(lka) + cl <- exp(lcl) + vc <- exp(lvc) + kel <- cl / vc + d/dt(depot) <- -ka*depot + d/dt(central) <- ka*depot-kel*central + cp <- central / vc + }) + } + + mod2 <- mod |> ini( ~ lka) + + expect_equal(mod2$omega, lotri(lcl ~ 1, lvc ~ 3.45)) + + # remove correlated eta + + mod <- function() { + ini({ + lka + lcl + lvc ~ + c(0.45, + 0.01, 1, + 0.01, -0.01, 3.45) + }) + model({ + ka <- exp(lka) + cl <- exp(lcl) + vc <- exp(lvc) + kel <- cl / vc + d/dt(depot) <- -ka*depot + d/dt(central) <- ka*depot-kel*central + cp <- central / vc + }) + } + + mod2 <- mod |> ini( ~ lka) + + expect_equal(mod2$omega, lotri(lcl + lvc ~ c(1, + -0.01, 3.45))) + + + # negative and zero + + mod <- function() { + ini({ + lka <- 0.45 + lcl <- -1 + lvc <- 0 + }) + model({ + ka <- exp(lka) + cl <- exp(lcl) + vc <- exp(lvc) + kel <- cl / vc + d/dt(depot) <- -ka*depot + d/dt(central) <- ka*depot-kel*central + cp <- central / vc + }) + } + + mod2 <- mod |> ini( ~ lcl) + + expect_equal(mod2$omega, lotri(lcl ~ 1)) + + mod2 <- mod |> ini( ~ lvc) + + expect_equal(mod2$omega, lotri(lvc ~ 1)) + + mod3 <- mod2 |> ini( ~ lvc) + + expect_equal(mod3$omega, lotri(lcl ~ 1, lvc ~ 1)) + + mod4 <- mod3 |> ini( ~ lvc) + + expect_equal(mod4$omega, lotri(lcl ~ 1)) + +}) From 5b90211a6a0414a5aefade2be18e50ce2c117e70 Mon Sep 17 00:00:00 2001 From: "Matthew L. Fidler" Date: Sat, 2 Dec 2023 22:10:22 -0600 Subject: [PATCH 2/5] Handle demote to covariate and add recalc on switch/demote --- R/piping-ini.R | 56 ++++++++++++++++++++++++++++++++ R/piping-model.R | 4 +++ tests/testthat/test-piping-ini.R | 51 +++++++++++++++++++++++++++++ 3 files changed, 111 insertions(+) diff --git a/R/piping-ini.R b/R/piping-ini.R index 797a38624..d0e7b14bf 100644 --- a/R/piping-ini.R +++ b/R/piping-ini.R @@ -408,6 +408,17 @@ assign("iniDf", ret, envir=rxui) invisible() } + +.iniHandleRecalc <- function(rxui) { + .fun <- rxUiDecompress(rxui$fun()) + for (.i in ls(.fun, all=TRUE)) { + if (.i != "meta") { + assign(.i, get(.i, envir=.fun), envir=rxui) + } + } + invisible() +} + #' Handle switching theta to eta and vice versa #' #' This is coded as model |> ini(~par) @@ -470,6 +481,49 @@ } .ini <- rbind(.theta, .eta) assign("iniDf", .ini, envir=rxui) + # This may change mu-referencing, recalculate everything + .iniHandleRecalc(rxui) + invisible() +} + +#' Handle dropping parameter and treating as if it is a covariate +#' +#' This is coded as model |> ini(-par) +#' +#' @param expr Expression, this would be the ~par expression +#' @param rxui rxui uncompressed environment +#' @param envir Environment for evaluation (if needed) +#' @return Nothing, called for side effects +#' @noRd +#' @author Matthew L. Fidler +.iniHandleDropType <- function(expr, rxui, envir=parent.frame()) { + .var <- as.character(expr[[2]]) + .iniDf <- rxui$iniDf + .w <- which(.iniDf$name == .var) + if (length(.w) != 1L) stop("no initial estimates for '", .var, "', cannot change to covariate", call.=FALSE) + .theta <- .iniDf[!is.na(.iniDf$ntheta),, drop = FALSE] + .eta <- .iniDf[is.na(.iniDf$ntheta),, drop = FALSE] + if (is.na(.iniDf$ntheta[.w])) { + .minfo(paste0("changing between subject variability parameter '", .var, "' to covariate parameter")) + .neta <- .iniDf$neta1[.w] + .eta <- .eta[.eta$neta1 != .neta,, drop = FALSE] + .eta <- .eta[.eta$neta2 != .neta,, drop = FALSE] + .eta$neta1 <- .eta$neta1 - ifelse(.eta$neta1 < .neta, 0L, 1L) + .eta$neta2 <- .eta$neta2 - ifelse(.eta$neta2 < .neta, 0L, 1L) + } else { + if (!is.na(.iniDf$err[.w])){ + stop("cannot switch error parameter '", .var, + "' to a covariate", call. = FALSE) + } + .minfo(paste0("changing population parameter '", .var, "' to covariate parameter")) + .ntheta <- .iniDf$ntheta[.w] + .theta <- .theta[.theta$ntheta != .ntheta,, drop = FALSE] + .theta$ntheta <- .theta$ntheta - ifelse(.theta$ntheta < .ntheta, 0L, 1L) + } + .ini <- rbind(.theta, .eta) + assign("iniDf", .ini, envir=rxui) + # This will change covariates, recalculate everything + .iniHandleRecalc(rxui) invisible() } @@ -530,6 +584,8 @@ .iniHandleFixOrUnfixEqual(expr=expr, rxui=rxui, envir=envir, maxLen=1L) } else if (.isTildeExpr(expr)) { .iniHandleSwitchType(expr=expr, rxui=rxui, envir=envir) + } else if (.isIniDropExpression(expr)){ + .iniHandleDropType(expr=expr, rxui=rxui, envir=envir) } else { # Can this error be improved to clarify what is the expression causing the # issue? It needs a single character string representation of something diff --git a/R/piping-model.R b/R/piping-model.R index 531b955fc..df4f62812 100644 --- a/R/piping-model.R +++ b/R/piping-model.R @@ -184,6 +184,10 @@ model.rxModelVars <- model.rxode2 .matchesLangTemplate(expr, str2lang("~ .")) } +.isIniDropExpression <- function(expr) { + .matchesLangTemplate(expr, str2lang("- .")) +} + # get the left hand side of an assignment or endpoint; returns NULL if the input # is not an assignment or endpoint .getLhs <- function(expr) { diff --git a/tests/testthat/test-piping-ini.R b/tests/testthat/test-piping-ini.R index 6eae63850..bf27620f7 100644 --- a/tests/testthat/test-piping-ini.R +++ b/tests/testthat/test-piping-ini.R @@ -686,3 +686,54 @@ test_that("change ini type with ~", { expect_equal(mod4$omega, lotri(lcl ~ 1)) }) + + + +test_that("change ini variable to covariate with -", { + + mod <- function() { + ini({ + lka + lcl + lvc ~ + c(0.45, + 0.01, 1, + 0.01, -0.01, 3.45) + }) + model({ + ka <- exp(lka) + cl <- exp(lcl) + vc <- exp(lvc) + kel <- cl / vc + d/dt(depot) <- -ka*depot + d/dt(central) <- ka*depot-kel*central + cp <- central / vc + }) + } + + mod2 <- mod |> ini(-lka) + + expect_equal(mod2$allCovs, "lka") + expect_equal(mod2$omega, lotri(lcl + lvc ~ c(1, -0.01, 3.45))) + + mod <- function() { + ini({ + lka ~ 0.45 + lcl ~ 1 + lvc ~ 3.45 + }) + model({ + ka <- exp(lka) + cl <- exp(lcl) + vc <- exp(lvc) + kel <- cl / vc + d/dt(depot) <- -ka*depot + d/dt(central) <- ka*depot-kel*central + cp <- central / vc + }) + } + + mod2 <- mod |> ini(-lka) + + expect_equal(mod2$allCovs, "lka") + + +}) From a55f3114b2bc01e60789166f382752098cd52bce Mon Sep 17 00:00:00 2001 From: "Matthew L. Fidler" Date: Sat, 2 Dec 2023 22:13:18 -0600 Subject: [PATCH 3/5] Make sure eta values are not bounded --- R/piping-ini.R | 2 ++ 1 file changed, 2 insertions(+) diff --git a/R/piping-ini.R b/R/piping-ini.R index d0e7b14bf..9515ffc7e 100644 --- a/R/piping-ini.R +++ b/R/piping-ini.R @@ -477,6 +477,8 @@ .minfo("old initial estimate was negative, changing to positive") .newEta$est <- -.newEta$est } + .newEta$lower <- -Inf + .newEta$upper <- Inf .eta <- rbind(.eta, .newEta) } .ini <- rbind(.theta, .eta) From 53bb5e68d8bfeb66bd723d8c6fda599e6e409a47 Mon Sep 17 00:00:00 2001 From: "Matthew L. Fidler" Date: Sat, 2 Dec 2023 22:15:24 -0600 Subject: [PATCH 4/5] CF fixes --- R/piping-ini.R | 6 +++--- 1 file changed, 3 insertions(+), 3 deletions(-) diff --git a/R/piping-ini.R b/R/piping-ini.R index 9515ffc7e..16e6ec12f 100644 --- a/R/piping-ini.R +++ b/R/piping-ini.R @@ -455,7 +455,7 @@ .theta <- rbind(.theta, .newTheta) } else { # switch theta to eta - if (!is.na(.iniDf$err[.w])){ + if (!is.na(.iniDf$err[.w])) { stop("cannot switch error parameter '", .var, "' to a different type", call. = FALSE) } @@ -513,7 +513,7 @@ .eta$neta1 <- .eta$neta1 - ifelse(.eta$neta1 < .neta, 0L, 1L) .eta$neta2 <- .eta$neta2 - ifelse(.eta$neta2 < .neta, 0L, 1L) } else { - if (!is.na(.iniDf$err[.w])){ + if (!is.na(.iniDf$err[.w])) { stop("cannot switch error parameter '", .var, "' to a covariate", call. = FALSE) } @@ -586,7 +586,7 @@ .iniHandleFixOrUnfixEqual(expr=expr, rxui=rxui, envir=envir, maxLen=1L) } else if (.isTildeExpr(expr)) { .iniHandleSwitchType(expr=expr, rxui=rxui, envir=envir) - } else if (.isIniDropExpression(expr)){ + } else if (.isIniDropExpression(expr)) { .iniHandleDropType(expr=expr, rxui=rxui, envir=envir) } else { # Can this error be improved to clarify what is the expression causing the From 991e1f98ab6a811eb12cc4f02dd4a9f825d6f85d Mon Sep 17 00:00:00 2001 From: "Matthew L. Fidler" Date: Sun, 3 Dec 2023 10:02:06 -0600 Subject: [PATCH 5/5] Fix tests --- tests/testthat/test-piping-ini.R | 4 ++-- 1 file changed, 2 insertions(+), 2 deletions(-) diff --git a/tests/testthat/test-piping-ini.R b/tests/testthat/test-piping-ini.R index bf27620f7..b7f655d0a 100644 --- a/tests/testthat/test-piping-ini.R +++ b/tests/testthat/test-piping-ini.R @@ -679,11 +679,11 @@ test_that("change ini type with ~", { mod3 <- mod2 |> ini( ~ lvc) - expect_equal(mod3$omega, lotri(lcl ~ 1, lvc ~ 1)) + expect_equal(mod3$omega, NULL) mod4 <- mod3 |> ini( ~ lvc) - expect_equal(mod4$omega, lotri(lcl ~ 1)) + expect_equal(mod4$omega, lotri(lvc ~ 1)) })