Nothing
R/RcppExports.R
In trtswitch: Treatment Switching
Defines functions
match3
qrcpp
survsplit
hasVariable
findInterval3
tssim
tsesimpcpp
tsegestsim
tsegestcpp
residuals_phregcpp
survfit_phregcpp
phregcpp
residuals_liferegcpp
liferegcpp
rmdiff
rmest
lrtest
kmdiff
kmest
survQuantile
nscpp
bscpp
splineDesigncpp
rpsftmcpp
logisregcpp
ipecpp
ipcwcpp
Documented in
bscpp
findInterval3
kmdiff
kmest
lrtest
nscpp
qrcpp
rmdiff
rmest
splineDesigncpp
survQuantile
survsplit
tsegestsim
tssim
# Generated by using Rcpp::compileAttributes() -> do not edit by hand
# Generator token: 10BE3573-1514-4C36-9D1C-5A225CD40393
ipcwcpp <- function(data, id = "id", stratum = "", tstart = "tstart", tstop = "tstop", event = "event", treat = "treat", swtrt = "swtrt", swtrt_time = "swtrt_time", base_cov = "", numerator = "", denominator = "", logistic_switching_model = 0L, strata_main_effect_only = 1L, firth = 0L, flic = 0L, ns_df = 3L, stabilized_weights = 1L, trunc = 0, trunc_upper_only = 1L, swtrt_control_only = 1L, alpha = 0.05, ties = "efron", boot = 1L, n_boot = 1000L, seed = NA_integer_) {
.Call(`_trtswitch_ipcwcpp`, data, id, stratum, tstart, tstop, event, treat, swtrt, swtrt_time, base_cov, numerator, denominator, logistic_switching_model, strata_main_effect_only, firth, flic, ns_df, stabilized_weights, trunc, trunc_upper_only, swtrt_control_only, alpha, ties, boot, n_boot, seed)
}
ipecpp <- function(data, id = "id", stratum = "", time = "time", event = "event", treat = "treat", rx = "rx", censor_time = "censor_time", base_cov = "", aft_dist = "weibull", strata_main_effect_only = 1L, low_psi = -2, hi_psi = 2, treat_modifier = 1, recensor = 1L, admin_recensor_only = 1L, autoswitch = 1L, alpha = 0.05, ties = "efron", tol = 1.0e-6, boot = 0L, n_boot = 1000L, seed = NA_integer_) {
.Call(`_trtswitch_ipecpp`, data, id, stratum, time, event, treat, rx, censor_time, base_cov, aft_dist, strata_main_effect_only, low_psi, hi_psi, treat_modifier, recensor, admin_recensor_only, autoswitch, alpha, ties, tol, boot, n_boot, seed)
}
logisregcpp <- function(data, rep = "", event = "event", covariates = "", freq = "", weight = "", offset = "", id = "", link = "logit", robust = 0L, firth = 0L, flic = 0L, plci = 0L, alpha = 0.05, maxiter = 50L, eps = 1.0e-9) {
.Call(`_trtswitch_logisregcpp`, data, rep, event, covariates, freq, weight, offset, id, link, robust, firth, flic, plci, alpha, maxiter, eps)
}
rpsftmcpp <- function(data, id = "id", stratum = "", time = "time", event = "event", treat = "treat", rx = "rx", censor_time = "censor_time", base_cov = "", psi_test = "logrank", aft_dist = "weibull", strata_main_effect_only = 1L, low_psi = -2, hi_psi = 2, n_eval_z = 101L, treat_modifier = 1, recensor = 1L, admin_recensor_only = 1L, autoswitch = 1L, gridsearch = 0L, alpha = 0.05, ties = "efron", tol = 1.0e-6, boot = 0L, n_boot = 1000L, seed = NA_integer_) {
.Call(`_trtswitch_rpsftmcpp`, data, id, stratum, time, event, treat, rx, censor_time, base_cov, psi_test, aft_dist, strata_main_effect_only, low_psi, hi_psi, n_eval_z, treat_modifier, recensor, admin_recensor_only, autoswitch, gridsearch, alpha, ties, tol, boot, n_boot, seed)
}
#' @title B-Spline Design Matrix
#' @description Computes the design matrix for B-splines based on the
#' specified \code{knots} and evaluated at the values in \code{x}.
#'
#' @param knots A numeric vector specifying the positions of the knots,
#' including both boundary and internal knots.
#' @param x A numeric vector of values where the B-spline functions
#' or their derivatives will be evaluated. The values of \code{x}
#' must lie within the range of the "inner" knots, i.e., between
#' \code{knots[ord]} and \code{knots[length(knots) - (ord - 1)]}.
#' @param ord A positive integer indicating the order of the B-spline.
#' This corresponds to the number of coefficients in each piecewise
#' polynomial segment, where \code{ord = degree + 1}.
#' @param derivs An integer vector specifying the order of derivatives
#' to be evaluated at the corresponding \code{x} values. Each value
#' must be between \code{0} and \code{ord - 1}, and the vector is
#' conceptually recycled to match the length of \code{x}.
#' The default is \code{0}, meaning the B-spline functions themselves
#' are evaluated.
#'
#' @return A matrix with dimensions
#' \code{c(length(x), length(knots) - ord)}. Each row corresponds
#' to a value in \code{x} and contains the coefficients of the
#' B-splines, or the specified derivatives, as defined by the
#' \code{knots} and evaluated at that particular value of \code{x}.
#' The total number of B-splines is \code{length(knots) - ord},
#' with each B-spline defined by a set of \code{ord} consecutive knots.
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @examples
#'
#' splineDesigncpp(knots = 1:10, x = 4:7)
#' splineDesigncpp(knots = 1:10, x = 4:7, derivs = 1)
#'
#' @export
splineDesigncpp <- function(knots = NA_real_, x = NA_real_, ord = 4L, derivs = as.integer( c(0))) {
.Call(`_trtswitch_splineDesigncpp`, knots, x, ord, derivs)
}
#' @title B-Spline Basis for Polynomial Splines
#' @description Computes the B-spline basis matrix for a given polynomial
#' spline.
#'
#' @param x A numeric vector representing the predictor variable.
#' @param df Degrees of freedom, specifying the number of columns in the
#' basis matrix. If \code{df} is provided, the function automatically
#' selects \code{df - degree - intercept} internal knots based on
#' appropriate quantiles of \code{x}, ignoring any missing values.
#' @param knots A numeric vector specifying the internal breakpoints
#' that define the spline. If not provided, \code{df} must be specified.
#' @param degree An integer specifying the degree of the piecewise
#' polynomial. The default value is \code{3}, which corresponds to
#' cubic splines.
#' @param intercept A logical value indicating whether to include an
#' intercept in the basis. The default is \code{FALSE}.
#' @param boundary_knots A numeric vector of length 2 specifying the
#' boundary points where the B-spline basis should be anchored.
#' If not supplied, the default is the range of non-missing values
#' in \code{x}.
#' @param warn_outside A logical value indicating whether a warning
#' should be issued if any values of \code{x} fall outside the
#' specified boundary knots.
#'
#' @return A matrix with dimensions \code{c(length(x), df)}. If
#' \code{df} is provided, the matrix will have \code{df} columns.
#' Alternatively, if \code{knots} are supplied, the number of columns
#' will be \code{length(knots) + degree + intercept}. The matrix
#' contains attributes that correspond to the arguments
#' passed to the \code{bscpp} function.
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @examples
#'
#' bscpp(women$height, df = 5)
#'
#' @export
bscpp <- function(x = NA_real_, df = NA_integer_, knots = NA_real_, degree = 3L, intercept = 0L, boundary_knots = NA_real_, warn_outside = 1L) {
.Call(`_trtswitch_bscpp`, x, df, knots, degree, intercept, boundary_knots, warn_outside)
}
#' @title Natural Cubic Spline Basis
#' @description Computes the B-spline basis matrix for a natural cubic
#' spline.
#'
#' @param x A numeric vector representing the predictor variable.
#' Missing values are allowed.
#' @param df Degrees of freedom, specifying the number of columns in
#' the basis matrix. If \code{df} is provided, the function selects
#' \code{df - 1 - intercept} internal knots based on appropriate
#' quantiles of \code{x}, ignoring any missing values.
#' @param knots A numeric vector specifying the internal breakpoints
#' that define the spline. If provided, the number of degrees of
#' freedom will be determined by the length of \code{knots}.
#' @param intercept A logical value indicating whether to include an
#' intercept in the basis. The default is \code{FALSE}.
#' @param boundary_knots A numeric vector of length 2 specifying the
#' boundary points where the natural boundary conditions are applied
#' and the B-spline basis is anchored. If not supplied, the default
#' is the range of non-missing values in \code{x}.
#'
#' @return A matrix with dimensions \code{c(length(x), df)}, where
#' \code{df} is either provided directly or computed as
#' \code{length(knots) + 1 + intercept} when \code{knots} are supplied.
#' The matrix contains attributes that correspond to the arguments
#' passed to the \code{nscpp} function.
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @examples
#'
#' nscpp(women$height, df = 5)
#'
#' @export
nscpp <- function(x = NA_real_, df = NA_integer_, knots = NA_real_, intercept = 0L, boundary_knots = NA_real_) {
.Call(`_trtswitch_nscpp`, x, df, knots, intercept, boundary_knots)
}
#' @title Brookmeyer-Crowley Confidence Interval for Quantiles of
#' Right-Censored Time-to-Event Data
#' @description Obtains the Brookmeyer-Crowley confidence
#' interval for quantiles of right-censored time-to-event data.
#'
#' @param time The vector of possibly right-censored survival times.
#' @param event The vector of event indicators.
#' @param cilevel The confidence interval level. Defaults to 0.95.
#' @param transform The transformation of the survival function to use
#' to construct the confidence interval. Options include
#' "linear" (alternatively "plain"), "log",
#' "loglog" (alternatively "log-log" or "cloglog"),
#' "asinsqrt" (alternatively "asin" or "arcsin"), and "logit".
#' Defaults to "loglog".
#'
#' @param probs The vector of probabilities to calculate the quantiles.
#' Defaults to c(0.25, 0.5, 0.75).
#'
#' @return A data frame containing the estimated quantile and
#' confidence interval corresponding to each specified probability.
#' It includes the following variables:
#'
#' * \code{prob}: The probability to calculate the quantile.
#'
#' * \code{quantile}: The estimated quantile.
#'
#' * \code{lower}: The lower limit of the confidence interval.
#'
#' * \code{upper}: The upper limit of the confidence interval.
#'
#' * \code{cilevel}: The confidence interval level.
#'
#' * \code{transform}: The transformation of the survival function to use
#' to construct the confidence interval.
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @examples
#'
#' survQuantile(
#' time = c(33.7, 3.9, 10.5, 5.4, 19.5, 23.8, 7.9, 16.9, 16.6,
#' 33.7, 17.1, 7.9, 10.5, 38),
#' event = c(0, 1, 1, 1, 1, 0, 1, 0, 0, 0, 0, 0, 1, 1),
#' probs = c(0.25, 0.5, 0.75))
#'
#' @export
survQuantile <- function(time = NA_real_, event = NA_real_, cilevel = 0.95, transform = "loglog", probs = NA_real_) {
.Call(`_trtswitch_survQuantile`, time, event, cilevel, transform, probs)
}
#' @title Kaplan-Meier Estimates of Survival Curve
#' @description Obtains the Kaplan-Meier estimates of the survival curve.
#'
#' @param data The input data frame that contains the following variables:
#'
#' * \code{rep}: The replication for by-group processing.
#'
#' * \code{stratum}: The stratum.
#'
#' * \code{time}: The possibly right-censored survival time.
#'
#' * \code{event}: The event indicator.
#'
#' @param rep The name(s) of the replication variable(s) in the input data.
#' @param stratum The name(s) of the stratum variable(s) in the input data.
#' @param time The name of the time variable in the input data.
#' @param event The name of the event variable in the input data.
#' @param conftype The type of the confidence interval. One of "none",
#' "plain", "log", "log-log" (the default), or "arcsin".
#' The arcsin option bases the intervals on asin(sqrt(survival)).
#' @param conflev The level of the two-sided confidence interval for
#' the survival probabilities. Defaults to 0.95.
#' @param keep_censor Whether to retain the censoring time in the output
#' data frame.
#'
#' @return A data frame with the following variables:
#'
#' * \code{size}: The number of subjects in the stratum.
#'
#' * \code{time}: The event time.
#'
#' * \code{nrisk}: The number of subjects at risk.
#'
#' * \code{nevent}: The number of subjects having the event.
#'
#' * \code{ncensor}: The number of censored subjects.
#'
#' * \code{survival}: The Kaplan-Meier estimate of the survival probability.
#'
#' * \code{stderr}: The standard error of the estimated survival
#' probability based on the Greendwood formula.
#'
#' * \code{lower}: The lower bound of confidence interval if requested.
#'
#' * \code{upper}: The upper bound of confidence interval if requested.
#'
#' * \code{conflev}: The level of confidence interval if requested.
#'
#' * \code{conftype}: The type of confidence interval if requested.
#'
#' * \code{stratum}: The stratum.
#'
#' * \code{rep}: The replication.
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @examples
#'
#' kmest(data = aml, stratum = "x", time = "time", event = "status")
#'
#' @export
kmest <- function(data, rep = "", stratum = "", time = "time", event = "event", conftype = "log-log", conflev = 0.95, keep_censor = 0L) {
.Call(`_trtswitch_kmest`, data, rep, stratum, time, event, conftype, conflev, keep_censor)
}
#' @title Estimate of Milestone Survival Difference
#' @description Obtains the estimate of milestone survival difference
#' between two treatment groups.
#'
#' @param data The input data frame that contains the following variables:
#'
#' * \code{rep}: The replication for by-group processing.
#'
#' * \code{stratum}: The stratum.
#'
#' * \code{treat}: The treatment.
#'
#' * \code{time}: The possibly right-censored survival time.
#'
#' * \code{event}: The event indicator.
#'
#' @param rep The name of the replication variable in the input data.
#' @param stratum The name of the stratum variable in the input data.
#' @param treat The name of the treatment variable in the input data.
#' @param time The name of the time variable in the input data.
#' @param event The name of the event variable in the input data.
#' @param milestone The milestone time at which to calculate the
#' survival probability.
#' @param survDiffH0 The difference in milestone survival probabilities
#' under the null hypothesis. Defaults to 0 for superiority test.
#' @param conflev The level of the two-sided confidence interval for
#' the difference in milestone survival probabilities. Defaults to 0.95.
#'
#' @return A data frame with the following variables:
#'
#' * \code{rep}: The replication.
#'
#' * \code{milestone}: The milestone time relative to randomization.
#'
#' * \code{survDiffH0}: The difference in milestone survival probabilities
#' under the null hypothesis.
#'
#' * \code{surv1}: The estimated milestone survival probability for
#' the treatment group.
#'
#' * \code{surv2}: The estimated milestone survival probability for
#' the control group.
#'
#' * \code{survDiff}: The estimated difference in milestone survival
#' probabilities.
#'
#' * \code{vsurv1}: The variance for surv1.
#'
#' * \code{vsurv2}: The variance for surv2.
#'
#' * \code{vsurvDiff}: The variance for survDiff.
#'
#' * \code{survDiffZ}: The Z-statistic value.
#'
#' * \code{survDiffPValue}: The one-sided p-value.
#'
#' * \code{lower}: The lower bound of confidence interval.
#'
#' * \code{upper}: The upper bound of confidence interval.
#'
#' * \code{conflev}: The level of confidence interval.
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @examples
#'
#' df <- kmdiff(data = rawdata, rep = "iterationNumber",
#' stratum = "stratum", treat = "treatmentGroup",
#' time = "timeUnderObservation", event = "event",
#' milestone = 12)
#' head(df)
#'
#' @export
kmdiff <- function(data, rep = "", stratum = "", treat = "treat", time = "time", event = "event", milestone = NA_real_, survDiffH0 = 0, conflev = 0.95) {
.Call(`_trtswitch_kmdiff`, data, rep, stratum, treat, time, event, milestone, survDiffH0, conflev)
}
#' @title Log-Rank Test of Survival Curve Difference
#' @description Obtains the log-rank test using the Fleming-Harrington
#' family of weights.
#'
#' @param data The input data frame that contains the following variables:
#'
#' * \code{rep}: The replication for by-group processing.
#'
#' * \code{stratum}: The stratum.
#'
#' * \code{treat}: The treatment.
#'
#' * \code{time}: The possibly right-censored survival time.
#'
#' * \code{event}: The event indicator.
#'
#' @param rep The name of the replication variable in the input data.
#' @param stratum The name of the stratum variable in the input data.
#' @param treat The name of the treatment variable in the input data.
#' @param time The name of the time variable in the input data.
#' @param event The name of the event variable in the input data.
#' @param rho1 The first parameter of the Fleming-Harrington family of
#' weighted log-rank test. Defaults to 0 for conventional log-rank test.
#' @param rho2 The second parameter of the Fleming-Harrington family of
#' weighted log-rank test. Defaults to 0 for conventional log-rank test.
#'
#' @return A data frame with the following variables:
#'
#' * \code{uscore}: The numerator of the log-rank test statistic.
#'
#' * \code{vscore}: The variance of the log-rank score test statistic.
#'
#' * \code{logRankZ}: The Z-statistic value.
#'
#' * \code{logRankPValue}: The one-sided p-value.
#'
#' * \code{rho1}: The first parameter of the Fleming-Harrington weights.
#'
#' * \code{rho2}: The second parameter of the Fleming-Harrington weights.
#'
#' * \code{rep}: The replication.
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @examples
#'
#' df <- lrtest(data = rawdata, rep = "iterationNumber",
#' stratum = "stratum", treat = "treatmentGroup",
#' time = "timeUnderObservation", event = "event",
#' rho1 = 0.5, rho2 = 0)
#' head(df)
#'
#' @export
lrtest <- function(data, rep = "", stratum = "", treat = "treat", time = "time", event = "event", rho1 = 0, rho2 = 0) {
.Call(`_trtswitch_lrtest`, data, rep, stratum, treat, time, event, rho1, rho2)
}
#' @title Estimate of Restricted Mean Survival Time
#' @description Obtains the estimate of restricted means survival time
#' for each stratum.
#'
#' @param data The input data frame that contains the following variables:
#'
#' * \code{rep}: The replication for by-group processing.
#'
#' * \code{stratum}: The stratum.
#'
#' * \code{time}: The possibly right-censored survival time.
#'
#' * \code{event}: The event indicator.
#'
#' @param rep The name of the replication variable in the input data.
#' @param stratum The name of the stratum variable in the input data.
#' @param time The name of the time variable in the input data.
#' @param event The name of the event variable in the input data.
#' @param milestone The milestone time at which to calculate the
#' restricted mean survival time.
#' @param conflev The level of the two-sided confidence interval for
#' the survival probabilities. Defaults to 0.95.
#' @param biascorrection Whether to apply bias correction for the
#' variance estimate. Defaults to no bias correction.
#'
#' @return A data frame with the following variables:
#'
#' * \code{rep}: The replication.
#'
#' * \code{stratum}: The stratum variable.
#'
#' * \code{size}: The number of subjects in the stratum.
#'
#' * \code{milestone}: The milestone time relative to randomization.
#'
#' * \code{rmst}: The estimate of restricted mean survival time.
#'
#' * \code{stderr}: The standard error of the estimated rmst.
#'
#' * \code{lower}: The lower bound of confidence interval if requested.
#'
#' * \code{upper}: The upper bound of confidence interval if requested.
#'
#' * \code{conflev}: The level of confidence interval if requested.
#'
#' * \code{biascorrection}: Whether to apply bias correction for the
#' variance estimate.
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @examples
#'
#' rmest(data = aml, stratum = "x",
#' time = "time", event = "status", milestone = 24)
#'
#' @export
rmest <- function(data, rep = "", stratum = "", time = "time", event = "event", milestone = NA_real_, conflev = 0.95, biascorrection = 0L) {
.Call(`_trtswitch_rmest`, data, rep, stratum, time, event, milestone, conflev, biascorrection)
}
#' @title Estimate of Restricted Mean Survival Time Difference
#' @description Obtains the estimate of restricted mean survival time
#' difference between two treatment groups.
#'
#' @param data The input data frame that contains the following variables:
#'
#' * \code{rep}: The replication for by-group processing.
#'
#' * \code{stratum}: The stratum.
#'
#' * \code{treat}: The treatment.
#'
#' * \code{time}: The possibly right-censored survival time.
#'
#' * \code{event}: The event indicator.
#'
#' @param rep The name of the replication variable in the input data.
#' @param stratum The name of the stratum variable in the input data.
#' @param treat The name of the treatment variable in the input data.
#' @param time The name of the time variable in the input data.
#' @param event The name of the event variable in the input data.
#' @param milestone The milestone time at which to calculate the
#' restricted mean survival time.
#' @param rmstDiffH0 The difference in restricted mean survival times
#' under the null hypothesis. Defaults to 0 for superiority test.
#' @param conflev The level of the two-sided confidence interval for
#' the difference in restricted mean survival times. Defaults to 0.95.
#' @param biascorrection Whether to apply bias correction for the
#' variance estimate of individual restricted mean survival times.
#' Defaults to no bias correction.
#'
#' @return A data frame with the following variables:
#'
#' * \code{rep}: The replication number.
#'
#' * \code{milestone}: The milestone time relative to randomization.
#'
#' * \code{rmstDiffH0}: The difference in restricted mean survival times
#' under the null hypothesis.
#'
#' * \code{rmst1}: The estimated restricted mean survival time for
#' the treatment group.
#'
#' * \code{rmst2}: The estimated restricted mean survival time for
#' the control group.
#'
#' * \code{rmstDiff}: The estimated difference in restricted mean
#' survival times.
#'
#' * \code{vrmst1}: The variance for rmst1.
#'
#' * \code{vrmst2}: The variance for rmst2.
#'
#' * \code{vrmstDiff}: The variance for rmstDiff.
#'
#' * \code{rmstDiffZ}: The Z-statistic value.
#'
#' * \code{rmstDiffPValue}: The one-sided p-value.
#'
#' * \code{lower}: The lower bound of confidence interval.
#'
#' * \code{upper}: The upper bound of confidence interval.
#'
#' * \code{conflev}: The level of confidence interval.
#'
#' * \code{biascorrection}: Whether to apply bias correction for the
#' variance estimate of individual restricted mean survival times.
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @examples
#'
#' df <- rmdiff(data = rawdata, rep = "iterationNumber",
#' stratum = "stratum", treat = "treatmentGroup",
#' time = "timeUnderObservation", event = "event",
#' milestone = 12)
#' head(df)
#'
#' @export
rmdiff <- function(data, rep = "", stratum = "", treat = "treat", time = "time", event = "event", milestone = NA_real_, rmstDiffH0 = 0, conflev = 0.95, biascorrection = 0L) {
.Call(`_trtswitch_rmdiff`, data, rep, stratum, treat, time, event, milestone, rmstDiffH0, conflev, biascorrection)
}
liferegcpp <- function(data, rep = "", stratum = "", time = "time", time2 = "", event = "event", covariates = "", weight = "", offset = "", id = "", dist = "weibull", robust = 0L, plci = 0L, alpha = 0.05, maxiter = 50L, eps = 1.0e-9) {
.Call(`_trtswitch_liferegcpp`, data, rep, stratum, time, time2, event, covariates, weight, offset, id, dist, robust, plci, alpha, maxiter, eps)
}
residuals_liferegcpp <- function(beta, vbeta, data, stratum = "", time = "time", time2 = "", event = "event", covariates = "", weight = "", offset = "", id = "", dist = "weibull", type = "response", collapse = 0L, weighted = 0L) {
.Call(`_trtswitch_residuals_liferegcpp`, beta, vbeta, data, stratum, time, time2, event, covariates, weight, offset, id, dist, type, collapse, weighted)
}
phregcpp <- function(data, rep = "", stratum = "", time = "time", time2 = "", event = "event", covariates = "", weight = "", offset = "", id = "", ties = "efron", robust = 0L, est_basehaz = 1L, est_resid = 1L, firth = 0L, plci = 0L, alpha = 0.05, maxiter = 50L, eps = 1.0e-9) {
.Call(`_trtswitch_phregcpp`, data, rep, stratum, time, time2, event, covariates, weight, offset, id, ties, robust, est_basehaz, est_resid, firth, plci, alpha, maxiter, eps)
}
survfit_phregcpp <- function(p, beta, vbeta, basehaz, newdata, covariates = "", stratum = "", offset = "", id = "", tstart = "", tstop = "", sefit = 1L, conftype = "log-log", conflev = 0.95) {
.Call(`_trtswitch_survfit_phregcpp`, p, beta, vbeta, basehaz, newdata, covariates, stratum, offset, id, tstart, tstop, sefit, conftype, conflev)
}
residuals_phregcpp <- function(p, beta, vbeta, resmart, data, stratum = "", time = "time", time2 = "", event = "event", covariates = "", weight = "", offset = "", id = "", ties = "efron", type = "schoenfeld", collapse = 0L, weighted = 0L) {
.Call(`_trtswitch_residuals_phregcpp`, p, beta, vbeta, resmart, data, stratum, time, time2, event, covariates, weight, offset, id, ties, type, collapse, weighted)
}
tsegestcpp <- function(data, id = "id", stratum = "", tstart = "tstart", tstop = "tstop", event = "event", treat = "treat", censor_time = "censor_time", pd = "pd", pd_time = "pd_time", swtrt = "swtrt", swtrt_time = "swtrt_time", base_cov = "", conf_cov = "", low_psi = -2, hi_psi = 2, n_eval_z = 101L, strata_main_effect_only = 1L, firth = 0L, flic = 0L, recensor = 1L, admin_recensor_only = 1L, swtrt_control_only = 1L, gridsearch = 0L, alpha = 0.05, ties = "efron", tol = 1.0e-6, offset = 1, boot = 1L, n_boot = 1000L, seed = NA_integer_) {
.Call(`_trtswitch_tsegestcpp`, data, id, stratum, tstart, tstop, event, treat, censor_time, pd, pd_time, swtrt, swtrt_time, base_cov, conf_cov, low_psi, hi_psi, n_eval_z, strata_main_effect_only, firth, flic, recensor, admin_recensor_only, swtrt_control_only, gridsearch, alpha, ties, tol, offset, boot, n_boot, seed)
}
#' @title Simulate Survival Data for Two-Stage Estimation Method Using
#' g-estimation
#' @description Obtains the simulated data for baseline prognosis,
#' disease progression, treatment switching, death, and
#' time-dependent covariates.
#'
#' @param n The total sample size for two treatment arms combined.
#' @param allocation1 The number of subjects in the active treatment group
#' in a randomization block.
#' @param allocation2 The number of subjects in the control group in
#' a randomization block.
#' @param pbprog The probability of having poor prognosis at baseline.
#' @param trtlghr The treatment effect in terms of log hazard ratio.
#' @param bprogsl The poor prognosis effect in terms of log hazard ratio.
#' @param shape1 The shape parameter for the Weibull event distribution
#' for the first component.
#' @param scale1 The scale parameter for the Weibull event distribution
#' for the first component.
#' @param shape2 The shape parameter for the Weibull event distribution
#' for the second component.
#' @param scale2 The scale parameter for the Weibull event distribution
#' for the second component.
#' @param pmix The mixing probability of the first component Weibull
#' distribution.
#' @param admin The administrative censoring time.
#' @param pcatnotrtbprog The probability of developing metastatic disease
#' on control treatment with poor baseline prognosis.
#' @param pcattrtbprog The probability of developing metastatic disease
#' on active treatment with poor baseline prognosis.
#' @param pcatnotrt The probability of developing metastatic disease
#' on control treatment with good baseline prognosis.
#' @param pcattrt The probability of developing metastatic disease
#' on active treatment with good baseline prognosis.
#' @param catmult The impact of metastatic disease on shortening remaining
#' survival time.
#' @param tdxo Whether treatment crossover depends on time-dependent
#' covariates between disease progression and treatment switching.
#' @param ppoor The probability of switching for poor baseline prognosis
#' with no metastatic disease.
#' @param pgood The probability of switching for good baseline prognosis
#' with no metastatic disease.
#' @param ppoormet The probability of switching for poor baseline prognosis
#' after developing metastatic disease.
#' @param pgoodmet The probability of switching for good baseline prognosis
#' after developing metastatic disease.
#' @param xomult The direct effect of crossover on extending remaining
#' survival time.
#' @param milestone The milestone to calculate restricted mean survival
#' time.
#' @param outputRawDataset Whether to output the raw data set.
#' @param seed The seed to reproduce the simulation results.
#' The seed from the environment will be used if left unspecified.
#'
#' @return A list with two data frames.
#'
#' * \code{sumdata}: A summary data frame with the following variables:
#'
#' - \code{simtrueconstmean}: The true control group restricted mean
#' survival time (RMST).
#'
#' - \code{simtrueconstlb}: The lower bound for control group RMST.
#'
#' - \code{simtrueconstub}: The upper bound for control group RMST.
#'
#' - \code{simtrueconstse}: The standard error for control group RMST.
#'
#' - \code{simtrueexpstmean}: The true experimental group restricted
#' mean survival time (RMST).
#'
#' - \code{simtrueexpstlb}: The lower bound for experimental group RMST.
#'
#' - \code{simtrueexpstub}: The upper bound for experimental group RMST.
#'
#' - \code{simtrueexpstse}: The standard error for experimental group
#' RMST.
#'
#' - \code{simtrue_coxwbprog_hr}: The treatment hazard ratio from the
#' Cox model adjusting for baseline prognosis.
#'
#' - \code{simtrue_cox_hr}: The treatment hazard ratio from the Cox
#' model without adjusting for baseline prognosis.
#'
#' * \code{paneldata}: A counting process style subject-level data frame
#' with the following variables:
#'
#' - \code{id}: The subject ID.
#'
#' - \code{trtrand}: The randomized treatment arm.
#'
#' - \code{bprog}: Whether the patient had poor baseline prognosis.
#'
#' - \code{tstart}: The left end of time interval.
#'
#' - \code{tstop}: The right end of time interval.
#'
#' - \code{event}: Whether the patient died at the end of the interval.
#'
#' - \code{timeOS}: The observed survival time.
#'
#' - \code{died}: Whether the patient died during the study.
#'
#' - \code{progressed}: Whether the patient had disease progression.
#'
#' - \code{timePFSobs}: The observed time of disease progression at
#' regular scheduled visits.
#'
#' - \code{progtdc}: The time-dependent covariate for progression.
#'
#' - \code{catevent}: Whether the patient developed metastatic disease.
#'
#' - \code{cattime}: When the patient developed metastatic disease.
#'
#' - \code{cattdc}: The time-dependent covariate for cat event.
#'
#' - \code{catlag}: The lagged value of \code{cattdc}.
#'
#' - \code{xo}: Whether the patient switched treatment.
#'
#' - \code{xotime}: When the patient switched treatment.
#'
#' - \code{xotdc}: The time-dependent covariate for treatment
#' switching.
#'
#' - \code{xotime_upper}: The upper bound of treatment switching time.
#'
#' - \code{censor_time}: The administrative censoring time.
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @references
#' NR Latimer, IR White, K Tilling, and U Siebert.
#' Improved two-stage estimation to adjust for treatment switching in
#' randomised trials: g-estimation to address time-dependent confounding.
#' Statistical Methods in Medical Research. 2020;29(10):2900-2918.
#'
#' @examples
#'
#' sim1 <- tsegestsim(
#' n = 500, allocation1 = 2, allocation2 = 1, pbprog = 0.5,
#' trtlghr = -0.5, bprogsl = 0.3, shape1 = 1.8,
#' scale1 = 360, shape2 = 1.7, scale2 = 688,
#' pmix = 0.5, admin = 5000, pcatnotrtbprog = 0.5,
#' pcattrtbprog = 0.25, pcatnotrt = 0.2, pcattrt = 0.1,
#' catmult = 0.5, tdxo = 1, ppoor = 0.1, pgood = 0.04,
#' ppoormet = 0.4, pgoodmet = 0.2, xomult = 1.4188308,
#' milestone = 546, outputRawDataset = 1, seed = 2000)
#'
#' @export
tsegestsim <- function(n = 500L, allocation1 = 2L, allocation2 = 1L, pbprog = 0.5, trtlghr = -0.5, bprogsl = 0.3, shape1 = 1.8, scale1 = 360, shape2 = 1.7, scale2 = 688, pmix = 0.5, admin = 5000, pcatnotrtbprog = 0.5, pcattrtbprog = 0.25, pcatnotrt = 0.2, pcattrt = 0.1, catmult = 0.5, tdxo = 1, ppoor = 0.1, pgood = 0.04, ppoormet = 0.4, pgoodmet = 0.2, xomult = 1.4188308, milestone = 546, outputRawDataset = 1L, seed = NA_integer_) {
.Call(`_trtswitch_tsegestsim`, n, allocation1, allocation2, pbprog, trtlghr, bprogsl, shape1, scale1, shape2, scale2, pmix, admin, pcatnotrtbprog, pcattrtbprog, pcatnotrt, pcattrt, catmult, tdxo, ppoor, pgood, ppoormet, pgoodmet, xomult, milestone, outputRawDataset, seed)
}
tsesimpcpp <- function(data, id = "id", stratum = "", time = "time", event = "event", treat = "treat", censor_time = "censor_time", pd = "pd", pd_time = "pd_time", swtrt = "swtrt", swtrt_time = "swtrt_time", base_cov = "", base2_cov = "", aft_dist = "weibull", strata_main_effect_only = 1L, recensor = 1L, admin_recensor_only = 1L, swtrt_control_only = 1L, alpha = 0.05, ties = "efron", offset = 1, boot = 1L, n_boot = 1000L, seed = NA_integer_) {
.Call(`_trtswitch_tsesimpcpp`, data, id, stratum, time, event, treat, censor_time, pd, pd_time, swtrt, swtrt_time, base_cov, base2_cov, aft_dist, strata_main_effect_only, recensor, admin_recensor_only, swtrt_control_only, alpha, ties, offset, boot, n_boot, seed)
}
#' @title Simulate Data for Treatment Switching
#' @description Simulates data for studies involving treatment switching,
#' incorporating time-dependent confounding. The generated data can be used
#' to evaluate methods for handling treatment switching in survival analysis.
#'
#' @param tdxo Logical indicator for timing of treatment switching:
#' \itemize{
#' \item 1: Treatment switching can occur at or after disease progression.
#' \item 0: Treatment switching is restricted to the time of disease
#' progression.
#' }
#' @param coxo Logical indicator for arm-specific treatment switching:
#' \itemize{
#' \item 1: Treatment switching occurs only in the control arm.
#' \item 0: Treatment switching is allowed in both arms.
#' }
#' @param p_R Probability of randomization to the experimental arm.
#' @param p_X_1 Probability of poor baseline prognosis in the experimental
#' arm.
#' @param p_X_0 Probability of poor baseline prognosis in the control arm.
#' @param rate_T Baseline hazard rate for time to death.
#' @param beta1 Log hazard ratio for randomized treatment (\code{R}).
#' @param beta2 Log hazard ratio for baseline covariate (\code{X}).
#' @param gamma0 Intercept for the time-dependent covariate model (\code{L}).
#' @param gamma1 Coefficient for lagged treatment switching (\code{Alag})
#' in the \code{L} model.
#' @param gamma2 Coefficient for lagged \code{L} in the \code{L} model.
#' @param gamma3 Coefficient for baseline covariate (\code{X}) in the
#' \code{L} model.
#' @param gamma4 Coefficient for randomized treatment (\code{R}) in the
#' \code{L} model.
#' @param zeta0 Intercept for the disease progression model (\code{Z}).
#' @param zeta1 Coefficient for \code{L} in the \code{Z} model.
#' @param zeta2 Coefficient for baseline covariate (\code{X}) in the
#' \code{Z} model.
#' @param zeta3 Coefficient for randomized treatment (\code{R}) in the
#' \code{Z} model.
#' @param alpha0 Intercept for the treatment switching model (\code{A}).
#' @param alpha1 Coefficient for \code{L} in the \code{A} model.
#' @param alpha2 Coefficient for baseline covariate (\code{X}) in the
#' \code{A} model.
#' @param theta1_1 Negative log time ratio for \code{A} (experimental arm).
#' @param theta1_0 Negative log time ratio for \code{A} (control arm).
#' @param theta2 Negative log time ratio for \code{L}.
#' @param rate_C Hazard rate for random (dropout) censoring.
#' @param followup Number of treatment cycles per subject.
#' @param days Number of days in each treatment cycle.
#' @param n Number of subjects per simulation.
#' @param NSim Number of simulated datasets.
#' @param seed Random seed for reproducibility.
#'
#' @return
#' A list of data frames, each containing simulated longitudinal and
#' event history data with the following variables:
#'
#' * \code{id}: Subject identifier.
#'
#' * \code{trtrand}: Randomized treatment assignment (0 = control,
#' 1 = experimental)
#'
#' * \code{bprog}: Baseline prognosis (0 = good, 1 = poor).
#'
#' * \code{tpoint}: Treatment cycle index.
#'
#' * \code{tstart}: Start day of the treatment cycle.
#'
#' * \code{tstop}: End day of the treatment cycle.
#'
#' * \code{L}: Time-dependent covariate predicting survival and switching;
#' affected by treatment switching.
#'
#' * \code{Llag}: Lagged value of \code{L}.
#'
#' * \code{Z}: Disease progression status at \code{tstop}.
#'
#' * \code{A}: Treatment switching status at \code{tstop}.
#'
#' * \code{Alag}: Lagged value of \code{A}.
#'
#' * \code{Y}: Death indicator at \code{tstop}.
#'
#' * \code{timeOS}: Observed time to death or censoring.
#'
#' * \code{died}: Indicator of death by end of follow-up.
#'
#' * \code{progressed}: Indicator of disease progression by end of follow-up.
#'
#' * \code{timePD}: Observed time to progression or censoring.
#'
#' * \code{xo}: Indicator for whether treatment switching occurred.
#'
#' * \code{xotime}: Time of treatment switching (if applicable).
#'
#' * \code{censor_time}: Administrative censoring time.
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @references
#'
#' Jessica G. Young, and Eric J. Tchetgen Tchetgen.
#' Simulation from a known Cox MSM using standard parametric models
#' for the g-formula.
#' Statistics in Medicine. 2014;33(6):1001-1014.
#'
#' NR Latimer, IR White, K Tilling, and U Siebert.
#' Improved two-stage estimation to adjust for treatment switching in
#' randomised trials: g-estimation to address time-dependent confounding.
#' Statistical Methods in Medical Research. 2020;29(10):2900-2918.
#'
#' Jing Xu, Guohui Liu, and Bingxia Wang.
#' Bias and type I error control in correcting treatment effect
#' for treatment switching using marginal structural models
#' in Phse III oncology trials.
#' Journal of Biopharmaceutical Statistics. 2022;32(6):897-914.
#'
#' @examples
#'
#' simulated.data <- tssim(
#' tdxo = 0, coxo = 0, p_R = 0.5, p_X_1 = 0.3, p_X_0 = 0.3,
#' rate_T = 0.002, beta1 = -0.5, beta2 = 0.3,
#' gamma0 = 0.3, gamma1 = -0.9, gamma2 = 0.7, gamma3 = 1.1, gamma4 = -0.8,
#' zeta0 = -3.5, zeta1 = 0.5, zeta2 = 0.2, zeta3 = -0.4,
#' alpha0 = 0.5, alpha1 = 0.5, alpha2 = 0.4,
#' theta1_1 = -0.4, theta1_0 = -0.4, theta2 = 0.2,
#' rate_C = 0.0000855, followup = 20, days = 30,
#' n = 500, NSim = 100, seed = 314159)
#'
#' @export
tssim <- function(tdxo = 0L, coxo = 1L, p_R = 0.5, p_X_1 = NA_real_, p_X_0 = NA_real_, rate_T = NA_real_, beta1 = NA_real_, beta2 = NA_real_, gamma0 = NA_real_, gamma1 = NA_real_, gamma2 = NA_real_, gamma3 = NA_real_, gamma4 = NA_real_, zeta0 = NA_real_, zeta1 = NA_real_, zeta2 = NA_real_, zeta3 = NA_real_, alpha0 = NA_real_, alpha1 = NA_real_, alpha2 = NA_real_, theta1_1 = NA_real_, theta1_0 = NA_real_, theta2 = NA_real_, rate_C = NA_real_, followup = NA_integer_, days = NA_integer_, n = NA_integer_, NSim = 1000L, seed = NA_integer_) {
.Call(`_trtswitch_tssim`, tdxo, coxo, p_R, p_X_1, p_X_0, rate_T, beta1, beta2, gamma0, gamma1, gamma2, gamma3, gamma4, zeta0, zeta1, zeta2, zeta3, alpha0, alpha1, alpha2, theta1_1, theta1_0, theta2, rate_C, followup, days, n, NSim, seed)
}
#' @title Find Interval Numbers of Indices
#' @description The implementation of \code{findInterval()} in R from
#' Advanced R by Hadley Wickham. Given a vector of non-decreasing
#' breakpoints in v, find the interval containing each element of x; i.e.,
#' if \code{i <- findInterval3(x,v)}, for each index \code{j} in \code{x},
#' \code{v[i[j]] <= x[j] < v[i[j] + 1]}, where \code{v[0] := -Inf},
#' \code{v[N+1] := +Inf}, and \code{N = length(v)}.
#'
#' @param x The numeric vector of interest.
#' @param v The vector of break points.
#' @return A vector of \code{length(x)} with values in \code{0:N} where
#' \code{N = length(v)}.
#'
#' @keywords internal
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @examples
#' x <- 2:18
#' v <- c(5, 10, 15) # create two bins [5,10) and [10,15)
#' cbind(x, findInterval3(x, v))
#'
#' @export
findInterval3 <- function(x, v) {
.Call(`_trtswitch_findInterval3`, x, v)
}
hasVariable <- function(df, varName) {
.Call(`_trtswitch_hasVariable`, df, varName)
}
#' @title Split a survival data set at specified cut points
#' @description For a given survival dataset and specified cut times,
#' each record is split into multiple subrecords at each cut time.
#' The resulting dataset is in counting process format, with each
#' subrecord containing a start time, stop time, and event status.
#' This is adapted from the survsplit.c function from the survival package.
#'
#' @param tstart The starting time of the time interval for
#' counting-process data.
#' @param tstop The stopping time of the time interval for
#' counting-process data.
#' @param cut The vector of cut points.
#'
#' @return A data frame with the following variables:
#'
#' * \code{row}: The row number of the observation in the input data
#' (starting from 0).
#'
#' * \code{start}: The starting time of the resulting subrecord.
#'
#' * \code{end}: The ending time of the resulting subrecord.
#'
#' * \code{censor}: Whether the subrecord lies strictly within a record
#' in the input data (1 for all but the last interval and 0 for the
#' last interval with cutpoint set equal to tstop).
#'
#' * \code{interval}: The interval number derived from cut (starting
#' from 0 if the interval lies to the left of the first cutpoint).
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @keywords internal
#'
#' @examples
#'
#' survsplit(15, 60, c(10, 30, 40))
#'
#' @export
survsplit <- function(tstart, tstop, cut) {
.Call(`_trtswitch_survsplit`, tstart, tstop, cut)
}
#' @title QR Decomposition of a Matrix
#' @description Computes the QR decomposition of a matrix.
#'
#' @param X A numeric matrix whose QR decomposition is to be computed.
#' @param tol The tolerance for detecting linear dependencies in the
#' columns of \code{X}.
#'
#' @details
#' This function performs Householder QR with column pivoting:
#' Given an \eqn{m}-by-\eqn{n} matrix \eqn{A} with \eqn{m \geq n},
#' the following algorithm computes \eqn{r = \textrm{rank}(A)} and
#' the factorization \eqn{Q^T A P} equal to
#' \tabular{ccccc}{
#' | \tab \eqn{R_{11}} \tab \eqn{R_{12}} \tab | \tab \eqn{r} \cr
#' | \tab 0 \tab 0 \tab | \tab \eqn{m-r} \cr
#' \tab \eqn{r} \tab \eqn{n-r} \tab \tab
#' }
#' with \eqn{Q = H_1 \cdots H_r} and \eqn{P = P_1 \cdots P_r}.
#' The upper triangular part of \eqn{A}
#' is overwritten by the upper triangular part of \eqn{R} and
#' components \eqn{(j+1):m} of
#' the \eqn{j}th Householder vector are stored in \eqn{A((j+1):m, j)}.
#' The permutation \eqn{P} is encoded in an integer vector \code{pivot}.
#'
#' @return A list with the following components:
#'
#' * \code{qr}: A matrix with the same dimensions as \code{X}. The upper
#' triangle contains the \code{R} of the decomposition and the lower
#' triangle contains Householder vectors (stored in compact form).
#'
#' * \code{rank}: The rank of \code{X} as computed by the decomposition.
#'
#' * \code{pivot}: The column permutation for the pivoting strategy used
#' during the decomposition.
#'
#' * \code{Q}: The complete \eqn{m}-by-\eqn{m} orthogonal matrix \eqn{Q}.
#'
#' * \code{R}: The complete \eqn{m}-by-\eqn{n} upper triangular
#' matrix \eqn{R}.
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @references
#' Gene N. Golub and Charles F. Van Loan.
#' Matrix Computations, second edition. Baltimore, Maryland:
#' The John Hopkins University Press, 1989, p.235.
#'
#' @examples
#'
#' hilbert <- function(n) { i <- 1:n; 1 / outer(i - 1, i, `+`) }
#' h9 <- hilbert(9)
#' qrcpp(h9)
#'
#' @export
qrcpp <- function(X, tol = 1e-12) {
.Call(`_trtswitch_qrcpp`, X, tol)
}
match3 <- function(id1, v1, id2, v2) {
.Call(`_trtswitch_match3`, id1, v1, id2, v2)
}
Try the trtswitch package in your browser
Any scripts or data that you put into this service are public.
trtswitch documentation built on June 8, 2025, 1:45 p.m.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Defines functions match3 qrcpp survsplit hasVariable findInterval3 tssim tsesimpcpp tsegestsim tsegestcpp residuals_phregcpp survfit_phregcpp phregcpp residuals_liferegcpp liferegcpp rmdiff rmest lrtest kmdiff kmest survQuantile nscpp bscpp splineDesigncpp rpsftmcpp logisregcpp ipecpp ipcwcpp
Documented in bscpp findInterval3 kmdiff kmest lrtest nscpp qrcpp rmdiff rmest splineDesigncpp survQuantile survsplit tsegestsim tssim
# Generated by using Rcpp::compileAttributes() -> do not edit by hand
# Generator token: 10BE3573-1514-4C36-9D1C-5A225CD40393
ipcwcpp <- function(data, id = "id", stratum = "", tstart = "tstart", tstop = "tstop", event = "event", treat = "treat", swtrt = "swtrt", swtrt_time = "swtrt_time", base_cov = "", numerator = "", denominator = "", logistic_switching_model = 0L, strata_main_effect_only = 1L, firth = 0L, flic = 0L, ns_df = 3L, stabilized_weights = 1L, trunc = 0, trunc_upper_only = 1L, swtrt_control_only = 1L, alpha = 0.05, ties = "efron", boot = 1L, n_boot = 1000L, seed = NA_integer_) {
.Call(`_trtswitch_ipcwcpp`, data, id, stratum, tstart, tstop, event, treat, swtrt, swtrt_time, base_cov, numerator, denominator, logistic_switching_model, strata_main_effect_only, firth, flic, ns_df, stabilized_weights, trunc, trunc_upper_only, swtrt_control_only, alpha, ties, boot, n_boot, seed)
}
ipecpp <- function(data, id = "id", stratum = "", time = "time", event = "event", treat = "treat", rx = "rx", censor_time = "censor_time", base_cov = "", aft_dist = "weibull", strata_main_effect_only = 1L, low_psi = -2, hi_psi = 2, treat_modifier = 1, recensor = 1L, admin_recensor_only = 1L, autoswitch = 1L, alpha = 0.05, ties = "efron", tol = 1.0e-6, boot = 0L, n_boot = 1000L, seed = NA_integer_) {
.Call(`_trtswitch_ipecpp`, data, id, stratum, time, event, treat, rx, censor_time, base_cov, aft_dist, strata_main_effect_only, low_psi, hi_psi, treat_modifier, recensor, admin_recensor_only, autoswitch, alpha, ties, tol, boot, n_boot, seed)
}
logisregcpp <- function(data, rep = "", event = "event", covariates = "", freq = "", weight = "", offset = "", id = "", link = "logit", robust = 0L, firth = 0L, flic = 0L, plci = 0L, alpha = 0.05, maxiter = 50L, eps = 1.0e-9) {
.Call(`_trtswitch_logisregcpp`, data, rep, event, covariates, freq, weight, offset, id, link, robust, firth, flic, plci, alpha, maxiter, eps)
}
rpsftmcpp <- function(data, id = "id", stratum = "", time = "time", event = "event", treat = "treat", rx = "rx", censor_time = "censor_time", base_cov = "", psi_test = "logrank", aft_dist = "weibull", strata_main_effect_only = 1L, low_psi = -2, hi_psi = 2, n_eval_z = 101L, treat_modifier = 1, recensor = 1L, admin_recensor_only = 1L, autoswitch = 1L, gridsearch = 0L, alpha = 0.05, ties = "efron", tol = 1.0e-6, boot = 0L, n_boot = 1000L, seed = NA_integer_) {
.Call(`_trtswitch_rpsftmcpp`, data, id, stratum, time, event, treat, rx, censor_time, base_cov, psi_test, aft_dist, strata_main_effect_only, low_psi, hi_psi, n_eval_z, treat_modifier, recensor, admin_recensor_only, autoswitch, gridsearch, alpha, ties, tol, boot, n_boot, seed)
}
#' @title B-Spline Design Matrix
#' @description Computes the design matrix for B-splines based on the
#' specified \code{knots} and evaluated at the values in \code{x}.
#'
#' @param knots A numeric vector specifying the positions of the knots,
#' including both boundary and internal knots.
#' @param x A numeric vector of values where the B-spline functions
#' or their derivatives will be evaluated. The values of \code{x}
#' must lie within the range of the "inner" knots, i.e., between
#' \code{knots[ord]} and \code{knots[length(knots) - (ord - 1)]}.
#' @param ord A positive integer indicating the order of the B-spline.
#' This corresponds to the number of coefficients in each piecewise
#' polynomial segment, where \code{ord = degree + 1}.
#' @param derivs An integer vector specifying the order of derivatives
#' to be evaluated at the corresponding \code{x} values. Each value
#' must be between \code{0} and \code{ord - 1}, and the vector is
#' conceptually recycled to match the length of \code{x}.
#' The default is \code{0}, meaning the B-spline functions themselves
#' are evaluated.
#'
#' @return A matrix with dimensions
#' \code{c(length(x), length(knots) - ord)}. Each row corresponds
#' to a value in \code{x} and contains the coefficients of the
#' B-splines, or the specified derivatives, as defined by the
#' \code{knots} and evaluated at that particular value of \code{x}.
#' The total number of B-splines is \code{length(knots) - ord},
#' with each B-spline defined by a set of \code{ord} consecutive knots.
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @examples
#'
#' splineDesigncpp(knots = 1:10, x = 4:7)
#' splineDesigncpp(knots = 1:10, x = 4:7, derivs = 1)
#'
#' @export
splineDesigncpp <- function(knots = NA_real_, x = NA_real_, ord = 4L, derivs = as.integer( c(0))) {
.Call(`_trtswitch_splineDesigncpp`, knots, x, ord, derivs)
}
#' @title B-Spline Basis for Polynomial Splines
#' @description Computes the B-spline basis matrix for a given polynomial
#' spline.
#'
#' @param x A numeric vector representing the predictor variable.
#' @param df Degrees of freedom, specifying the number of columns in the
#' basis matrix. If \code{df} is provided, the function automatically
#' selects \code{df - degree - intercept} internal knots based on
#' appropriate quantiles of \code{x}, ignoring any missing values.
#' @param knots A numeric vector specifying the internal breakpoints
#' that define the spline. If not provided, \code{df} must be specified.
#' @param degree An integer specifying the degree of the piecewise
#' polynomial. The default value is \code{3}, which corresponds to
#' cubic splines.
#' @param intercept A logical value indicating whether to include an
#' intercept in the basis. The default is \code{FALSE}.
#' @param boundary_knots A numeric vector of length 2 specifying the
#' boundary points where the B-spline basis should be anchored.
#' If not supplied, the default is the range of non-missing values
#' in \code{x}.
#' @param warn_outside A logical value indicating whether a warning
#' should be issued if any values of \code{x} fall outside the
#' specified boundary knots.
#'
#' @return A matrix with dimensions \code{c(length(x), df)}. If
#' \code{df} is provided, the matrix will have \code{df} columns.
#' Alternatively, if \code{knots} are supplied, the number of columns
#' will be \code{length(knots) + degree + intercept}. The matrix
#' contains attributes that correspond to the arguments
#' passed to the \code{bscpp} function.
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @examples
#'
#' bscpp(women$height, df = 5)
#'
#' @export
bscpp <- function(x = NA_real_, df = NA_integer_, knots = NA_real_, degree = 3L, intercept = 0L, boundary_knots = NA_real_, warn_outside = 1L) {
.Call(`_trtswitch_bscpp`, x, df, knots, degree, intercept, boundary_knots, warn_outside)
}
#' @title Natural Cubic Spline Basis
#' @description Computes the B-spline basis matrix for a natural cubic
#' spline.
#'
#' @param x A numeric vector representing the predictor variable.
#' Missing values are allowed.
#' @param df Degrees of freedom, specifying the number of columns in
#' the basis matrix. If \code{df} is provided, the function selects
#' \code{df - 1 - intercept} internal knots based on appropriate
#' quantiles of \code{x}, ignoring any missing values.
#' @param knots A numeric vector specifying the internal breakpoints
#' that define the spline. If provided, the number of degrees of
#' freedom will be determined by the length of \code{knots}.
#' @param intercept A logical value indicating whether to include an
#' intercept in the basis. The default is \code{FALSE}.
#' @param boundary_knots A numeric vector of length 2 specifying the
#' boundary points where the natural boundary conditions are applied
#' and the B-spline basis is anchored. If not supplied, the default
#' is the range of non-missing values in \code{x}.
#'
#' @return A matrix with dimensions \code{c(length(x), df)}, where
#' \code{df} is either provided directly or computed as
#' \code{length(knots) + 1 + intercept} when \code{knots} are supplied.
#' The matrix contains attributes that correspond to the arguments
#' passed to the \code{nscpp} function.
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @examples
#'
#' nscpp(women$height, df = 5)
#'
#' @export
nscpp <- function(x = NA_real_, df = NA_integer_, knots = NA_real_, intercept = 0L, boundary_knots = NA_real_) {
.Call(`_trtswitch_nscpp`, x, df, knots, intercept, boundary_knots)
}
#' @title Brookmeyer-Crowley Confidence Interval for Quantiles of
#' Right-Censored Time-to-Event Data
#' @description Obtains the Brookmeyer-Crowley confidence
#' interval for quantiles of right-censored time-to-event data.
#'
#' @param time The vector of possibly right-censored survival times.
#' @param event The vector of event indicators.
#' @param cilevel The confidence interval level. Defaults to 0.95.
#' @param transform The transformation of the survival function to use
#' to construct the confidence interval. Options include
#' "linear" (alternatively "plain"), "log",
#' "loglog" (alternatively "log-log" or "cloglog"),
#' "asinsqrt" (alternatively "asin" or "arcsin"), and "logit".
#' Defaults to "loglog".
#'
#' @param probs The vector of probabilities to calculate the quantiles.
#' Defaults to c(0.25, 0.5, 0.75).
#'
#' @return A data frame containing the estimated quantile and
#' confidence interval corresponding to each specified probability.
#' It includes the following variables:
#'
#' * \code{prob}: The probability to calculate the quantile.
#'
#' * \code{quantile}: The estimated quantile.
#'
#' * \code{lower}: The lower limit of the confidence interval.
#'
#' * \code{upper}: The upper limit of the confidence interval.
#'
#' * \code{cilevel}: The confidence interval level.
#'
#' * \code{transform}: The transformation of the survival function to use
#' to construct the confidence interval.
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @examples
#'
#' survQuantile(
#' time = c(33.7, 3.9, 10.5, 5.4, 19.5, 23.8, 7.9, 16.9, 16.6,
#' 33.7, 17.1, 7.9, 10.5, 38),
#' event = c(0, 1, 1, 1, 1, 0, 1, 0, 0, 0, 0, 0, 1, 1),
#' probs = c(0.25, 0.5, 0.75))
#'
#' @export
survQuantile <- function(time = NA_real_, event = NA_real_, cilevel = 0.95, transform = "loglog", probs = NA_real_) {
.Call(`_trtswitch_survQuantile`, time, event, cilevel, transform, probs)
}
#' @title Kaplan-Meier Estimates of Survival Curve
#' @description Obtains the Kaplan-Meier estimates of the survival curve.
#'
#' @param data The input data frame that contains the following variables:
#'
#' * \code{rep}: The replication for by-group processing.
#'
#' * \code{stratum}: The stratum.
#'
#' * \code{time}: The possibly right-censored survival time.
#'
#' * \code{event}: The event indicator.
#'
#' @param rep The name(s) of the replication variable(s) in the input data.
#' @param stratum The name(s) of the stratum variable(s) in the input data.
#' @param time The name of the time variable in the input data.
#' @param event The name of the event variable in the input data.
#' @param conftype The type of the confidence interval. One of "none",
#' "plain", "log", "log-log" (the default), or "arcsin".
#' The arcsin option bases the intervals on asin(sqrt(survival)).
#' @param conflev The level of the two-sided confidence interval for
#' the survival probabilities. Defaults to 0.95.
#' @param keep_censor Whether to retain the censoring time in the output
#' data frame.
#'
#' @return A data frame with the following variables:
#'
#' * \code{size}: The number of subjects in the stratum.
#'
#' * \code{time}: The event time.
#'
#' * \code{nrisk}: The number of subjects at risk.
#'
#' * \code{nevent}: The number of subjects having the event.
#'
#' * \code{ncensor}: The number of censored subjects.
#'
#' * \code{survival}: The Kaplan-Meier estimate of the survival probability.
#'
#' * \code{stderr}: The standard error of the estimated survival
#' probability based on the Greendwood formula.
#'
#' * \code{lower}: The lower bound of confidence interval if requested.
#'
#' * \code{upper}: The upper bound of confidence interval if requested.
#'
#' * \code{conflev}: The level of confidence interval if requested.
#'
#' * \code{conftype}: The type of confidence interval if requested.
#'
#' * \code{stratum}: The stratum.
#'
#' * \code{rep}: The replication.
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @examples
#'
#' kmest(data = aml, stratum = "x", time = "time", event = "status")
#'
#' @export
kmest <- function(data, rep = "", stratum = "", time = "time", event = "event", conftype = "log-log", conflev = 0.95, keep_censor = 0L) {
.Call(`_trtswitch_kmest`, data, rep, stratum, time, event, conftype, conflev, keep_censor)
}
#' @title Estimate of Milestone Survival Difference
#' @description Obtains the estimate of milestone survival difference
#' between two treatment groups.
#'
#' @param data The input data frame that contains the following variables:
#'
#' * \code{rep}: The replication for by-group processing.
#'
#' * \code{stratum}: The stratum.
#'
#' * \code{treat}: The treatment.
#'
#' * \code{time}: The possibly right-censored survival time.
#'
#' * \code{event}: The event indicator.
#'
#' @param rep The name of the replication variable in the input data.
#' @param stratum The name of the stratum variable in the input data.
#' @param treat The name of the treatment variable in the input data.
#' @param time The name of the time variable in the input data.
#' @param event The name of the event variable in the input data.
#' @param milestone The milestone time at which to calculate the
#' survival probability.
#' @param survDiffH0 The difference in milestone survival probabilities
#' under the null hypothesis. Defaults to 0 for superiority test.
#' @param conflev The level of the two-sided confidence interval for
#' the difference in milestone survival probabilities. Defaults to 0.95.
#'
#' @return A data frame with the following variables:
#'
#' * \code{rep}: The replication.
#'
#' * \code{milestone}: The milestone time relative to randomization.
#'
#' * \code{survDiffH0}: The difference in milestone survival probabilities
#' under the null hypothesis.
#'
#' * \code{surv1}: The estimated milestone survival probability for
#' the treatment group.
#'
#' * \code{surv2}: The estimated milestone survival probability for
#' the control group.
#'
#' * \code{survDiff}: The estimated difference in milestone survival
#' probabilities.
#'
#' * \code{vsurv1}: The variance for surv1.
#'
#' * \code{vsurv2}: The variance for surv2.
#'
#' * \code{vsurvDiff}: The variance for survDiff.
#'
#' * \code{survDiffZ}: The Z-statistic value.
#'
#' * \code{survDiffPValue}: The one-sided p-value.
#'
#' * \code{lower}: The lower bound of confidence interval.
#'
#' * \code{upper}: The upper bound of confidence interval.
#'
#' * \code{conflev}: The level of confidence interval.
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @examples
#'
#' df <- kmdiff(data = rawdata, rep = "iterationNumber",
#' stratum = "stratum", treat = "treatmentGroup",
#' time = "timeUnderObservation", event = "event",
#' milestone = 12)
#' head(df)
#'
#' @export
kmdiff <- function(data, rep = "", stratum = "", treat = "treat", time = "time", event = "event", milestone = NA_real_, survDiffH0 = 0, conflev = 0.95) {
.Call(`_trtswitch_kmdiff`, data, rep, stratum, treat, time, event, milestone, survDiffH0, conflev)
}
#' @title Log-Rank Test of Survival Curve Difference
#' @description Obtains the log-rank test using the Fleming-Harrington
#' family of weights.
#'
#' @param data The input data frame that contains the following variables:
#'
#' * \code{rep}: The replication for by-group processing.
#'
#' * \code{stratum}: The stratum.
#'
#' * \code{treat}: The treatment.
#'
#' * \code{time}: The possibly right-censored survival time.
#'
#' * \code{event}: The event indicator.
#'
#' @param rep The name of the replication variable in the input data.
#' @param stratum The name of the stratum variable in the input data.
#' @param treat The name of the treatment variable in the input data.
#' @param time The name of the time variable in the input data.
#' @param event The name of the event variable in the input data.
#' @param rho1 The first parameter of the Fleming-Harrington family of
#' weighted log-rank test. Defaults to 0 for conventional log-rank test.
#' @param rho2 The second parameter of the Fleming-Harrington family of
#' weighted log-rank test. Defaults to 0 for conventional log-rank test.
#'
#' @return A data frame with the following variables:
#'
#' * \code{uscore}: The numerator of the log-rank test statistic.
#'
#' * \code{vscore}: The variance of the log-rank score test statistic.
#'
#' * \code{logRankZ}: The Z-statistic value.
#'
#' * \code{logRankPValue}: The one-sided p-value.
#'
#' * \code{rho1}: The first parameter of the Fleming-Harrington weights.
#'
#' * \code{rho2}: The second parameter of the Fleming-Harrington weights.
#'
#' * \code{rep}: The replication.
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @examples
#'
#' df <- lrtest(data = rawdata, rep = "iterationNumber",
#' stratum = "stratum", treat = "treatmentGroup",
#' time = "timeUnderObservation", event = "event",
#' rho1 = 0.5, rho2 = 0)
#' head(df)
#'
#' @export
lrtest <- function(data, rep = "", stratum = "", treat = "treat", time = "time", event = "event", rho1 = 0, rho2 = 0) {
.Call(`_trtswitch_lrtest`, data, rep, stratum, treat, time, event, rho1, rho2)
}
#' @title Estimate of Restricted Mean Survival Time
#' @description Obtains the estimate of restricted means survival time
#' for each stratum.
#'
#' @param data The input data frame that contains the following variables:
#'
#' * \code{rep}: The replication for by-group processing.
#'
#' * \code{stratum}: The stratum.
#'
#' * \code{time}: The possibly right-censored survival time.
#'
#' * \code{event}: The event indicator.
#'
#' @param rep The name of the replication variable in the input data.
#' @param stratum The name of the stratum variable in the input data.
#' @param time The name of the time variable in the input data.
#' @param event The name of the event variable in the input data.
#' @param milestone The milestone time at which to calculate the
#' restricted mean survival time.
#' @param conflev The level of the two-sided confidence interval for
#' the survival probabilities. Defaults to 0.95.
#' @param biascorrection Whether to apply bias correction for the
#' variance estimate. Defaults to no bias correction.
#'
#' @return A data frame with the following variables:
#'
#' * \code{rep}: The replication.
#'
#' * \code{stratum}: The stratum variable.
#'
#' * \code{size}: The number of subjects in the stratum.
#'
#' * \code{milestone}: The milestone time relative to randomization.
#'
#' * \code{rmst}: The estimate of restricted mean survival time.
#'
#' * \code{stderr}: The standard error of the estimated rmst.
#'
#' * \code{lower}: The lower bound of confidence interval if requested.
#'
#' * \code{upper}: The upper bound of confidence interval if requested.
#'
#' * \code{conflev}: The level of confidence interval if requested.
#'
#' * \code{biascorrection}: Whether to apply bias correction for the
#' variance estimate.
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @examples
#'
#' rmest(data = aml, stratum = "x",
#' time = "time", event = "status", milestone = 24)
#'
#' @export
rmest <- function(data, rep = "", stratum = "", time = "time", event = "event", milestone = NA_real_, conflev = 0.95, biascorrection = 0L) {
.Call(`_trtswitch_rmest`, data, rep, stratum, time, event, milestone, conflev, biascorrection)
}
#' @title Estimate of Restricted Mean Survival Time Difference
#' @description Obtains the estimate of restricted mean survival time
#' difference between two treatment groups.
#'
#' @param data The input data frame that contains the following variables:
#'
#' * \code{rep}: The replication for by-group processing.
#'
#' * \code{stratum}: The stratum.
#'
#' * \code{treat}: The treatment.
#'
#' * \code{time}: The possibly right-censored survival time.
#'
#' * \code{event}: The event indicator.
#'
#' @param rep The name of the replication variable in the input data.
#' @param stratum The name of the stratum variable in the input data.
#' @param treat The name of the treatment variable in the input data.
#' @param time The name of the time variable in the input data.
#' @param event The name of the event variable in the input data.
#' @param milestone The milestone time at which to calculate the
#' restricted mean survival time.
#' @param rmstDiffH0 The difference in restricted mean survival times
#' under the null hypothesis. Defaults to 0 for superiority test.
#' @param conflev The level of the two-sided confidence interval for
#' the difference in restricted mean survival times. Defaults to 0.95.
#' @param biascorrection Whether to apply bias correction for the
#' variance estimate of individual restricted mean survival times.
#' Defaults to no bias correction.
#'
#' @return A data frame with the following variables:
#'
#' * \code{rep}: The replication number.
#'
#' * \code{milestone}: The milestone time relative to randomization.
#'
#' * \code{rmstDiffH0}: The difference in restricted mean survival times
#' under the null hypothesis.
#'
#' * \code{rmst1}: The estimated restricted mean survival time for
#' the treatment group.
#'
#' * \code{rmst2}: The estimated restricted mean survival time for
#' the control group.
#'
#' * \code{rmstDiff}: The estimated difference in restricted mean
#' survival times.
#'
#' * \code{vrmst1}: The variance for rmst1.
#'
#' * \code{vrmst2}: The variance for rmst2.
#'
#' * \code{vrmstDiff}: The variance for rmstDiff.
#'
#' * \code{rmstDiffZ}: The Z-statistic value.
#'
#' * \code{rmstDiffPValue}: The one-sided p-value.
#'
#' * \code{lower}: The lower bound of confidence interval.
#'
#' * \code{upper}: The upper bound of confidence interval.
#'
#' * \code{conflev}: The level of confidence interval.
#'
#' * \code{biascorrection}: Whether to apply bias correction for the
#' variance estimate of individual restricted mean survival times.
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @examples
#'
#' df <- rmdiff(data = rawdata, rep = "iterationNumber",
#' stratum = "stratum", treat = "treatmentGroup",
#' time = "timeUnderObservation", event = "event",
#' milestone = 12)
#' head(df)
#'
#' @export
rmdiff <- function(data, rep = "", stratum = "", treat = "treat", time = "time", event = "event", milestone = NA_real_, rmstDiffH0 = 0, conflev = 0.95, biascorrection = 0L) {
.Call(`_trtswitch_rmdiff`, data, rep, stratum, treat, time, event, milestone, rmstDiffH0, conflev, biascorrection)
}
liferegcpp <- function(data, rep = "", stratum = "", time = "time", time2 = "", event = "event", covariates = "", weight = "", offset = "", id = "", dist = "weibull", robust = 0L, plci = 0L, alpha = 0.05, maxiter = 50L, eps = 1.0e-9) {
.Call(`_trtswitch_liferegcpp`, data, rep, stratum, time, time2, event, covariates, weight, offset, id, dist, robust, plci, alpha, maxiter, eps)
}
residuals_liferegcpp <- function(beta, vbeta, data, stratum = "", time = "time", time2 = "", event = "event", covariates = "", weight = "", offset = "", id = "", dist = "weibull", type = "response", collapse = 0L, weighted = 0L) {
.Call(`_trtswitch_residuals_liferegcpp`, beta, vbeta, data, stratum, time, time2, event, covariates, weight, offset, id, dist, type, collapse, weighted)
}
phregcpp <- function(data, rep = "", stratum = "", time = "time", time2 = "", event = "event", covariates = "", weight = "", offset = "", id = "", ties = "efron", robust = 0L, est_basehaz = 1L, est_resid = 1L, firth = 0L, plci = 0L, alpha = 0.05, maxiter = 50L, eps = 1.0e-9) {
.Call(`_trtswitch_phregcpp`, data, rep, stratum, time, time2, event, covariates, weight, offset, id, ties, robust, est_basehaz, est_resid, firth, plci, alpha, maxiter, eps)
}
survfit_phregcpp <- function(p, beta, vbeta, basehaz, newdata, covariates = "", stratum = "", offset = "", id = "", tstart = "", tstop = "", sefit = 1L, conftype = "log-log", conflev = 0.95) {
.Call(`_trtswitch_survfit_phregcpp`, p, beta, vbeta, basehaz, newdata, covariates, stratum, offset, id, tstart, tstop, sefit, conftype, conflev)
}
residuals_phregcpp <- function(p, beta, vbeta, resmart, data, stratum = "", time = "time", time2 = "", event = "event", covariates = "", weight = "", offset = "", id = "", ties = "efron", type = "schoenfeld", collapse = 0L, weighted = 0L) {
.Call(`_trtswitch_residuals_phregcpp`, p, beta, vbeta, resmart, data, stratum, time, time2, event, covariates, weight, offset, id, ties, type, collapse, weighted)
}
tsegestcpp <- function(data, id = "id", stratum = "", tstart = "tstart", tstop = "tstop", event = "event", treat = "treat", censor_time = "censor_time", pd = "pd", pd_time = "pd_time", swtrt = "swtrt", swtrt_time = "swtrt_time", base_cov = "", conf_cov = "", low_psi = -2, hi_psi = 2, n_eval_z = 101L, strata_main_effect_only = 1L, firth = 0L, flic = 0L, recensor = 1L, admin_recensor_only = 1L, swtrt_control_only = 1L, gridsearch = 0L, alpha = 0.05, ties = "efron", tol = 1.0e-6, offset = 1, boot = 1L, n_boot = 1000L, seed = NA_integer_) {
.Call(`_trtswitch_tsegestcpp`, data, id, stratum, tstart, tstop, event, treat, censor_time, pd, pd_time, swtrt, swtrt_time, base_cov, conf_cov, low_psi, hi_psi, n_eval_z, strata_main_effect_only, firth, flic, recensor, admin_recensor_only, swtrt_control_only, gridsearch, alpha, ties, tol, offset, boot, n_boot, seed)
}
#' @title Simulate Survival Data for Two-Stage Estimation Method Using
#' g-estimation
#' @description Obtains the simulated data for baseline prognosis,
#' disease progression, treatment switching, death, and
#' time-dependent covariates.
#'
#' @param n The total sample size for two treatment arms combined.
#' @param allocation1 The number of subjects in the active treatment group
#' in a randomization block.
#' @param allocation2 The number of subjects in the control group in
#' a randomization block.
#' @param pbprog The probability of having poor prognosis at baseline.
#' @param trtlghr The treatment effect in terms of log hazard ratio.
#' @param bprogsl The poor prognosis effect in terms of log hazard ratio.
#' @param shape1 The shape parameter for the Weibull event distribution
#' for the first component.
#' @param scale1 The scale parameter for the Weibull event distribution
#' for the first component.
#' @param shape2 The shape parameter for the Weibull event distribution
#' for the second component.
#' @param scale2 The scale parameter for the Weibull event distribution
#' for the second component.
#' @param pmix The mixing probability of the first component Weibull
#' distribution.
#' @param admin The administrative censoring time.
#' @param pcatnotrtbprog The probability of developing metastatic disease
#' on control treatment with poor baseline prognosis.
#' @param pcattrtbprog The probability of developing metastatic disease
#' on active treatment with poor baseline prognosis.
#' @param pcatnotrt The probability of developing metastatic disease
#' on control treatment with good baseline prognosis.
#' @param pcattrt The probability of developing metastatic disease
#' on active treatment with good baseline prognosis.
#' @param catmult The impact of metastatic disease on shortening remaining
#' survival time.
#' @param tdxo Whether treatment crossover depends on time-dependent
#' covariates between disease progression and treatment switching.
#' @param ppoor The probability of switching for poor baseline prognosis
#' with no metastatic disease.
#' @param pgood The probability of switching for good baseline prognosis
#' with no metastatic disease.
#' @param ppoormet The probability of switching for poor baseline prognosis
#' after developing metastatic disease.
#' @param pgoodmet The probability of switching for good baseline prognosis
#' after developing metastatic disease.
#' @param xomult The direct effect of crossover on extending remaining
#' survival time.
#' @param milestone The milestone to calculate restricted mean survival
#' time.
#' @param outputRawDataset Whether to output the raw data set.
#' @param seed The seed to reproduce the simulation results.
#' The seed from the environment will be used if left unspecified.
#'
#' @return A list with two data frames.
#'
#' * \code{sumdata}: A summary data frame with the following variables:
#'
#' - \code{simtrueconstmean}: The true control group restricted mean
#' survival time (RMST).
#'
#' - \code{simtrueconstlb}: The lower bound for control group RMST.
#'
#' - \code{simtrueconstub}: The upper bound for control group RMST.
#'
#' - \code{simtrueconstse}: The standard error for control group RMST.
#'
#' - \code{simtrueexpstmean}: The true experimental group restricted
#' mean survival time (RMST).
#'
#' - \code{simtrueexpstlb}: The lower bound for experimental group RMST.
#'
#' - \code{simtrueexpstub}: The upper bound for experimental group RMST.
#'
#' - \code{simtrueexpstse}: The standard error for experimental group
#' RMST.
#'
#' - \code{simtrue_coxwbprog_hr}: The treatment hazard ratio from the
#' Cox model adjusting for baseline prognosis.
#'
#' - \code{simtrue_cox_hr}: The treatment hazard ratio from the Cox
#' model without adjusting for baseline prognosis.
#'
#' * \code{paneldata}: A counting process style subject-level data frame
#' with the following variables:
#'
#' - \code{id}: The subject ID.
#'
#' - \code{trtrand}: The randomized treatment arm.
#'
#' - \code{bprog}: Whether the patient had poor baseline prognosis.
#'
#' - \code{tstart}: The left end of time interval.
#'
#' - \code{tstop}: The right end of time interval.
#'
#' - \code{event}: Whether the patient died at the end of the interval.
#'
#' - \code{timeOS}: The observed survival time.
#'
#' - \code{died}: Whether the patient died during the study.
#'
#' - \code{progressed}: Whether the patient had disease progression.
#'
#' - \code{timePFSobs}: The observed time of disease progression at
#' regular scheduled visits.
#'
#' - \code{progtdc}: The time-dependent covariate for progression.
#'
#' - \code{catevent}: Whether the patient developed metastatic disease.
#'
#' - \code{cattime}: When the patient developed metastatic disease.
#'
#' - \code{cattdc}: The time-dependent covariate for cat event.
#'
#' - \code{catlag}: The lagged value of \code{cattdc}.
#'
#' - \code{xo}: Whether the patient switched treatment.
#'
#' - \code{xotime}: When the patient switched treatment.
#'
#' - \code{xotdc}: The time-dependent covariate for treatment
#' switching.
#'
#' - \code{xotime_upper}: The upper bound of treatment switching time.
#'
#' - \code{censor_time}: The administrative censoring time.
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @references
#' NR Latimer, IR White, K Tilling, and U Siebert.
#' Improved two-stage estimation to adjust for treatment switching in
#' randomised trials: g-estimation to address time-dependent confounding.
#' Statistical Methods in Medical Research. 2020;29(10):2900-2918.
#'
#' @examples
#'
#' sim1 <- tsegestsim(
#' n = 500, allocation1 = 2, allocation2 = 1, pbprog = 0.5,
#' trtlghr = -0.5, bprogsl = 0.3, shape1 = 1.8,
#' scale1 = 360, shape2 = 1.7, scale2 = 688,
#' pmix = 0.5, admin = 5000, pcatnotrtbprog = 0.5,
#' pcattrtbprog = 0.25, pcatnotrt = 0.2, pcattrt = 0.1,
#' catmult = 0.5, tdxo = 1, ppoor = 0.1, pgood = 0.04,
#' ppoormet = 0.4, pgoodmet = 0.2, xomult = 1.4188308,
#' milestone = 546, outputRawDataset = 1, seed = 2000)
#'
#' @export
tsegestsim <- function(n = 500L, allocation1 = 2L, allocation2 = 1L, pbprog = 0.5, trtlghr = -0.5, bprogsl = 0.3, shape1 = 1.8, scale1 = 360, shape2 = 1.7, scale2 = 688, pmix = 0.5, admin = 5000, pcatnotrtbprog = 0.5, pcattrtbprog = 0.25, pcatnotrt = 0.2, pcattrt = 0.1, catmult = 0.5, tdxo = 1, ppoor = 0.1, pgood = 0.04, ppoormet = 0.4, pgoodmet = 0.2, xomult = 1.4188308, milestone = 546, outputRawDataset = 1L, seed = NA_integer_) {
.Call(`_trtswitch_tsegestsim`, n, allocation1, allocation2, pbprog, trtlghr, bprogsl, shape1, scale1, shape2, scale2, pmix, admin, pcatnotrtbprog, pcattrtbprog, pcatnotrt, pcattrt, catmult, tdxo, ppoor, pgood, ppoormet, pgoodmet, xomult, milestone, outputRawDataset, seed)
}
tsesimpcpp <- function(data, id = "id", stratum = "", time = "time", event = "event", treat = "treat", censor_time = "censor_time", pd = "pd", pd_time = "pd_time", swtrt = "swtrt", swtrt_time = "swtrt_time", base_cov = "", base2_cov = "", aft_dist = "weibull", strata_main_effect_only = 1L, recensor = 1L, admin_recensor_only = 1L, swtrt_control_only = 1L, alpha = 0.05, ties = "efron", offset = 1, boot = 1L, n_boot = 1000L, seed = NA_integer_) {
.Call(`_trtswitch_tsesimpcpp`, data, id, stratum, time, event, treat, censor_time, pd, pd_time, swtrt, swtrt_time, base_cov, base2_cov, aft_dist, strata_main_effect_only, recensor, admin_recensor_only, swtrt_control_only, alpha, ties, offset, boot, n_boot, seed)
}
#' @title Simulate Data for Treatment Switching
#' @description Simulates data for studies involving treatment switching,
#' incorporating time-dependent confounding. The generated data can be used
#' to evaluate methods for handling treatment switching in survival analysis.
#'
#' @param tdxo Logical indicator for timing of treatment switching:
#' \itemize{
#' \item 1: Treatment switching can occur at or after disease progression.
#' \item 0: Treatment switching is restricted to the time of disease
#' progression.
#' }
#' @param coxo Logical indicator for arm-specific treatment switching:
#' \itemize{
#' \item 1: Treatment switching occurs only in the control arm.
#' \item 0: Treatment switching is allowed in both arms.
#' }
#' @param p_R Probability of randomization to the experimental arm.
#' @param p_X_1 Probability of poor baseline prognosis in the experimental
#' arm.
#' @param p_X_0 Probability of poor baseline prognosis in the control arm.
#' @param rate_T Baseline hazard rate for time to death.
#' @param beta1 Log hazard ratio for randomized treatment (\code{R}).
#' @param beta2 Log hazard ratio for baseline covariate (\code{X}).
#' @param gamma0 Intercept for the time-dependent covariate model (\code{L}).
#' @param gamma1 Coefficient for lagged treatment switching (\code{Alag})
#' in the \code{L} model.
#' @param gamma2 Coefficient for lagged \code{L} in the \code{L} model.
#' @param gamma3 Coefficient for baseline covariate (\code{X}) in the
#' \code{L} model.
#' @param gamma4 Coefficient for randomized treatment (\code{R}) in the
#' \code{L} model.
#' @param zeta0 Intercept for the disease progression model (\code{Z}).
#' @param zeta1 Coefficient for \code{L} in the \code{Z} model.
#' @param zeta2 Coefficient for baseline covariate (\code{X}) in the
#' \code{Z} model.
#' @param zeta3 Coefficient for randomized treatment (\code{R}) in the
#' \code{Z} model.
#' @param alpha0 Intercept for the treatment switching model (\code{A}).
#' @param alpha1 Coefficient for \code{L} in the \code{A} model.
#' @param alpha2 Coefficient for baseline covariate (\code{X}) in the
#' \code{A} model.
#' @param theta1_1 Negative log time ratio for \code{A} (experimental arm).
#' @param theta1_0 Negative log time ratio for \code{A} (control arm).
#' @param theta2 Negative log time ratio for \code{L}.
#' @param rate_C Hazard rate for random (dropout) censoring.
#' @param followup Number of treatment cycles per subject.
#' @param days Number of days in each treatment cycle.
#' @param n Number of subjects per simulation.
#' @param NSim Number of simulated datasets.
#' @param seed Random seed for reproducibility.
#'
#' @return
#' A list of data frames, each containing simulated longitudinal and
#' event history data with the following variables:
#'
#' * \code{id}: Subject identifier.
#'
#' * \code{trtrand}: Randomized treatment assignment (0 = control,
#' 1 = experimental)
#'
#' * \code{bprog}: Baseline prognosis (0 = good, 1 = poor).
#'
#' * \code{tpoint}: Treatment cycle index.
#'
#' * \code{tstart}: Start day of the treatment cycle.
#'
#' * \code{tstop}: End day of the treatment cycle.
#'
#' * \code{L}: Time-dependent covariate predicting survival and switching;
#' affected by treatment switching.
#'
#' * \code{Llag}: Lagged value of \code{L}.
#'
#' * \code{Z}: Disease progression status at \code{tstop}.
#'
#' * \code{A}: Treatment switching status at \code{tstop}.
#'
#' * \code{Alag}: Lagged value of \code{A}.
#'
#' * \code{Y}: Death indicator at \code{tstop}.
#'
#' * \code{timeOS}: Observed time to death or censoring.
#'
#' * \code{died}: Indicator of death by end of follow-up.
#'
#' * \code{progressed}: Indicator of disease progression by end of follow-up.
#'
#' * \code{timePD}: Observed time to progression or censoring.
#'
#' * \code{xo}: Indicator for whether treatment switching occurred.
#'
#' * \code{xotime}: Time of treatment switching (if applicable).
#'
#' * \code{censor_time}: Administrative censoring time.
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @references
#'
#' Jessica G. Young, and Eric J. Tchetgen Tchetgen.
#' Simulation from a known Cox MSM using standard parametric models
#' for the g-formula.
#' Statistics in Medicine. 2014;33(6):1001-1014.
#'
#' NR Latimer, IR White, K Tilling, and U Siebert.
#' Improved two-stage estimation to adjust for treatment switching in
#' randomised trials: g-estimation to address time-dependent confounding.
#' Statistical Methods in Medical Research. 2020;29(10):2900-2918.
#'
#' Jing Xu, Guohui Liu, and Bingxia Wang.
#' Bias and type I error control in correcting treatment effect
#' for treatment switching using marginal structural models
#' in Phse III oncology trials.
#' Journal of Biopharmaceutical Statistics. 2022;32(6):897-914.
#'
#' @examples
#'
#' simulated.data <- tssim(
#' tdxo = 0, coxo = 0, p_R = 0.5, p_X_1 = 0.3, p_X_0 = 0.3,
#' rate_T = 0.002, beta1 = -0.5, beta2 = 0.3,
#' gamma0 = 0.3, gamma1 = -0.9, gamma2 = 0.7, gamma3 = 1.1, gamma4 = -0.8,
#' zeta0 = -3.5, zeta1 = 0.5, zeta2 = 0.2, zeta3 = -0.4,
#' alpha0 = 0.5, alpha1 = 0.5, alpha2 = 0.4,
#' theta1_1 = -0.4, theta1_0 = -0.4, theta2 = 0.2,
#' rate_C = 0.0000855, followup = 20, days = 30,
#' n = 500, NSim = 100, seed = 314159)
#'
#' @export
tssim <- function(tdxo = 0L, coxo = 1L, p_R = 0.5, p_X_1 = NA_real_, p_X_0 = NA_real_, rate_T = NA_real_, beta1 = NA_real_, beta2 = NA_real_, gamma0 = NA_real_, gamma1 = NA_real_, gamma2 = NA_real_, gamma3 = NA_real_, gamma4 = NA_real_, zeta0 = NA_real_, zeta1 = NA_real_, zeta2 = NA_real_, zeta3 = NA_real_, alpha0 = NA_real_, alpha1 = NA_real_, alpha2 = NA_real_, theta1_1 = NA_real_, theta1_0 = NA_real_, theta2 = NA_real_, rate_C = NA_real_, followup = NA_integer_, days = NA_integer_, n = NA_integer_, NSim = 1000L, seed = NA_integer_) {
.Call(`_trtswitch_tssim`, tdxo, coxo, p_R, p_X_1, p_X_0, rate_T, beta1, beta2, gamma0, gamma1, gamma2, gamma3, gamma4, zeta0, zeta1, zeta2, zeta3, alpha0, alpha1, alpha2, theta1_1, theta1_0, theta2, rate_C, followup, days, n, NSim, seed)
}
#' @title Find Interval Numbers of Indices
#' @description The implementation of \code{findInterval()} in R from
#' Advanced R by Hadley Wickham. Given a vector of non-decreasing
#' breakpoints in v, find the interval containing each element of x; i.e.,
#' if \code{i <- findInterval3(x,v)}, for each index \code{j} in \code{x},
#' \code{v[i[j]] <= x[j] < v[i[j] + 1]}, where \code{v[0] := -Inf},
#' \code{v[N+1] := +Inf}, and \code{N = length(v)}.
#'
#' @param x The numeric vector of interest.
#' @param v The vector of break points.
#' @return A vector of \code{length(x)} with values in \code{0:N} where
#' \code{N = length(v)}.
#'
#' @keywords internal
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @examples
#' x <- 2:18
#' v <- c(5, 10, 15) # create two bins [5,10) and [10,15)
#' cbind(x, findInterval3(x, v))
#'
#' @export
findInterval3 <- function(x, v) {
.Call(`_trtswitch_findInterval3`, x, v)
}
hasVariable <- function(df, varName) {
.Call(`_trtswitch_hasVariable`, df, varName)
}
#' @title Split a survival data set at specified cut points
#' @description For a given survival dataset and specified cut times,
#' each record is split into multiple subrecords at each cut time.
#' The resulting dataset is in counting process format, with each
#' subrecord containing a start time, stop time, and event status.
#' This is adapted from the survsplit.c function from the survival package.
#'
#' @param tstart The starting time of the time interval for
#' counting-process data.
#' @param tstop The stopping time of the time interval for
#' counting-process data.
#' @param cut The vector of cut points.
#'
#' @return A data frame with the following variables:
#'
#' * \code{row}: The row number of the observation in the input data
#' (starting from 0).
#'
#' * \code{start}: The starting time of the resulting subrecord.
#'
#' * \code{end}: The ending time of the resulting subrecord.
#'
#' * \code{censor}: Whether the subrecord lies strictly within a record
#' in the input data (1 for all but the last interval and 0 for the
#' last interval with cutpoint set equal to tstop).
#'
#' * \code{interval}: The interval number derived from cut (starting
#' from 0 if the interval lies to the left of the first cutpoint).
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @keywords internal
#'
#' @examples
#'
#' survsplit(15, 60, c(10, 30, 40))
#'
#' @export
survsplit <- function(tstart, tstop, cut) {
.Call(`_trtswitch_survsplit`, tstart, tstop, cut)
}
#' @title QR Decomposition of a Matrix
#' @description Computes the QR decomposition of a matrix.
#'
#' @param X A numeric matrix whose QR decomposition is to be computed.
#' @param tol The tolerance for detecting linear dependencies in the
#' columns of \code{X}.
#'
#' @details
#' This function performs Householder QR with column pivoting:
#' Given an \eqn{m}-by-\eqn{n} matrix \eqn{A} with \eqn{m \geq n},
#' the following algorithm computes \eqn{r = \textrm{rank}(A)} and
#' the factorization \eqn{Q^T A P} equal to
#' \tabular{ccccc}{
#' | \tab \eqn{R_{11}} \tab \eqn{R_{12}} \tab | \tab \eqn{r} \cr
#' | \tab 0 \tab 0 \tab | \tab \eqn{m-r} \cr
#' \tab \eqn{r} \tab \eqn{n-r} \tab \tab
#' }
#' with \eqn{Q = H_1 \cdots H_r} and \eqn{P = P_1 \cdots P_r}.
#' The upper triangular part of \eqn{A}
#' is overwritten by the upper triangular part of \eqn{R} and
#' components \eqn{(j+1):m} of
#' the \eqn{j}th Householder vector are stored in \eqn{A((j+1):m, j)}.
#' The permutation \eqn{P} is encoded in an integer vector \code{pivot}.
#'
#' @return A list with the following components:
#'
#' * \code{qr}: A matrix with the same dimensions as \code{X}. The upper
#' triangle contains the \code{R} of the decomposition and the lower
#' triangle contains Householder vectors (stored in compact form).
#'
#' * \code{rank}: The rank of \code{X} as computed by the decomposition.
#'
#' * \code{pivot}: The column permutation for the pivoting strategy used
#' during the decomposition.
#'
#' * \code{Q}: The complete \eqn{m}-by-\eqn{m} orthogonal matrix \eqn{Q}.
#'
#' * \code{R}: The complete \eqn{m}-by-\eqn{n} upper triangular
#' matrix \eqn{R}.
#'
#' @author Kaifeng Lu, \email{kaifenglu@@gmail.com}
#'
#' @references
#' Gene N. Golub and Charles F. Van Loan.
#' Matrix Computations, second edition. Baltimore, Maryland:
#' The John Hopkins University Press, 1989, p.235.
#'
#' @examples
#'
#' hilbert <- function(n) { i <- 1:n; 1 / outer(i - 1, i, `+`) }
#' h9 <- hilbert(9)
#' qrcpp(h9)
#'
#' @export
qrcpp <- function(X, tol = 1e-12) {
.Call(`_trtswitch_qrcpp`, X, tol)
}
match3 <- function(id1, v1, id2, v2) {
.Call(`_trtswitch_match3`, id1, v1, id2, v2)
}
Try the trtswitch package in your browser
Any scripts or data that you put into this service are public.
R Package Documentation
Browse R Packages
We want your feedback!
Note that we can't provide technical support on individual packages. You should contact the package authors for that.
Embedding an R snippet on your website
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.