elimFit: Calculate the elimination rate(s) for a concentration-time...

In shirewoman2/LaurasHelpers: Miscellaneous Helper Functions for Common PK or Analytical Chemistry Tasks

Calculate the elimination rate(s) for a concentration-time curve

Description

elimFit fits an exponential decay equation to concentration-time data. The equation, which is of the form f(t) = concentration = A * exp(-kt) where A is ~Cmax, k is the elimination rate constant, and t is time. The equation can be monoexponential, biexponential, or triexponential.

Usage

elimFit(
  DF,
  concentration = Concentration,
  time = Time,
  startValues = NA,
  omit = NA,
  tmax = NA,
  modelType = "monoexponential",
  returnDataUsed = FALSE,
  weights = NULL,
  returnRSS = FALSE,
  useNLS_outnames = TRUE,
  maxiter = 50,
  graph = FALSE
)

Arguments

DF	The data.frame with the concentration-time data
concentration	The name of the column that contains concentration data
time	The name of the column that contains time data
startValues	The starting values to be used in the fit. Options: NA If left as NA, the starting values will be determined automatically. a list A list of starting values for the fit. Use the same values that you'd supply for nls. For a monoexponential fit, this will be a list of A and k. For a biexponential fit: A, alpha, B, beta. And for a triexponential fit: A, alpha, B, beta, G, gamma. a data.frame Sometimes, especially with sparsely sampled time points, the nls algorithm can fail to converge if you don't have perfect estimates for the starting values. When that is the case, set startValues to a two row data.frame with columns for each coefficient. The first row should contain the minimum of the range to search, and the second row should contain the maximum of the range to search for that coefficient. This option uses nls2, which does a more rigorous search for starting values than nls, to perform the regression. (For more information, see the documentation on nls2 and the algorithm "random-search".) A warning: Because nls2 searches more possible starting values, it can be appreciably slower than nls. A piece of advice: Regardless of whether you choose to use nls (supply no starting values or supply a list) or nls2 (supply a data.frame), you truly will benefit by supplying reasonable start values. Even if you're supplying a data.frame for nls2 to search, those values should still be reasonable or you just won't randomly select enough decent starting places to come up with regressions that will converge.
omit	An index of which, if any, samples to omit from the regression. These samples will be depicted as red open circles in the graph, if you choose to make one, but will not be included in the regression. Note that only points after tmax are used in the regression anyway, so there's no need to omit, e.g., t0.
tmax	The putative time at which the maximum concentration is observed, the time at which drug elimination becomes the dominant process. If left as NA, this will be whatever time is tmax. If you only want to fit a subset of times that don't necessarily include the actual tmax, set this to the time you want to start fitting. Time points before tmax will be omitted from the fit.
modelType	The mathematical model to use for fitting the data; options are "monoexponential", "biexponential", or "triexponential".
returnDataUsed	Should the data used be returned? I wrote this script initially for bootstrapping, where it can be useful to see what data were used as input. For that reason, I'm including the option of returning the data that were used.
weights	Weighting scheme to use for the regression. User may supply a numeric vector of weights to use or choose from "1/x", "1/x^2", "1/y" or "1/y^2". If left as NULL, no weighting scheme will be used. Be careful that you don't have any infinite values or this will fail!
returnRSS	TRUE or FALSE for whether to resturn the residual sum of squares. If set to TRUE, this will be the last column of the output data.frame where all rows = the residual sum of squares. (I wanted the output to still be a data.frame, so that's the place I could think of to put it.)
useNLS_outnames	TRUE or FALSE for whether to use the standard output coeffecient names that come with the nls or nls2 functions, e.g., "Estimate", "Std. Error", "t value", and "Pr(>|t|)". These names are annoying to work with for output data.frames b/c they don't follow standard column-naming practices (they contain spaces and symbols). If set to FALSE, the names of the output coefficient data.frame will be "Estimate", "SE", "tvalue" and "pvalue".
maxiter	Maximum number of iterations of start values to use; default is 50, just like nls. (See also nls.control.) Using more iterations means more random sampling of starting values and thus a higher likelihood of the fit converging. However, it also means more processing time.
graph	TRUE or FALSE for whether to create a graph of the data

Value

Returns a data.frame of the coefficients or returns a list containing whatever combination the user has specified of:

DataUsed

A data.frame of the input data

Estimates

A data.frame of the estimated coefficients

Graph

A ggplot2 graph of the input data with a line showing the fit to the terminal phase data

Examples

# Example data to work with:
data(ConcTime)
IV1 <- dplyr::filter(ConcTime, SubjectID == 101 & DoseRoute == "IV" &
                                     Drug == "A")
IV1 <- dplyr::select(IV1, SubjectID, TimeHr, Concentration)

# Automatically select the start values for the regression and start
# fitting at whatever time is tmax
elimFit(IV1, concentration = Concentration, time = TimeHr,
            modelType = "monoexponential")

# Set the start values yourself using a list (algorithm uses nls
# to fit the data).
elimFit(IV1, concentration = Concentration, time = TimeHr,
            startValues = list(A = 30, k = 0.01),
            modelType = "monoexponential")

# Set the start values yourself but use the more robust nls2
# function to do the regression. Provide a range of values to search.
elimFit(IV1, concentration = Concentration, time = TimeHr,
            startValues = data.frame(A = c(5, 50), k = c(0.0001, 0.05)),
            modelType = "monoexponential")

# Omit a point. In this case, omit the point at t = 8.
elimFit(IV1, concentration = Concentration, time = TimeHr,
            modelType = "monoexponential",
            omit = which(IV1$TimeHr == 8))

# Don't start fitting until a later time than tmax, e.g., t = 4.
elimFit(IV1, concentration = Concentration, time = TimeHr,
            modelType = "monoexponential", tmax = 4))

# Weight by 1/y.
elimFit(IV1, concentration = Concentration, time = TimeHr,
            weight = 1/IV1$Concentration,
            modelType = "monoexponential")

# Another way to weight by 1/y
elimFit(IV1, concentration = Concentration, time = TimeHr,
            weight = "1/y",
            modelType = "monoexponential")

# Get the residual sum of squares
elimFit(IV1, concentration = Concentration, time = TimeHr,
            modelType = "monoexponential", returnRSS = TRUE)

# Use better names for the columns in the output
elimFit(IV1, concentration = Concentration, time = TimeHr,
            modelType = "monoexponential", useNLS_outnames = FALSE)

# Graph the data; good for visually inspecting the fit.
elimFit(IV1, concentration = Concentration, time = TimeHr,
            modelType = "monoexponential", tmax = 4, graph = TRUE)

# Some data for using the triexponential fit
TriExpDF <- data.frame(Time_min = c(0, 18, 33, 48, 63, 95, 123, 153, 183,
                                    213, 243, 303, 483),
                       Concentration = c(0, 420, 228, 143, 90, 48, 28,
                                         18, 13, 10, 7, 4, 2))

# Some starting values to use for the triexponential fit
Start_tri <- data.frame(A = c(100, 500),
                        alpha = c(0.01, 0.5),
                        B = c(50, 200),
                        beta = c(0.001, 0.05),
                        G = c(1, 100),
                        gamma = c(1e-04, 0.01))

elimFit(TriExpDF, startValues = Start_tri, concentration = Concentration,
        time = Time_min, tmax = 33, modelType = "triexponential",
        graph = TRUE)

shirewoman2/LaurasHelpers documentation built on Oct. 22, 2023, 2:07 p.m.