From 70f8f86d62a5eacfb004ca39f87af6e3d501ba89 Mon Sep 17 00:00:00 2001 From: "Matthew L. Fidler" Date: Fri, 29 Sep 2023 21:20:53 -0500 Subject: [PATCH 1/5] na.rm more often to help avoid situations like #59 --- R/computingutil.R | 77 ++++++++++++++++++++++++----------------------- 1 file changed, 39 insertions(+), 38 deletions(-) diff --git a/R/computingutil.R b/R/computingutil.R index 888fbd1..4a29a3b 100644 --- a/R/computingutil.R +++ b/R/computingutil.R @@ -29,7 +29,7 @@ # mean by groups (Individual) groupMeans <- with(data, ave(get(covariate),get(uidCol), FUN = function(x) mean(x, na.rm = TRUE))) # pop mean - popMean <- mean(groupMeans) + popMean <- mean(groupMeans, na.rm=TRUE) # pop std popStd <- .sd.p (groupMeans) @@ -171,7 +171,8 @@ foldgen <- function(data,nfold=5,stratVar=NULL){ unique( quantile(y, probs = - seq(0, 1, length = cuts))), + seq(0, 1, length = cuts), + na.rm=TRUE)), include.lowest = TRUE) } @@ -235,8 +236,8 @@ optimUnisampling <- function(xvec,N=1000,medValue,floorT=TRUE) { if (floorT){ x <- floor(stats::runif(N, xmin, xmax))} else{ - x <- stats::runif(N, xmin, xmax) - } + x <- stats::runif(N, xmin, xmax) + } xdist <- (median(x)-medValue)^2 xdist } @@ -518,7 +519,7 @@ bootstrapFit <- function(fit, xPosthoc <- nlmixr2(x, data = origData, est = "posthoc", control = list(calcTables = FALSE, print = 1, compress=FALSE) - ) + ) saveRDS(xPosthoc, .path) } xPosthoc$objf - fit$objf @@ -617,7 +618,7 @@ sampling <- function(data, len = 1, any.missing = FALSE, lower = 2 - ) + ) } if (performStrat && missing(stratVar)) { @@ -629,7 +630,7 @@ sampling <- function(data, lower = 2, len = 1, any.missing = FALSE - ) + ) if (missing(uid_colname)) { # search the dataframe for a column name of 'ID' @@ -761,7 +762,7 @@ modelBootstrap <- function(fit, len = 1, any.missing = FALSE, lower = 1 - ) + ) if (missing(nSampIndiv)) { nSampIndiv <- length(unique(data[, uidCol])) @@ -806,7 +807,7 @@ modelBootstrap <- function(fit, fnameBootDataPattern <- paste0("boot_data_", "[0-9]+", ".rds", sep = "" - ) + ) fileExists <- list.files(paste0("./", output_dir), pattern = fnameBootDataPattern) @@ -886,7 +887,7 @@ modelBootstrap <- function(fit, if (!restart) { if (length(modFileExists) > 0 && - (length(fileExists) > 0)) { + (length(fileExists) > 0)) { # read bootData and modelsEnsemble files from disk cli::cli_alert_success( @@ -960,28 +961,28 @@ modelBootstrap <- function(fit, modIdx)) fit <- tryCatch( - { - fit <- suppressWarnings(nlmixr2(ui, - boot_data, - est = fitMeth, - control = .ctl)) - - .env$multipleFits <- list( - # objf = fit$OBJF, - # aic = fit$AIC, - omega = fit$omega, - parFixedDf = fit$parFixedDf[, c("Estimate", "Back-transformed")], - message = fit$message, - warnings = fit$warnings) - - fit # to return 'fit' - }, - error = function(error_message) { - message("error fitting the model") - message(error_message) - message("storing the models as NA ...") - return(NA) # return NA otherwise (instead of NULL) - }) + { + fit <- suppressWarnings(nlmixr2(ui, + boot_data, + est = fitMeth, + control = .ctl)) + + .env$multipleFits <- list( + # objf = fit$OBJF, + # aic = fit$AIC, + omega = fit$omega, + parFixedDf = fit$parFixedDf[, c("Estimate", "Back-transformed")], + message = fit$message, + warnings = fit$warnings) + + fit # to return 'fit' + }, + error = function(error_message) { + message("error fitting the model") + message(error_message) + message("storing the models as NA ...") + return(NA) # return NA otherwise (instead of NULL) + }) saveRDS( .env$multipleFits, @@ -1070,7 +1071,7 @@ extractVars <- function(fitlist, id = "method") { if (!(id == "omega" || - id == "parFixedDf")) { + id == "parFixedDf")) { # check if all message strings are empty if (id == "message") { prev <- TRUE @@ -1136,11 +1137,11 @@ getBootstrapSummary <- function(fitList, # omega estimates omegaMatlist <- extractVars(fitList, id) varVec <- simplify2array(omegaMatlist) - mn <- apply(varVec, 1:2, mean) - sd <- apply(varVec, 1:2, sd) + mn <- apply(varVec, 1:2, mean, na.rm=TRUE) + sd <- apply(varVec, 1:2, sd, na.rm=TRUE) quants <- apply(varVec, 1:2, function(x) { - unname(quantile(x, quantLevels)) + unname(quantile(x, quantLevels, na.rm=TRUE)) }) median <- quants[1, , ] confLower <- quants[2, , ] @@ -1412,11 +1413,11 @@ bootplot.nlmixr2FitCore <- function(x, ...) { } else { stop("this nlmixr2 object does not include boostrap distribution statics for comparison", call. = FALSE - ) + ) } } else { stop("this is not a nlmixr2 object", call. = FALSE - ) + ) } } From 7092131e4b6bf67d2b9ad339af5582eee2118083 Mon Sep 17 00:00:00 2001 From: "Matthew L. Fidler" Date: Fri, 29 Sep 2023 21:25:47 -0500 Subject: [PATCH 2/5] add news about na.rm --- NEWS.md | 11 ++++++++--- 1 file changed, 8 insertions(+), 3 deletions(-) diff --git a/NEWS.md b/NEWS.md index c0a9ef8..58993da 100644 --- a/NEWS.md +++ b/NEWS.md @@ -1,3 +1,8 @@ +# nlmixr2extra development version + +* `bootstrapFit()` now will be more careful handling `NA` values so + they do not completely affect results (Issue #59) + # nlmixr2extra 2.0.9 * New method for `knit_print()` will generate model equations for LaTeX @@ -7,7 +12,7 @@ * Use `assignInMyNamespace()` instead of using the global assignment operator for the horseshoe prior - + * Be specific in version requirements (as requested by CRAN checks) * Move the `theoFitOde.rda` data build to `devtools::document()` to @@ -19,10 +24,10 @@ * Fix `cli` issues with the new `cli` 3.4+ release that will allow bootstrapping to run again (before `cli` would error, this fixes the `donttest` issues on CRAN). - + * Fixed step-wise covariate selection to work a bit better with the updated UI, thanks to Vishal Sarsani - + * Added lasso covariate selection (thanks to Vishal Sarsani) * Added horseshoe prior covarite selecion (thanks to Vishal Sarsani) From 5686d521326705aad5ed338dcd803e1fda2187a1 Mon Sep 17 00:00:00 2001 From: "Matthew L. Fidler" Date: Fri, 29 Sep 2023 21:27:14 -0500 Subject: [PATCH 3/5] CF fix; changed to seq_along() --- R/computingutil.R | 2 +- 1 file changed, 1 insertion(+), 1 deletion(-) diff --git a/R/computingutil.R b/R/computingutil.R index 4a29a3b..37da6df 100644 --- a/R/computingutil.R +++ b/R/computingutil.R @@ -1163,7 +1163,7 @@ getBootstrapSummary <- function(fitList, parFixedlistVec <- do.call("rbind", parFixedlistVec) omgVecBoot <- list() - omegaIdx <- seq(length(omegaMatlist)) + omegaIdx <- seq_along(omegaMatlist) omgVecBoot <- lapply(omegaIdx, function(idx) { omgMat <- omegaMatlist[[idx]] From ce18e3ba01dcb18f0a5b3b0e48c73d454d4abdc8 Mon Sep 17 00:00:00 2001 From: "Matthew L. Fidler" Date: Fri, 29 Sep 2023 21:49:27 -0500 Subject: [PATCH 4/5] Only take cov2cor with non-zero diagonals --- R/computingutil.R | 10 +++++++++- 1 file changed, 9 insertions(+), 1 deletion(-) diff --git a/R/computingutil.R b/R/computingutil.R index 37da6df..256dc3c 100644 --- a/R/computingutil.R +++ b/R/computingutil.R @@ -1202,8 +1202,16 @@ getBootstrapSummary <- function(fitList, parFixedOmegaCombined <- cbind(parFixedlistVec, omgVecBoot) covMatrix <- cov(parFixedOmegaCombined) - corMatrix <- cov2cor(covMatrix) + w <- which(diag(covMatrix) == 0) + if (length(w) > 0) { + d <- dim(covMatrix)[1] + corMatrix <- matrix(rep(0,d * d), d, d) + corMatrix[-w, -w] <- cov2cor(covMatrix[-w, -w]) + } else { + corMatrix <- cov2cor(covMatrix) + } diag(corMatrix) <- sqrt(diag(covMatrix)) + dimnames(corMatrix) <- dimnames(covMatrix) lst <- list( mean = mn, median = median, From 91d57724a97fbd6c881a6895d3b25cd826332c6d Mon Sep 17 00:00:00 2001 From: "Matthew L. Fidler" Date: Fri, 29 Sep 2023 21:50:43 -0500 Subject: [PATCH 5/5] Add to news --- NEWS.md | 3 +++ 1 file changed, 3 insertions(+) diff --git a/NEWS.md b/NEWS.md index 58993da..c1c88ac 100644 --- a/NEWS.md +++ b/NEWS.md @@ -3,6 +3,9 @@ * `bootstrapFit()` now will be more careful handling `NA` values so they do not completely affect results (Issue #59) +* `bootstrapFit()` will now only take the correlation of the non-zero + diagonals (Issue #59). + # nlmixr2extra 2.0.9 * New method for `knit_print()` will generate model equations for LaTeX