Skip to content

Commit

Permalink
Merge pull request #105 from keaven/almostIntegerIssue101
Browse files Browse the repository at this point in the history
Add `gsSurvCalendar()`, fix `toInteger()`, and fix `toBinomialExact()`
  • Loading branch information
nanxstats authored Nov 12, 2023
2 parents 3305b2f + 4bc7749 commit 6de6583
Show file tree
Hide file tree
Showing 23 changed files with 877 additions and 203 deletions.
2 changes: 1 addition & 1 deletion DESCRIPTION
Original file line number Diff line number Diff line change
@@ -1,5 +1,5 @@
Package: gsDesign
Version: 3.5.0
Version: 3.6.0
Title: Group Sequential Design
Authors@R: person(given = "Keaven", family = "Anderson", email =
"[email protected]", role = c('aut','cre'))
Expand Down
1 change: 1 addition & 0 deletions NAMESPACE
Original file line number Diff line number Diff line change
Expand Up @@ -51,6 +51,7 @@ export(gsPosterior)
export(gsProbability)
export(gsRR)
export(gsSurv)
export(gsSurvCalendar)
export(hGraph)
export(hrn2z)
export(hrz2n)
Expand Down
27 changes: 27 additions & 0 deletions NEWS.md
Original file line number Diff line number Diff line change
@@ -1,3 +1,30 @@
# gsDesign 3.6.0 (November, 2023)

## Breaking changes

