R/get_hf.R

In nlmixr2autoinit: Automatic Generation of Initial Estimates for Population Pharmacokinetic Modeling

#' Estimate half-life from pooled pharmacokinetic data
#'
#' Estimates the terminal half-life of a drug using pooled pharmacokinetic data.
#' The method supports analysis based on first-dose, repeated-dose, or combined
#' dosing profiles.
#'
#' @param dat A data frame containing raw time–concentration data in the
#'   standard nlmixr2 format.
#'
#' @param dose_type Specifies the dosing context of the pharmacokinetic
#'   observations. Classified as:
#'   \itemize{
#'     \item first_dose: observations occur only after the initial administration
#'     \item repeated_doses: observations occur only after multiple
#'           administrations or at steady state
#'     \item combined_doses: observations include both first-dose and
#'           repeated-dose intervals
#'   }
#'
#' @param pooled Optional pooled data object generated by get_pooled_data. If
#'   not provided, pooled data are automatically created using the input dataset.
#'
#' @param verbose Logical; if TRUE (default), progress and completion
#'   messages are printed to the console. Set to FALSE to suppress all
#'   informational messages, for example when running automated scripts or tests.

#' @param ... Additional arguments passed to bin.time for pooling operations or
#'   to find_best_lambdaz for elimination slope estimation.
#'
#' @details
#' The function estimates terminal half-life using the following procedure:
#'   - Generate or use existing pooled pharmacokinetic data
#'   - Identify the terminal phase using elimination slope estimation
#'   - Calculate the terminal elimination rate constant (lambda_z)
#'   - Derive half-life using:
#'
#' \deqn{
#'   t_{1/2} = \frac{\log(2)}{\lambda_z}
#' }
#'
#' @return A list containing individual and median half-life estimates for
#'   first-dose, repeated-dose, and combined dosing profiles.
#'
#' @author Zhonghui Huang
#'
#' @examples
#' dat <- Bolus_1CPT
#' dat <- processData(dat)$dat
#' get_hf(dat, dose_type = "combined_doses")
#'
#' @seealso \link{get_pooled_data}, \link{bin.time}, \link{find_best_lambdaz}
#' @export
#'
get_hf <- function(dat,
                   dose_type = "first_dose",
                   pooled = NULL,
                   verbose=TRUE,
                   ...) {

  if (verbose){
  message(crayon::black(paste0("Estimating half-life", strrep(".", 20))))
  }
  dots <- list(...)
  bin_args <- dots[names(dots) %in% names(formals(bin.time))]
  slope_args <-
    dots[names(dots) %in% names(formals(find_best_lambdaz))]

  if (is.null(pooled)) {
    pooled <- do.call(get_pooled_data,
                      c(list(
                        dat = dat, dose_type = dose_type
                      ), bin_args))
  }

  # Initialize results
  half_life_fd <- NA
  half_life_md <- NA
  half_life_all <- NA
  slope_results_fd <- NA
  slope_results_md <- NA
  slope_results_all <- NA

  # First dose
  if (!is.null(pooled$datpooled_fd) &&
      "binned.df" %in% names(pooled$datpooled_fd)) {
    slope_results_fd <- do.call(force_find_lambdaz, c(
      list(
        time = pooled$datpooled_fd$binned.df$Time,
        conc = pooled$datpooled_fd$binned.df$Conc
      ),
      slope_args
    ))
    ke <- slope_results_fd$lambdaz
    half_life_fd <- ifelse(ke > 0, log(2) / ke, NA)
  }

  # Repeated doses
  if (!is.null(pooled$datpooled_efd) &&
      "binned.df" %in% names(pooled$datpooled_efd)) {
    slope_results_md <- do.call(force_find_lambdaz, c(
      list(
        time = pooled$datpooled_efd$binned.df$Time,
        conc = pooled$datpooled_efd$binned.df$Conc
      ),
      slope_args
    ))
    ke <- slope_results_md$lambdaz
    half_life_md <- ifelse(ke > 0, log(2) / ke, NA)
  }

  # Combined
  if (!is.null(pooled$datpooled_all) &&
      "binned.df" %in% names(pooled$datpooled_all)) {
    slope_results_all <-  do.call(force_find_lambdaz, c(
      list(
        time = pooled$datpooled_all$binned.df$Time,
        conc = pooled$datpooled_all$binned.df$Conc
      ),
      slope_args
    ))
    ke <- slope_results_all$lambdaz
    half_life_all <- ifelse(ke > 0, log(2) / ke, NA)
  }

  # Summarize
  half_life_values <- c(half_life_fd, half_life_md, half_life_all)
  positive_values <-
    half_life_values[half_life_values > 0 &
                       !is.na(half_life_values)]
  half_life_median <-
    if (length(positive_values) > 0)
      round(median(positive_values), 2)
  else
    NA

  # Final message
  message_text <- paste0(
    crayon::green("Half-life estimation complete: "),
    "Estimated t1/2 = ", half_life_median, " h"
  )

  if (verbose){
  message(message_text)
  }

  return(
    list(
      half_life_median = half_life_median,
      half_life_fd = half_life_fd,
      half_life_md = half_life_md,
      half_life_all = half_life_all,
      slope_results_fd = slope_results_fd,
      slope_results_md = slope_results_md,
      slope_results_all = slope_results_all
    )
  )
}