Skip to content

Commit

Permalink
Drop prediction interval for binom
Browse files Browse the repository at this point in the history
  • Loading branch information
mattfidler committed Oct 5, 2023
1 parent 5ec949b commit 0e5f784
Show file tree
Hide file tree
Showing 11 changed files with 135 additions and 116 deletions.
4 changes: 2 additions & 2 deletions R/RcppExports.R
Original file line number Diff line number Diff line change
Expand Up @@ -589,8 +589,8 @@ rxErf <- function(v) {
.Call(`_rxode2_rxErf`, v)
}

binomProbs_ <- function(x, probs, naRm, pred, nIn, mIn) {
.Call(`_rxode2_binomProbs_`, x, probs, naRm, pred, nIn, mIn)
binomProbs_ <- function(x, probs, naRm, nIn, cont) {
.Call(`_rxode2_binomProbs_`, x, probs, naRm, nIn, cont)
}

meanProbs_ <- function(x, probs, naRm, useT, pred, nIn) {
Expand Down
8 changes: 6 additions & 2 deletions R/confint.R
Original file line number Diff line number Diff line change
Expand Up @@ -63,6 +63,10 @@ confint.rxSolve <- function(object, parm = NULL, level = 0.95, ...) {
if (any(names(.args) == "m")) {
.m <- unique(.args$m)
}
.method <- "wald"
if (any(names(.args) == "method")) {
.method <- .args$method
}
.stk <- rxStack(object, parm, doSim=.doSim)
if (!any(names(.stk) == "id") &&
any(names(.stk) == "sim.id")) {
Expand Down Expand Up @@ -112,7 +116,7 @@ confint.rxSolve <- function(object, parm = NULL, level = 0.95, ...) {
} else if (.binom) {
.stk <- .stk[, list(
p1 = .p, eff = rxode2::binomProbs(.SD$value, probs = .p, na.rm = TRUE,
n=.nC, m=.m, pred=.pred),
n=.nC, method=.method),
Percentile = sprintf("%s%%", .p * 100)
),
by = c("time", "trt", .by)
Expand Down Expand Up @@ -154,7 +158,7 @@ confint.rxSolve <- function(object, parm = NULL, level = 0.95, ...) {
print(.nC)
.ret <- .ret[, list(p1 = .p,
eff = rxode2::binomProbs(.SD$value, probs = .p, na.rm = TRUE,
n=.nC, m=.m, pred=.pred)),
n=.nC, method=.method)),
by = c("id", "time", "trt", .by)]
} else {
.ret <- .ret[, list(p1 = .p,
Expand Down
54 changes: 34 additions & 20 deletions R/utils.R
Original file line number Diff line number Diff line change
Expand Up @@ -879,14 +879,6 @@ meanProbs.default <- function(x, probs=seq(0, 1, 0.25), na.rm=FALSE,

#' Calculate expected confidence bands with binomial sampling distribution
#'
#' The generic function `binomProbs` produces expected confidence bands
#' using the Wilson score interval with continuity correction
#' (described in Newcombe 1998).
#'
#' The expected prediction intervals use the Frequentist Nelson's
#' predictional interval (described by Lu 2020).
#'
#'
#' This is meant to perform in the same way as `quantile()` so it can
#' be a drop in replacement for code using `quantile()` but using
#' distributional assumptions.
Expand All @@ -896,33 +888,55 @@ meanProbs.default <- function(x, probs=seq(0, 1, 0.25), na.rm=FALSE,
#'
#' @param x numeric vector whose mean and probability based confidence
#' values are wanted, NA and NaN values are not allowed in numeric
#' vectors unless ‘na.rm’ is ‘TRUE’.
#' vectors unless `na.rm` is `TRUE`.
#'
#' @param probs numeric vector of probabilities with values in the
#' interval 0 to 1, inclusive. When 0, it represents the maximum
#' observed, when 1, it represents the maximum observed. When 0.5 it
#' represents the expected probability (mean)
#' represents the expected probability (mean).
#'
#' @param na.rm logical; if true, any NA and NaN's are removed from
#' `x` before the quantiles are computed.
#'
#' @param names logical; if true, the result has a names attribute.
#'
#' @param onlyProbs logical; if true, only return the probability
#' based confidence interval/prediction interval estimates,
#' otherwise return extra statistics
#' otherwise return extra statistics.
#'
#' @param n integer/integerish; this is the n used to calculate the
#' prediction or confidence interval. When `n=0` (default) use the
#' number of non-`NA` observations. When calculating the prediction
#' interval, this represents the number of observations used in the
#' input ("true") distribution.
#' @param m integer/integerish; this is the m used to calculate the
#' prediction interval. It represents the expected number of
#' samples in the new statistic. When `1=0` (default) use the
#' number of non-`NA` observations in the input `x`.
#'
#' @param method gives the method for calculating the confidence
#' interval.
#'
#' Can be:
#'
#' - "argestiCoull" -- Agresti-Coull method. For a 95\% confidence
#' interval, this method does not use the concept of "adding 2
#' successes and 2 failures," but rather uses the formulas explicitly
#' described in the following link:
#'
#' https://en.wikipedia.org/wiki/Binomial_proportion_confidence_interval#Agresti-Coull_Interval.
#'
#' - "wilson" -- Wilson Method
#'
#' - "wilsonCorrect" -- Wilson method with continuity correction
#'
#' - "wald" -- Wald confidence interval or standard z approximation.
#'
#' @param ... Arguments passed to default method, allows many
#' different methods to be applied.
#'
#' @return By default the return has the probabilities as names (if
#' named) with the points where the expected distribution are
#' located given the sampling mean and standard deviation. If
#' `onlyProbs=FALSE` then it would prepend mean, variance, standard
#' deviation, minimum, maximum and number of non-NA observations.
#'
#' @export
#' @author Matthew L. Fidler
#' @references
Expand Down Expand Up @@ -962,19 +976,19 @@ binomProbs <- function(x, ...) {
#' @rdname binomProbs
#' @export
binomProbs.default <- function(x, probs=c(0.025, 0.05, 0.5, 0.95, 0.975), na.rm=FALSE,
names=TRUE, onlyProbs=TRUE, pred=FALSE, n=0L, m=0L) {
names=TRUE, onlyProbs=TRUE, n=0L,
method=c("wilson", "wilsonCorrect", "agrestiCoull", "wald")) {
method <- match.arg(method)
method <- setNames(c("wilson"=1L, "wilsonCorrect"=0L, "agrestiCoull"=3, "wald"=2L)[method], NULL)
checkmate::assertNumeric(x, min.len=1, lower=0.0, upper=1.0)
x <- as.double(x)
checkmate::assertIntegerish(n, min.len=1, lower=0, any.missing=FALSE)
checkmate::assertIntegerish(m, min.len=1, lower=0, any.missing=FALSE)
n <- as.integer(n)
m <- as.integer(m)
checkmate::assertNumeric(probs, min.len=1, any.missing = FALSE, lower=0.0, upper=1.0)
checkmate::assertLogical(na.rm, any.missing=FALSE, len=1)
checkmate::assertLogical(names, any.missing=FALSE, len=1)
checkmate::assertLogical(onlyProbs, any.missing=FALSE, len=1)
checkmate::assertLogical(pred, any.missing=FALSE, len=1)
.ret <- .Call(`_rxode2_binomProbs_`, x, probs, na.rm, pred, n, m)
.ret <- .Call(`_rxode2_binomProbs_`, x, probs, na.rm, n, method)
.names <- NULL
if (names) {
.names <- paste0(probs*100, "%")
Expand Down
8 changes: 4 additions & 4 deletions inst/include/rxode2_RcppExports.h
Original file line number Diff line number Diff line change
Expand Up @@ -1073,17 +1073,17 @@ namespace rxode2 {
return Rcpp::as<NumericVector >(rcpp_result_gen);
}

inline NumericVector binomProbs_(NumericVector x, NumericVector probs, bool naRm, bool pred, int nIn, int mIn) {
typedef SEXP(*Ptr_binomProbs_)(SEXP,SEXP,SEXP,SEXP,SEXP,SEXP);
inline NumericVector binomProbs_(NumericVector x, NumericVector probs, bool naRm, int nIn, int cont) {
typedef SEXP(*Ptr_binomProbs_)(SEXP,SEXP,SEXP,SEXP,SEXP);
static Ptr_binomProbs_ p_binomProbs_ = NULL;
if (p_binomProbs_ == NULL) {
validateSignature("NumericVector(*binomProbs_)(NumericVector,NumericVector,bool,bool,int,int)");
validateSignature("NumericVector(*binomProbs_)(NumericVector,NumericVector,bool,int,int)");
p_binomProbs_ = (Ptr_binomProbs_)R_GetCCallable("rxode2", "_rxode2_binomProbs_");
}
RObject rcpp_result_gen;
{
RNGScope RCPP_rngScope_gen;
rcpp_result_gen = p_binomProbs_(Shield<SEXP>(Rcpp::wrap(x)), Shield<SEXP>(Rcpp::wrap(probs)), Shield<SEXP>(Rcpp::wrap(naRm)), Shield<SEXP>(Rcpp::wrap(pred)), Shield<SEXP>(Rcpp::wrap(nIn)), Shield<SEXP>(Rcpp::wrap(mIn)));
rcpp_result_gen = p_binomProbs_(Shield<SEXP>(Rcpp::wrap(x)), Shield<SEXP>(Rcpp::wrap(probs)), Shield<SEXP>(Rcpp::wrap(naRm)), Shield<SEXP>(Rcpp::wrap(nIn)), Shield<SEXP>(Rcpp::wrap(cont)));
}
if (rcpp_result_gen.inherits("interrupted-error"))
throw Rcpp::internal::InterruptedException();
Expand Down
41 changes: 23 additions & 18 deletions man/binomProbs.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

2 changes: 1 addition & 1 deletion man/meanProbs.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

4 changes: 2 additions & 2 deletions man/rxode2.Rd

Some generated files are not rendered by default. Learn more about how customized files appear on GitHub.

15 changes: 7 additions & 8 deletions src/RcppExports.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -1825,25 +1825,24 @@ RcppExport SEXP _rxode2_rxErf(SEXP vSEXP) {
return rcpp_result_gen;
}
// binomProbs_
NumericVector binomProbs_(NumericVector x, NumericVector probs, bool naRm, bool pred, int nIn, int mIn);
static SEXP _rxode2_binomProbs__try(SEXP xSEXP, SEXP probsSEXP, SEXP naRmSEXP, SEXP predSEXP, SEXP nInSEXP, SEXP mInSEXP) {
NumericVector binomProbs_(NumericVector x, NumericVector probs, bool naRm, int nIn, int cont);
static SEXP _rxode2_binomProbs__try(SEXP xSEXP, SEXP probsSEXP, SEXP naRmSEXP, SEXP nInSEXP, SEXP contSEXP) {
BEGIN_RCPP
Rcpp::RObject rcpp_result_gen;
Rcpp::traits::input_parameter< NumericVector >::type x(xSEXP);
Rcpp::traits::input_parameter< NumericVector >::type probs(probsSEXP);
Rcpp::traits::input_parameter< bool >::type naRm(naRmSEXP);
Rcpp::traits::input_parameter< bool >::type pred(predSEXP);
Rcpp::traits::input_parameter< int >::type nIn(nInSEXP);
Rcpp::traits::input_parameter< int >::type mIn(mInSEXP);
rcpp_result_gen = Rcpp::wrap(binomProbs_(x, probs, naRm, pred, nIn, mIn));
Rcpp::traits::input_parameter< int >::type cont(contSEXP);
rcpp_result_gen = Rcpp::wrap(binomProbs_(x, probs, naRm, nIn, cont));
return rcpp_result_gen;
END_RCPP_RETURN_ERROR
}
RcppExport SEXP _rxode2_binomProbs_(SEXP xSEXP, SEXP probsSEXP, SEXP naRmSEXP, SEXP predSEXP, SEXP nInSEXP, SEXP mInSEXP) {
RcppExport SEXP _rxode2_binomProbs_(SEXP xSEXP, SEXP probsSEXP, SEXP naRmSEXP, SEXP nInSEXP, SEXP contSEXP) {
SEXP rcpp_result_gen;
{
Rcpp::RNGScope rcpp_rngScope_gen;
rcpp_result_gen = PROTECT(_rxode2_binomProbs__try(xSEXP, probsSEXP, naRmSEXP, predSEXP, nInSEXP, mInSEXP));
rcpp_result_gen = PROTECT(_rxode2_binomProbs__try(xSEXP, probsSEXP, naRmSEXP, nInSEXP, contSEXP));
}
Rboolean rcpp_isInterrupt_gen = Rf_inherits(rcpp_result_gen, "interrupted-error");
if (rcpp_isInterrupt_gen) {
Expand Down Expand Up @@ -1957,7 +1956,7 @@ static int _rxode2_RcppExport_validate(const char* sig) {
signatures.insert("RObject(*rxSymInvCholEnvCalculate)(List,std::string,Nullable<NumericVector>)");
signatures.insert("LogicalVector(*isNullZero)(RObject)");
signatures.insert("NumericVector(*rxErf)(NumericVector)");
signatures.insert("NumericVector(*binomProbs_)(NumericVector,NumericVector,bool,bool,int,int)");
signatures.insert("NumericVector(*binomProbs_)(NumericVector,NumericVector,bool,int,int)");
signatures.insert("NumericVector(*meanProbs_)(NumericVector,NumericVector,bool,bool,bool,int)");
}
return signatures.find(sig) != signatures.end();
Expand Down
4 changes: 2 additions & 2 deletions src/init.c
Original file line number Diff line number Diff line change
Expand Up @@ -283,7 +283,7 @@ SEXP _rxode2_trans(SEXP parse_file, SEXP prefix, SEXP model_md5, SEXP parseStr,
SEXP _rxode2_assignSeedInfo(void);
SEXP _rxSetSeed(SEXP);
SEXP _rxode2_meanProbs_(SEXP x, SEXP probs, SEXP naRm, SEXP useT, SEXP pred, SEXP inN);
SEXP _rxode2_binomProbs_(SEXP x, SEXP probs, SEXP naRm, SEXP pred, SEXP inN, SEXP inM);
SEXP _rxode2_binomProbs_(SEXP x, SEXP probs, SEXP naRm, SEXP inN, SEXP cont);

typedef SEXP (*lotriMat_type) (SEXP, SEXP, SEXP);
typedef SEXP (*asLotriMat_type) (SEXP, SEXP, SEXP);
Expand Down Expand Up @@ -361,7 +361,7 @@ extern SEXP chin(SEXP x, SEXP table);
SEXP _rxode2_RcppExport_registerCCallable(void);
void R_init_rxode2(DllInfo *info){
R_CallMethodDef callMethods[] = {
{"_rxode2_binomProbs_", (DL_FUNC) &_rxode2_binomProbs_, 6},
{"_rxode2_binomProbs_", (DL_FUNC) &_rxode2_binomProbs_, 5},
{"_rxode2_meanProbs_", (DL_FUNC) &_rxode2_meanProbs_, 6},
{"_rxode2_getEtRxsolve", (DL_FUNC) &_rxode2_getEtRxsolve, 1},
{"_rxode2_assignSeedInfo", (DL_FUNC) &_rxode2_assignSeedInfo, 0},
Expand Down
Loading

0 comments on commit 0e5f784

Please sign in to comment.