- `gsSurv()` and `nSurv()` have updated default values for `T` and `minfup`
so that function calls with no arguments will run. Legacy code with `T`
or `minfup` not explicitly specified could break (#105).

## New features

- `gsSurvCalendar()` function added to enable group sequential design for
time-to-event outcomes using calendar timing of interim analysis
specification (#105).
- `as_rtf()` method for `gsBinomialExact()` objects added,
enable RTF table outputs for standard word processing software (#102).

## Improvements

- `toBinomialExact()` and `gsBinomialExact()`: fix error checking in bound
computations, improve documentation and error messages (#105).
- `print.gsSurv()`: Improve the display of targeted events (very minor).
The boundary crossing probability computations did not change.
The need is made evident by the addition of the `toInteger()` function (#105).
- `toInteger()`: Fix the documentation and execution based on the `ratio` argument (#105).
- Update the vaccine efficacy, Poisson mixture model, and toInteger vignettes (#105).
- Standardize and improve roxygen2 documentation (#104).

# gsDesign 3.5.0 (July, 2023)

- `sfPower()` now allows a wider parameter range (0, 15].
Expand Down
2 changes: 1 addition & 1 deletion R/as_rtf.R
Original file line number Diff line number Diff line change
Expand Up @@ -21,7 +21,7 @@
#'
#' @examples
#' zz <- gsBinomialExact(
#' k = 3, theta = seq(0, 1, 0.1), n.I = c(12, 24, 36),
#' k = 3, theta = seq(0.1, 0.9, 0.1), n.I = c(12, 24, 36),
#' a = c(-1, 0, 11), b = c(5, 9, 12)
#' )
#' zz %>%
Expand Down
24 changes: 13 additions & 11 deletions R/gsBinomialExact.R
Original file line number Diff line number Diff line change
Expand Up @@ -81,10 +81,10 @@ utils::globalVariables(c("N", "EN", "Bound", "rr", "Percent", "Outcome"))
#' length k.
#' @param a Number of "successes" required to cross lower bound cutoffs to
#' reject \code{p1} in favor of \code{p0} at each analysis; vector of length k;
#' -1 means no lower bound.
#' -1 (minimum allowed) means no lower bound.
#' @param b Number of "successes" required to cross upper bound cutoffs for
#' rejecting \code{p0} in favor of \code{p1} at each analysis; vector of length
#' k.
#' k; value > n.I means no upper bound.
#' @param p0 Lower of the two response (event) rates hypothesized.
#' @param p1 Higher of the two response (event) rates hypothesized.
#' @param alpha Nominal probability of rejecting response (event) rate
Expand Down Expand Up @@ -167,7 +167,7 @@ utils::globalVariables(c("N", "EN", "Bound", "rr", "Percent", "Outcome"))
#' library(ggplot2)
#'
#' zz <- gsBinomialExact(
#' k = 3, theta = seq(0, 1, 0.1), n.I = c(12, 24, 36),
#' k = 3, theta = seq(0.1, 0.9, 0.1), n.I = c(12, 24, 36),
#' a = c(-1, 0, 11), b = c(5, 9, 12)
#' )
#'
Expand Down Expand Up @@ -245,22 +245,24 @@ utils::globalVariables(c("N", "EN", "Bound", "rr", "Percent", "Outcome"))
# gsBinomialExact function [sinew] ----
gsBinomialExact <- function(k = 2, theta = c(.1, .2), n.I = c(50, 100), a = c(3, 7), b = c(20, 30)) {
checkScalar(k, "integer", c(2, Inf), inclusion = c(TRUE, FALSE))
checkVector(theta, "numeric", interval = 0:1, inclusion = c(TRUE, TRUE))
checkVector(theta, "numeric", interval = 0:1, inclusion = c(FALSE, FALSE))
checkVector(n.I, "integer", interval = c(1, Inf), inclusion = c(TRUE, FALSE))
checkVector(a, "integer", interval = c(-Inf, Inf), inclusion = c(FALSE, FALSE))
checkVector(a, "integer", interval = c(-1, Inf), inclusion = c(TRUE, FALSE))
checkVector(b, "integer", interval = c(1, Inf), inclusion = c(FALSE, FALSE))
ntheta <- as.integer(length(theta))
theta <- as.double(theta)
if (k != length(n.I) || k != length(a) || k != length(b)) {
stop("Lengths of n.I, a, and b must equal k on input")
}
m <- c(n.I[1], diff(n.I))
if (min(m) < 1) stop("n.I must must contain an increasing sequence of positive integers")
if (min(n.I - a) < 0) stop("Input a-vector must be less than n.I")
if (min(b - a) <= 0) stop("Input b-vector must be strictly greater than a")
if (min(diff(a)) < 0) stop("a must contain a non-decreasing sequence of integers")
if (min(diff(b)) < 0) stop("b must contain a non-decreasing sequence of integers")
if (min(diff(n.I - b)) < 0) stop("n.I - b must be non-decreasing")
if (min(m) < 1) stop(paste("n.I must contain an increasing sequence of positive integers\n n.I = ", paste(n.I, collapse = ", ")))
if (min(n.I - a) <= 0) stop(paste("Input a-vector must be strictly less than n.I\n a = ", paste(a, collapse = ", "), "\n n.I = ", paste(n.I, collapse = ", ")))
if (min(b - a) <= 0) stop(paste("Input b-vector must be strictly greater than a\n a = ", paste(a, collapse = ", "), "\n b = ", paste(b, collapse = ", ")))
if (min(diff(a)) < 0) stop(paste("a must contain a non-decreasing sequence of non-negative integers\n a =", paste(a, collapse = ", ")))
if (min(diff(b)) < 0) stop(paste("b must contain a non-decreasing sequence of non-negative integers\n b =", paste(b, collapse = ", ")))
if (min(diff(n.I - b)) < 0) stop(paste("n.I - b must be non-decreasing\n n.I =",
paste(n.I, collapse=", "), "\n b = ", paste(b, collapse=", "),
"\n n.I - b = ", paste(n.I - b, collapse=", ")))

plo <- matrix(nrow = k, ncol = ntheta)
rownames(plo) <- paste(rep("Analysis ", k), 1:k)
Expand Down
2 changes: 1 addition & 1 deletion R/gsMethods.R
Original file line number Diff line number Diff line change
Expand Up @@ -690,7 +690,7 @@ gsBoundSummary <- function(x, deltaname = NULL, logdelta = FALSE, Nname = NULL,
nstat <- 2
} else {
nstat <- 4
statframe[statframe$Value == statframe$Value[3], ]$Analysis <- paste("Events:", ceiling(rowSums(x$eDC + x$eDE)))
statframe[statframe$Value == statframe$Value[3], ]$Analysis <- paste("Events:", ceiling(x$n.I))
if (x$ratio == 1) N <- 2 * ceiling(rowSums(x$eNE)) else N <- ceiling(rowSums(x$eNE)) + ceiling(rowSums(x$eNC))
Time <- round(x$T, tdigits)
statframe[statframe$Value == statframe$Value[4], ]$Analysis <- paste(timename, ": ", as.character(Time), sep = "")
Expand Down
14 changes: 6 additions & 8 deletions R/gsSurv.R
Original file line number Diff line number Diff line change
Expand Up @@ -865,7 +865,7 @@ KT <- function(alpha = .025, sided = 1, beta = .1,
#' \code{r} will not be changed by the user.
#' @param usTime Default is NULL in which case upper bound spending time is
#' determined by \code{timing}. Otherwise, this should be a vector of length
#' code{k} with the spending time at each analysis (see Details in help for \code{gsDesign}).
#' \code{k} with the spending time at each analysis (see Details in help for \code{gsDesign}).
#' @param lsTime Default is NULL in which case lower bound spending time is
#' determined by \code{timing}. Otherwise, this should be a vector of length
#' \code{k} with the spending time at each analysis (see Details in help for \code{gsDesign}).
Expand Down Expand Up @@ -989,7 +989,7 @@ KT <- function(alpha = .025, sided = 1, beta = .1,
#' # piecewise constant enrollment rates (vary accrual duration)
#' nSurv(
#' lambdaC = log(2) / 6, hr = .5, eta = log(2) / 40, gamma = c(1, 3, 6),
#' R = c(3, 6, 9), minfup = 12
#' R = c(3, 6, 9), T = NULL, minfup = 12
#' )
#'
#' # stratified population (vary accrual duration)
Expand Down Expand Up @@ -1030,7 +1030,7 @@ KT <- function(alpha = .025, sided = 1, beta = .1,
#' @export
# nSurv function [sinew] ----
nSurv <- function(lambdaC = log(2) / 6, hr = .6, hr0 = 1, eta = 0, etaE = NULL,
gamma = 1, R = 12, S = NULL, T = NULL, minfup = NULL, ratio = 1,
gamma = 1, R = 12, S = NULL, T = 18, minfup = 6, ratio = 1,
alpha = 0.025, beta = 0.10, sided = 1, tol = .Machine$double.eps^0.25) {
if (is.null(etaE)) etaE <- eta
# set up rates as matrices with row and column names
Expand Down Expand Up @@ -1189,7 +1189,7 @@ gsSurv <- function(k = 3, test.type = 4, alpha = 0.025, sided = 1,
beta = 0.1, astar = 0, timing = 1, sfu = sfHSD, sfupar = -4,
sfl = sfHSD, sflpar = -2, r = 18,
lambdaC = log(2) / 6, hr = .6, hr0 = 1, eta = 0, etaE = NULL,
gamma = 1, R = 12, S = NULL, T = NULL, minfup = NULL, ratio = 1,
gamma = 1, R = 12, S = NULL, T = 18, minfup = 6, ratio = 1,
tol = .Machine$double.eps^0.25,
usTime = NULL, lsTime = NULL) # KA: last 2 arguments added 10/8/2017
{
Expand Down Expand Up @@ -1283,16 +1283,14 @@ print.gsSurv <- function(x, digits = 2, ...) {
print.gsDesign(x)
if (x$test.type != 1) {
y <- cbind(
x$T, (x$eNC + x$eNE) %*% rep(1, ncol(x$eNE)),
(x$eDC + x$eDE) %*% rep(1, ncol(x$eNE)),
x$T, (x$eNC + x$eNE) %*% rep(1, ncol(x$eNE)), x$n.I,
round(zn2hr(x$lower$bound, x$n.I, x$ratio, hr0 = x$hr0, hr1 = x$hr), 3),
round(zn2hr(x$upper$bound, x$n.I, x$ratio, hr0 = x$hr0, hr1 = x$hr), 3)
)
colnames(y) <- c("T", "n", "Events", "HR futility", "HR efficacy")
} else {
y <- cbind(
x$T, (x$eNC + x$eNE) %*% rep(1, ncol(x$eNE)),
(x$eDC + x$eDE) %*% rep(1, ncol(x$eNE)),
x$T, (x$eNC + x$eNE) %*% rep(1, ncol(x$eNE)), x$n.I,
round(zn2hr(x$upper$bound, x$n.I, x$ratio, hr0 = x$hr0, hr1 = x$hr), 3)
)
colnames(y) <- c("T", "n", "Events", "HR efficacy")
Expand Down
205 changes: 205 additions & 0 deletions R/gsSurvCalendar.R
Original file line number Diff line number Diff line change
@@ -0,0 +1,205 @@
#' Time-to-event endpoint design with calendar timing of analyses
#'
#' @param test.type Test type. See \code{\link{gsSurv}}.
#' @param alpha Type I error rate. Default is 0.025 since 1-sided
#' testing is default.
#' @param sided \code{1} for 1-sided testing, \code{2} for 2-sided testing.
#' @param beta Type II error rate. Default is 0.10
#' (90\% power); \code{NULL} if power is to be computed based on
#' other input values.
#' @param astar Normally not specified. If \code{test.type = 5}
#' or \code{6}, \code{astar} specifies the total probability
#' of crossing a lower bound at all analyses combined. This
#' will be changed to \code{1 - alpha} when default value of
#' \code{0} is used. Since this is the expected usage,
#' normally \code{astar} is not specified by the user.
#' @param sfu A spending function or a character string
#' indicating a boundary type (that is, \code{"WT"} for
#' Wang-Tsiatis bounds, \code{"OF"} for O'Brien-Fleming bounds and
#' \code{"Pocock"} for Pocock bounds). For one-sided and symmetric
#' two-sided testing is used to completely specify spending
#' (\code{test.type = 1}, \code{2}), \code{sfu}. The default value is
#' \code{sfHSD} which is a Hwang-Shih-DeCani spending function.
#' @param sfupar Real value, default is \code{-4} which is an
#' O'Brien-Fleming-like conservative bound when used with the
#' default Hwang-Shih-DeCani spending function. This is a
#' real-vector for many spending functions. The parameter
#' \code{sfupar} specifies any parameters needed for the spending
#' function specified by \code{sfu}; this will be ignored for
#' spending functions (\code{sfLDOF}, \code{sfLDPocock})
#' or bound types (\code{"OF"}, \code{"Pocock"})
#' that do not require parameters.
#' @param sfl Specifies the spending function for lower
#' boundary crossing probabilities when asymmetric,
#' two-sided testing is performed
#' (\code{test.type = 3}, \code{4}, \code{5}, or \code{6}).
#' Unlike the upper bound,
#' only spending functions are used to specify the lower bound.
#' The default value is \code{sfHSD} which is a
#' Hwang-Shih-DeCani spending function. The parameter
#' \code{sfl} is ignored for one-sided testing
#' (\code{test.type = 1}) or symmetric 2-sided testing
#' (\code{test.type = 2}).
#' @param sflpar Real value, default is \code{-2}, which, with the
#' default Hwang-Shih-DeCani spending function, specifies a
#' less conservative spending rate than the default for the
#' upper bound.
#' @param calendarTime Vector of increasing positive numbers
#' with calendar times of analyses. Time 0 is start of
#' randomization.
#' @param spending Select between calendar-based spending and
#' information-based spending.
#' @param lambdaC Scalar, vector or matrix of event hazard
#' rates for the control group; rows represent time periods while
#' columns represent strata; a vector implies a single stratum.
#' @param hr Hazard ratio (experimental/control) under the
#' alternate hypothesis (scalar).
#' @param hr0 Hazard ratio (experimental/control) under the null
#' hypothesis (scalar).
#' @param eta Scalar, vector or matrix of dropout hazard rates
#' for the control group; rows represent time periods while
#' columns represent strata; if entered as a scalar, rate is
#' constant across strata and time periods; if entered as a
#' vector, rates are constant across strata.
#' @param etaE Matrix dropout hazard rates for the experimental
#' group specified in like form as \code{eta}; if \code{NULL},
#' this is set equal to \code{eta}.
#' @param gamma A scalar, vector or matrix of rates of entry by
#' time period (rows) and strata (columns); if entered as a
#' scalar, rate is constant across strata and time periods;
#' if entered as a vector, rates are constant across strata.
#' @param R A scalar or vector of durations of time periods for
#' recruitment rates specified in rows of \code{gamma}. Length is the
#' same as number of rows in \code{gamma}. Note that when variable
#' enrollment duration is specified (input \code{T = NULL}), the final
#' enrollment period is extended as long as needed.
#' @param S A scalar or vector of durations of piecewise constant
#' event rates specified in rows of \code{lambda}, \code{eta} and \code{etaE};
#' this is \code{NULL} if there is a single event rate per stratum
#' (exponential failure) or length of the number of rows in \code{lambda}
#' minus 1, otherwise.
#' @param minfup A non-negative scalar less than the maximum value
#' in \code{calendarTime}. Enrollment will be cut off at the
#' difference between the maximum value in \code{calendarTime}
#' and \code{minfup}.
#' @param ratio Randomization ratio of experimental treatment
#' divided by control; normally a scalar, but may be a vector with
#' length equal to number of strata.
#' @param r Integer value controlling grid for numerical
#' integration as in Jennison and Turnbull (2000); default is 18,
#' range is 1 to 80. Larger values provide larger number of grid
#' points and greater accuracy. Normally \code{r} will not be changed by
#' the user.
#' @param tol Tolerance for error passed to the \code{\link{gsDesign}} function.
#'
#' @export
#'
#' @rdname gsSurvCalendar
#'
#' @examples
#' # First example: while timing is calendar-based, spending is event-based
#' x <- gsSurvCalendar() %>% toInteger()
#' gsBoundSummary(x)
#'
#' # Second example: both timing and spending are calendar-based
#' # This results in less spending at interims and leaves more for final analysis
#' y <- gsSurvCalendar(spending = "calendar") %>% toInteger()
#' gsBoundSummary(y)
#'
#' # Note that calendar timing for spending relates to planned timing for y
#' # rather than timing in y after toInteger() conversion
#'
#' # Values plugged into spending function for calendar time
#' y$usTime
#' # Actual calendar fraction from design after toInteger() conversion
#' y$T / max(y$T)
gsSurvCalendar <- function(
test.type = 4, alpha = 0.025, sided = 1, beta = 0.1, astar = 0,
sfu = gsDesign::sfHSD, sfupar = -4,
sfl = gsDesign::sfHSD, sflpar = -2,
calendarTime = c(12, 24, 36),
spending = c("information", "calendar"),
lambdaC = log(2) / 6, hr = .6, hr0 = 1, eta = 0, etaE = NULL,
gamma = 1, R = 12, S = NULL, minfup = 18, ratio = 1,
r = 18, tol = 1e-06) {
checkVector(calendarTime, isType = "numeric")
if (min(diff(calendarTime)) < 0) stop(paste("calendarTime must be an increasing vector of positive numbers, including final analysis\n input calendarTime = ",
paste(calendarTime, collapse = ", ")))
x <- nSurv(
lambdaC = lambdaC, hr = hr, hr0 = hr0, eta = eta, etaE = etaE,
gamma = gamma, R = R, S = S, T = max(calendarTime),
minfup = minfup, ratio = ratio,
alpha = alpha, beta = beta, sided = sided
)

# Get interim expected event counts and sample size based on
# input gamma, eta, lambdaC, R, S, minfup
eDC <- NULL
eDE <- NULL
eNC <- NULL
eNE <- NULL
T <- NULL
k <- length(calendarTime)
for (i in 1:k) {
xx <- nEventsIA(tIA = calendarTime[i], x = x, simple = FALSE)
eDC <- rbind(eDC, xx$eDC)
eDE <- rbind(eDE, xx$eDE)
eNC <- rbind(eNC, xx$eNC)
eNE <- rbind(eNE, xx$eNE)
}
timing <- rowSums(eDC) + rowSums(eDE)
timing <- timing / max(timing)
# if calendar spending, set usTime, lsTime
if (spending[1] == "calendar") {
lsTime <- calendarTime / max(calendarTime)
} else {
lsTime <- NULL
}
usTime <- lsTime

# Now inflate events to get targeted power
y <- gsDesign::gsDesign(
k = k, test.type = test.type, alpha = alpha / sided,
beta = beta, astar = astar, n.fix = x$d, timing = timing,
sfu = sfu, sfupar = sfupar, sfl = sfl, sflpar = sflpar, tol = tol,
delta1 = log(hr), delta0 = log(hr0),
usTime = usTime, lsTime = lsTime
)
y$hr <- hr
y$hr0 <- hr0
y$R <- x$R
y$S <- x$S
y$minfup <- x$minfup
# Inflate fixed design enrollment to get targeted events
inflate <- max(y$n.I) / x$d
y$gamma <- x$gamma * inflate
y$eDC <- inflate * eDC
y$eDE <- inflate * eDE
y$eNC <- inflate * eNC
y$eNE <- inflate * eNE
y$ratio <- ratio
y$lambdaC <- x$lambdaC
y$etaC <- x$etaC
y$etaE <- x$etaE
y$variable <- x$variable
y$tol <- tol
y$T <- calendarTime
class(y) <- c("gsSurv", "gsDesign")

nameR <- nameperiod(cumsum(y$R))
stratnames <- paste("Stratum", 1:ncol(y$lambdaC))
if (is.null(y$S)) {
nameS <- "0-Inf"
} else {
nameS <- nameperiod(cumsum(c(y$S, Inf)))
}
rownames(y$lambdaC) <- nameS
colnames(y$lambdaC) <- stratnames
rownames(y$etaC) <- nameS
colnames(y$etaC) <- stratnames
rownames(y$etaE) <- nameS
colnames(y$etaE) <- stratnames
rownames(y$gamma) <- nameR
colnames(y$gamma) <- stratnames
return(y)
}
Loading

0 comments on commit 6de6583

Please sign in to comment.