version 0.0.1.0

DESCRIPTION

@@ -0,0 +1,18 @@
+Package: pmxcv
+Title: Integration-Based Coefficients of Variance
+Version: 0.0.1.0
+Authors@R: 
+    person("John", "Prybylski", , "john.prybylski@pfizer.com", role = c("aut", "cre"),
+           comment = c(ORCID = "0000-0001-5802-0539"))
+Description: Estimate coefficient of variance percent (CV%) for any arbitrary distribution, including some built-in estimates for commonly-used transformations in pharmacometrics. Methods are described in various sources, but applied here as summarized in: Prybylski, (2024) <doi:10.1007/s40262-023-01343-2>.
+License: MIT + file LICENSE
+Encoding: UTF-8
+RoxygenNote: 7.2.3
+Suggests: testthat (>= 3.0.0)
+Config/testthat/edition: 3
+NeedsCompilation: no
+Packaged: 2024-03-26 12:47:11 UTC; PRYBYJ
+Author: John Prybylski [aut, cre] (<https://orcid.org/0000-0001-5802-0539>)
+Maintainer: John Prybylski <john.prybylski@pfizer.com>
+Repository: CRAN
+Date/Publication: 2024-03-26 17:00:09 UTC

LICENSE

@@ -0,0 +1,2 @@
+YEAR: 2023
+COPYRIGHT HOLDER: pmxcv authors

MD5

@@ -0,0 +1,29 @@
+99aa8cb58d1e9d2b03fc0a77dc69547d *DESCRIPTION
+d6fde29bfc7eb42a6f8fa78f8ebfe486 *LICENSE
+6ab1c8faa1822829fd9ec29bc2a50492 *NAMESPACE
+fc4587d41b69c8762b67dc1b696d853c *R/dist.intcv.R
+5f54d3f9a3255f57c29b4c3bb1dba45f *R/dist.moment.R
+b3532e7a51fe2002c12bfecae44e0c36 *R/globals.R
+cef3fcaebd681a208eb6196b5730a066 *R/intcv.R
+b169640e90ff07c6ce2585d2ef525893 *R/invcv.R
+1e17e1da888585d104beab100f3d73a3 *R/moment.R
+07a5fc8c601450baefb9621a0fa1c676 *R/moment_f.R
+2b9247bad63afc49e58d6de060143866 *R/nonmemboxcox.R
+b59e7565c0d430be5d4acb0156623467 *R/numcv.R
+6f409bf021205b41698c2b4db0bffeb9 *README.md
+272c50dc213e362ed6a29be99ac57a3f *inst/CITATION
+8551d90b2b1ec8d5d7e0f36bcaf61140 *man/dist.intcv.Rd
+3fb2b2f0743dad71ef4783612c86e6b2 *man/dist.moment.Rd
+49203ca2c610fd4fa48baf42f29c980e *man/intcv.Rd
+25f4ba7db82d36735eae573397ffb802 *man/invcv.Rd
+4c2ee51123cf3caf5d6a4f7f15a02c35 *man/moment.Rd
+1b3026a9d9ed4f42eb0e9a5c1f0f52d8 *man/moment_f.Rd
+815d1e13670ebe92b4955da0ce5f31dc *man/nonmemboxcox.Rd
+093929504244eb91cb5ac3e33f31d739 *man/numcv.Rd
+6749f03a8c4544d7e7db128ef0c295cd *tests/testthat.R
+6433aa8d62c88111eb37d854d4fb6197 *tests/testthat/test-dist_builtins.R
+c97970652416447c68623a6f23d5e650 *tests/testthat/test-integral_cv.R
+1c86f21189db1dc3a5a5142964bb9d03 *tests/testthat/test-inverse_cv.R
+5ff1a5e7896c6b1be00893b57e3e7ca7 *tests/testthat/test-moments.R
+1e24aa2ea361ec4eadad77e26dfbb205 *tests/testthat/test-nm_boxcox.R
+9d9cfbc315a3e8b054d68f62218df058 *tests/testthat/test-numeric_cv.R

NAMESPACE

@@ -0,0 +1,14 @@
+# Generated by roxygen2: do not edit by hand
+
+export(dist.intcv)
+export(dist.moment)
+export(intcv)
+export(invcv)
+export(moment)
+export(moment_f)
+export(nonmemboxcox)
+export(numcv)
+importFrom(stats,optim)
+importFrom(stats,plogis)
+importFrom(stats,qlogis)
+importFrom(stats,sd)

R/dist.intcv.R

@@ -0,0 +1,37 @@
+#' Built-in integration-based %CV functions
+#'
+#' @param dist Selection of built-in distributions.
+#' @param ... passed to moment()
+#' @param exact If there is an exact moment generating function, use that. Default TRUE only for log
+#' @param lambda shape parameter for nonmemboxcox()
+#' @param fun return function (for use in invcv())
+#'
+#' @return Percent CV
+#' @export
+#'
+#' @importFrom stats plogis qlogis
+#'
+#' @examples
+dist.intcv <- function(dist="log",...,exact=ifelse(dist=="log",TRUE,FALSE),lambda=NULL,fun=FALSE) {
+  if (fun) return(function(...) dist.intcv(dist=dist, ...))
+  out <- NULL
+  if (dist=="log") {
+    list2env(list(...), environment())
+    if (exact) out <- 100*sqrt(exp(v) - 1)
+    else if ("u" %in% names(list(...))) out <- intcv(...,pdist = exp,qdist= log)
+    else out <- intcv(...,u=1,pdist = exp,qdist= log)
+  } else if (dist == "logexp") {
+    out <- intcv(...,pdist = function(x) exp(x) - 1,qdist= function(x) log(x+1))
+  } else if (dist == "logit") {
+    out <- intcv(...,pdist = plogis,qdist= qlogis)
+  } else if (dist == "arcsin") {
+    out <- intcv(...,pdist = function(x) sin(x)^2,qdist= function(x) asin(sqrt(x)))
+  } else if (dist == "nmboxcox") {
+    list2env(list(...), environment())
+    qbxcxt <- function(x) nonmemboxcox(x,lambda=lambda, theta = u)
+    pbxcxt <- function(x) nonmemboxcox(x,lambda=lambda, theta = u,inv=TRUE)
+    out <- intcv(...,pdist = pbxcxt,qdist= qbxcxt)
+  }
+  if (is.null(out)) stop(sprintf("Distribution %s not built into package.", dist))
+  out
+}

R/dist.moment.R

@@ -0,0 +1,34 @@
+#' Built-in moment functions
+#'
+#' @param dist Selection of built-in distributions.
+#' @param ... passed to moment()
+#' @param exact If there is an exact moment generating function, use that. Default TRUE only for log
+#' @param lambda shape parameter for nonmemboxcox()
+#'
+#' @return moment
+#' @export
+#'
+#' @importFrom stats plogis qlogis
+#'
+#' @examples
+dist.moment <- function(dist="log",...,exact=ifelse(dist=="log",TRUE,FALSE),lambda=NULL) {
+  out <- NULL
+  if (dist=="log") {
+    list2env(list(...), environment())
+    if (exact) out <- exp(n*log(u) + (n^2)*v/2) # u expected for exact
+    else out <- moment(...,pdist = exp,qdist= log)
+  } else if (dist == "logexp") {
+    out <- moment(...,pdist = function(x) exp(x) - 1,qdist= function(x) log(x+1))
+  } else if (dist == "logit") {
+    out <- moment(...,pdist = plogis,qdist= qlogis)
+  } else if (dist == "arcsin") {
+    out <- moment(...,pdist = function(x) sin(x)^2,qdist= function(x) asin(sqrt(x)))
+  } else if (dist == "nmboxcox") {
+    list2env(list(...), environment())
+    qbxcxt <- function(x) nonmemboxcox(x,lambda=lambda, theta = u)
+    pbxcxt <- function(x) nonmemboxcox(x,lambda=lambda, theta = u,inv=TRUE)
+    out <- moment(...,pdist = pbxcxt,qdist= qbxcxt)
+  }
+  if (is.null(out)) stop(sprintf("Distribution %s not built into package.", dist))
+  out
+}

R/globals.R

@@ -0,0 +1 @@
+utils::globalVariables(c("n", "u", "v"))

R/intcv.R

@@ -0,0 +1,13 @@
+#' Integration-based CV%
+#'
+#' @param ... Arguments passed to moment()
+#'
+#' @return Percent CV
+#' @export
+#'
+#' @examples
+intcv <- function(...) {
+  mom1 <- moment(n=1,...)
+  mom2 <- moment(n=2,...)
+  100*sqrt( mom2/mom1^2 - 1 )
+}

R/invcv.R

@@ -0,0 +1,31 @@
+#' Variance from CV%
+#'
+#' @param cvfun intcv()-based function
+#' @param cv CV% generated from cvfun
+#' @param verbose extra output
+#' @param ... Other parameters to pass to cvfun
+#'
+#' @return Best-fit variance
+#' @export
+#'
+#' @importFrom stats optim
+#'
+#' @examples
+invcv <- function(cvfun,cv,verbose=FALSE,...) {
+  # given cv function and cv, find omega
+  ssq <- function(x) {
+    expcv <- cvfun(...,v=x)
+    out <- (expcv - cv)^2
+    out[!is.finite(out)] <- 1e6
+    out
+  }
+  # Starting value should assume lognormal
+  start <- log((cv/100)^2 + 1)
+  if (ssq(start)<.Machine$double.eps) { # check if "exact"
+    if (verbose) message("Solved numerically as lognormal")
+    return(start)
+  }
+  opt <- optim(start,ssq,lower=0,method="L-BFGS-B") # L-BFGS-needed for bounds
+  if (verbose) message(opt$message)
+  opt$par
+}

R/moment.R

@@ -0,0 +1,12 @@
+#' Moment function
+#'
+#' @param ... all arguments passed to moment_f()
+#'
+#' @return moment
+#' @export
+#'
+#' @examples
+moment <- function(...) {
+  int_f <- function(x) moment_f(x,...)
+  stats::integrate(int_f, -Inf, Inf, abs.tol = 0)$value
+}

R/moment_f.R

@@ -0,0 +1,23 @@
+#' Integratable moment function
+#'
+#' @param x numeric vector
+#' @param u mean
+#' @param v variance
+#' @param n moment number
+#' @param pdist un-transform function for transformed random variable (eg, exp())
+#' @param qdist transform function (eg, log())
+#'
+#' @return Point result of the moment function
+#' @export
+#'
+#' @examples
+moment_f <- function(x,u,v,n,pdist,qdist) {
+  mu <- qdist(u) # u provided in natural units
+  sigma <- sqrt(v)
+
+  dens <-  stats::dnorm(x, mean = mu, sd = sigma)
+  p <- pdist(x)
+  o <- (p^n)*dens
+  o[!is.finite(o)] <- 0
+  o
+}

R/nonmemboxcox.R

@@ -0,0 +1,34 @@
+#' Box-Cox transform typically used in NONMEM
+#'
+#' Parameters are typically treated as lognormally-distributed by NONMEM users.
+#' Box-Cox transforms are typically applied to the exponentiated individual ETA parameters;
+#' this means the parameter is neither Box-Cox distributed nor lognormally-distributed, but both.
+#' To get the "Box-Cox Transform" as it would be relevant for CV% calculation, these properties
+#' have to be considered.
+#'
+#' @param x random vector. Must be positive.
+#' @param lambda shape parameter
+#' @param theta centrality parameter
+#' @param inv inverse transform
+#'
+#' @return Box-Cox transformed or untransformed vector
+#' @export
+#'
+#' @examples
+nonmemboxcox <- function(x,lambda,theta=1,inv=FALSE) {
+  if (!inv) {
+    normalized <- x/theta
+    {
+      if (lambda==0) bc1 <- log(normalized)
+      else bc1 <- ((normalized^lambda) - 1)/lambda
+    }
+    return( log(theta) + bc1 )
+  } else {
+    normalized <- x - log(theta)
+    {
+      if (lambda==0) bc1 <- exp(normalized)
+      else bc1 <- (normalized*lambda + 1)^(1/lambda)
+    }
+    return( theta*bc1 )
+  }
+}

R/numcv.R

@@ -0,0 +1,17 @@
+#' Numeric CV% of a sample
+#'
+#' @param x numeric vector
+#' @param ... other arguments for sd() and mean()
+#'
+#' @return Percent cv
+#' @export
+#'
+#' @importFrom stats sd
+#'
+#' @examples
+#' test_x <- rnorm(1000, mean=50, sd=5)
+#' cv <- numcv(test_x)
+#' cv # expect ~ 10(%)
+numcv <- function(x,...) {
+  100*sd(x,...)/mean(x,...)
+}

README.md

@@ -0,0 +1,55 @@
+
+<!-- README.md is generated from README.Rmd. Please edit that file -->
+
+# pmxcv
+
+<!-- badges: start -->
+<!-- badges: end -->
+
+This package is intended to provide easy access to methods of reporting
+CV% in non-lognormal, non-normal distributions. Often in pharmacometric
+literature, CV% is only reported for lognormal variability (or treated
+as lognormal), favoring reporting of the harder-to-interpret variance
+parameter; this package attempts to provide an alternative approach.
+
+## Installation
+
+You can install the development version of pmxcv like so:
+
+``` r
+devtools::install_github("pfizer-rd/pmxcv")
+```
+
+## Example
+
+Example of a typical use case:
+
+``` r
+library(pmxcv)
+
+## Parameters from NONMEM (etc)
+theta_bio <- 0.689
+omega_bio <- 1.2
+
+# Representation of bioavailability calculation in NONMEM syntax:L
+# LOGITF1 = LOG( THETA(BIO) ) - LOG( 1 - THETA(BIO) )
+# F1 = 1/( 1 + 1/EXP( LOGITF1 + ETA(BIO) )
+
+## Numeric parameter variability
+eta_sample <- rnorm(10^9, sd=sqrt(omega_bio))
+logit_theta <- log(theta_bio) - log( 1 - theta_bio )
+indiv_bios <- 1 / ( 1 + 1/exp( logit_theta + eta_sample ) )
+expected_cv <- 100*sd(indiv_bios)/mean(indiv_bios)
+expected_cv
+#> [1] 31.39837
+
+## Lognormal reported CV% (erroneous)
+bio_cv_lnorm <- 100*sqrt(exp(omega_bio) - 1)
+bio_cv_lnorm
+#> [1] 152.3193
+
+## Logitnormal reported CV%
+bio_cv_lgtnorm <- dist.intcv("logit", u=theta_bio, v=omega_bio)
+bio_cv_lgtnorm # should be approximately equal to numeric
+#> [1] 31.39758
+```

inst/CITATION

@@ -0,0 +1,8 @@
+bibentry(
+  bibtype  = "Article",
+  title    = "Reporting Coefficient of Variation for Logit, Box-Cox and Other Non-log-normal Parameters",
+  author   = c(person("John", "Prybylski", , "john.prybylski@pfizer.com", role = c("aut", "cre"), comment = c(ORCID = "0000-0001-5802-0539"))),
+  journal  = "Clinical Pharmacokinetics",
+  year     = 2024,
+  doi      = "10.1007/s40262-023-01343-2"
+)