| Title: | Confidence Intervals for Exceedance Probability |
|---|---|
| Description: | Computes confidence intervals for the exceedance probability of normally distributed estimators. Currently supports general linear models, with a beta function for Cox models. Please see Segal (2019) <doi:10.1080/00031305.2019.1678521> for more information. |
| Authors: | Brian D. Segal [aut, cre] |
| Maintainer: | Brian D. Segal <[email protected]> |
| License: | GPL (>= 3) |
| Version: | 0.0.2 |
| Built: | 2025-03-07 06:35:33 UTC |
| Source: | https://github.com/bdsegal/exceedprob |
This function obtains point estimates for the exceedance probability for Cox model parameters.
CoxEP(data, surv_formula, j, cutoff, n, m, lower_tail)CoxEP(data, surv_formula, j, cutoff, n, m, lower_tail)
data |
survival data (data.frame) |
surv_formula |
Survival formula |
j |
Index of parameter for which the exceedance probability is obtained |
cutoff |
Cutoff values (scalar or vector) |
n |
Number of observations in the original study |
m |
Number of observations in the replication study (defaults to n if NULL) |
lower_tail |
If TRUE, reports lower tail probabilities |
point (scalar or vector) Point estimate of exceedance probability
This function obtains confidence intervals for exceedance probability
exceedProb(cutoff, theta_hat, sd_hat, alpha, d, n, m, interval = c(-100, 100), lower_tail = FALSE)exceedProb(cutoff, theta_hat, sd_hat, alpha, d, n, m, interval = c(-100, 100), lower_tail = FALSE)
cutoff |
Cutoff values (scalar or vector) |
theta_hat |
Point estimate for the parameter of interest |
sd_hat |
Estimated standard deviation for the parameter of interest (Note: not the standard error) |
alpha |
Significance level |
d |
Number of parameters in the general linear model |
n |
Number of observations in the initial study |
m |
Number of observations in the replication study |
interval |
Interval within which to search for roots |
lower_tail |
If TRUE, reports lower tail probabilities |
ep Exceedance probability with confidence intervals
library(exceedProb) # Sample mean ------------------------------------------------------- n <- 100 x <- rnorm(n = n) theta_hat <- mean(x) sd_hat <- sd(x) cutoff <- seq(from = theta_hat - 0.5, to = theta_hat + 0.5, by = 0.1) exceedProb(cutoff = cutoff, theta_hat = theta_hat, sd_hat = sd_hat, alpha = 0.05, d = 1, n = n, m = n) # Linear regression ------------------------------------------------- n <- 100 beta <- c(1, 2) x <-runif(n = n, min = 0, max = 10) y <- rnorm(n = n, mean = cbind(1, x) %*% beta, sd = 1) j <- 2 fit <- lm(y ~ x) theta_hat <- coef(fit)[j] sd_hat <- sqrt(n * vcov(fit)[j, j]) cutoff <- seq(from = theta_hat - 0.5, to = theta_hat + 0.5, by = 0.1) exceedProb(cutoff = cutoff, theta_hat = theta_hat, sd_hat = sd_hat, alpha = 0.05, d = length(beta), n = n, m = n)library(exceedProb) # Sample mean ------------------------------------------------------- n <- 100 x <- rnorm(n = n) theta_hat <- mean(x) sd_hat <- sd(x) cutoff <- seq(from = theta_hat - 0.5, to = theta_hat + 0.5, by = 0.1) exceedProb(cutoff = cutoff, theta_hat = theta_hat, sd_hat = sd_hat, alpha = 0.05, d = 1, n = n, m = n) # Linear regression ------------------------------------------------- n <- 100 beta <- c(1, 2) x <-runif(n = n, min = 0, max = 10) y <- rnorm(n = n, mean = cbind(1, x) %*% beta, sd = 1) j <- 2 fit <- lm(y ~ x) theta_hat <- coef(fit)[j] sd_hat <- sqrt(n * vcov(fit)[j, j]) cutoff <- seq(from = theta_hat - 0.5, to = theta_hat + 0.5, by = 0.1) exceedProb(cutoff = cutoff, theta_hat = theta_hat, sd_hat = sd_hat, alpha = 0.05, d = length(beta), n = n, m = n)
This function obtains nonparametric bootstrap percentile confidence for Cox model parameters. Beta version.
exceedProbCoxBoot(data, cox_fit, j, alpha, R, cutoff = NULL, m = NULL, lower_tail = FALSE, sim = "model")exceedProbCoxBoot(data, cox_fit, j, alpha, R, cutoff = NULL, m = NULL, lower_tail = FALSE, sim = "model")
data |
survival data (data.frame) |
cox_fit |
(coxph.object) A fitted Cox model |
j |
Index of parameter for which the exceedance probability is obtained |
alpha |
Significance level |
R |
Number of bootstrap resamples |
cutoff |
Cutoff values (scalar or vector if supplied, otherwise set to +/- 0.5 of theta_hat) |
m |
Number of observations in the replication study (defaults to n if NULL) |
lower_tail |
If TRUE, reports lower tail probabilities; otherwise reports upper tail probabilities |
sim |
type of simulation, input to boot::censboot |
ep Exceedance probability with confidence intervals
library(exceedProb) library(survival) # Cox model ------------------------------------------------------- # Simulate exponential data n <- 50 baseline_hazard <- 1 theta <- 0.4 p_censor <- 0.3 prop_tx <- 0.5 tx_indicator = rbinom(n = n, size = 1, prob = prop_tx) event_rate <- baseline_hazard * exp(theta * tx_indicator) censor_rate <- event_rate * p_censor / (1 - p_censor) event_time <- rexp(n = n, rate = event_rate) censor_time <- rexp(n = n, rate = censor_rate) time <- pmin(event_time, censor_time) event <- time == event_time surv_data = data.frame(time = time, event = event, group = tx_indicator) # Fit Cox model and get bootstrap percentile confidence intervals for the exceedance probability # with model-based resampling (see documentation for boot::censboot) cox_fit <- coxph(Surv(time, event) ~ group, data = surv_data) ep <- exceedProbCoxBoot(data = surv_data, cox_fit = cox_fit, j = 1, alpha = 0.05, R = 500) # Plot results with(ep, plot(cutoff, point, type = "l")) with(ep, lines(cutoff, lower, lty = 2)) with(ep, lines(cutoff, upper, lty = 2))library(exceedProb) library(survival) # Cox model ------------------------------------------------------- # Simulate exponential data n <- 50 baseline_hazard <- 1 theta <- 0.4 p_censor <- 0.3 prop_tx <- 0.5 tx_indicator = rbinom(n = n, size = 1, prob = prop_tx) event_rate <- baseline_hazard * exp(theta * tx_indicator) censor_rate <- event_rate * p_censor / (1 - p_censor) event_time <- rexp(n = n, rate = event_rate) censor_time <- rexp(n = n, rate = censor_rate) time <- pmin(event_time, censor_time) event <- time == event_time surv_data = data.frame(time = time, event = event, group = tx_indicator) # Fit Cox model and get bootstrap percentile confidence intervals for the exceedance probability # with model-based resampling (see documentation for boot::censboot) cox_fit <- coxph(Surv(time, event) ~ group, data = surv_data) ep <- exceedProbCoxBoot(data = surv_data, cox_fit = cox_fit, j = 1, alpha = 0.05, R = 500) # Plot results with(ep, plot(cutoff, point, type = "l")) with(ep, lines(cutoff, lower, lty = 2)) with(ep, lines(cutoff, upper, lty = 2))
This function obtains confidence intervals for the non-centrality parameter of a t-distribution.
getDeltaCI(test_stat, alpha, d, n, interval)getDeltaCI(test_stat, alpha, d, n, interval)
test_stat |
Test statistics |
alpha |
Significance level |
d |
Number of parameters in general linear model |
n |
Number of observations in initial study |
interval |
Interval within which to search for roots |
ep Exceedance probability with confidence intervals (vector if cutoff is scalar and matrix otherwise)
This function returns the cdf of a noncentral t-distribution. It is more accurate than stats::pt() for large ncp
pnct(x, df, ncp)pnct(x, df, ncp)
x |
Test statistic |
df |
Degrees of freedom |
ncp |
Noncentrality parameter |
Cumulative probability
This function returns the difference between the lower tail probability of a non-central t-distribution and a confidence level q where the t-distribution has df degrees of freedom and non-centrality parameter delta.
tRoot(delta, test_stat, df, conf_level)tRoot(delta, test_stat, df, conf_level)
delta |
Non-centrality parameter |
test_stat |
Test statistic at which to evaluate the t-distribution |
df |
Degrees of freedom |
conf_level |
Confidence level (usually alpha/2 or 1-alpha/2) |
dif Difference between t-distribution quantile and confidence level