| Title: | Joint Confidence Regions |
|---|---|
| Description: | Computing and plotting joint confidence regions and intervals. Regions include classical ellipsoids, minimum-volume or minimum-length regions, and an empirical Bayes region. Intervals include the TOST procedure with ordinary or expanded intervals and a fixed-sequence procedure. Such regions and intervals are useful e.g., for the assessment of multi-parameter (bio-)equivalence. Joint confidence regions for the mean and variance of a normal distribution are available as well. |
| Authors: | Philip Pallmann [aut, cre] |
| Maintainer: | Philip Pallmann <[email protected]> |
| License: | GPL-2 |
| Version: | 0.3.3 |
| Built: | 2025-03-04 03:46:38 UTC |
| Source: | https://github.com/philippallmann/jocre |
This package provides functions for computing and plotting joint confidence regions as well as (simultaneous) confidence intervals, with a focus on multivariate normal parameter vectors and parameters of the normal distribution i.e., mean and variance.
| Package: | jocre |
| Type: | Package |
| Version: | 0.3.3 |
| Date: | 2017-05-12 |
| License: | GPL-2 |
cset computes joint confidence regions and (simultaneous) confidence intervals around multivariate normal means. The parameter estimates and interval bounds (in the case of confidence regions their boundaries are projected onto the axes) are displayed with print or summary, and plot generates a graph of the estimate and region or intervals and allows for some fine tuning e.g., shading of an equivalence region.
csetMV computes joint confidence regions for the mean and variance (or standard deviation) of a normal distribution. plot, print, and summary produce the corresponding graphics and summaries.
Insights into the conservativeness of the two one-sided test procedure can be obtained with iutsize.
Confidence levels of a joint ellipse and marginal intervals can be "translated" back and forth with translate.
Philip Pallmann ([email protected])
Philip Pallmann (2017) Joint confidence regions with the R package jocre. In preparation.
## Not run: data("marzo") plot(cset(dat=marzo, method="limacon")) data("wires") plot(csetMV(dat=wires, method="mood")) ## End(Not run)## Not run: data("marzo") plot(cset(dat=marzo, method="limacon")) data("wires") plot(csetMV(dat=wires, method="mood")) ## End(Not run)
Computes boundaries of (simultaneous) confidence regions and intervals around multivariate normal means using different methods.
cset(dat, method, alpha=0.1, steps=NULL, nboot=1e4, TsengBrownA=1, TsengBrownB=1)cset(dat, method, alpha=0.1, steps=NULL, nboot=1e4, TsengBrownA=1, TsengBrownB=1)
dat |
A matrix or data.frame with independent units in rows and multivariate outcomes in columns. |
method |
A character string specifying the method to be used. See details for available |
alpha |
A numeric value giving the type I error level to be controlled. Default is |
steps |
An integer setting the initial number of steps for the search algorithm. Default is |
nboot |
A numeric giving the number of bootstrap replications to be used when |
TsengBrownA |
A numeric giving the parameter A to be used when |
TsengBrownB |
A numeric giving the parameter B to be used when |
Available methods for confidence regions are: boot.kern for the nonparametric bootstrap method using kernel density estimation described in Pallmann & Jaki (2017); emp.bayes for the empirical Bayes region described in Casella & Hwang (1983); hotelling for the Hotelling-type region described in Wang et al (1999); limacon.asy for the limacon-shaped mimimum expected volume region described in Brown et al (1995); limacon.fin for the finite-sample variant of the minimum expected volume region described in Berger & Hsu (1996); standard.cor for the standard region incorporating correlation between parameters described in Chew (1966); standard.ind for the standard region ignoring correlation between parameters; tost for the two one-sided test (TOST) intervals described in Schuirmann (1987); tseng for the mimimum expected interval length region described in Tseng (2002); tseng.brown for the pseudo-empirical Bayes region described in Tseng & Brown (1997).
Available methods for confidence intervals are: expanded for the two one-sided test (TOST) procedure (Schuirmann 1987) using the expanded intervals described e.g., in Bofinger (1992) and Hsu et al. (1994); fix.seq for the fixed sequence intervals described in Maurer et al (1995) and Hsu & Berger (1999); tost for the two one-sided test (TOST) intervals described in Schuirmann (1987).
See also an overview and comparison of all methods in Pallmann & Jaki (2017).
An object of class JOC.
Warning: please use with care! Some of the functionality has not yet been thoroughly tested.
Philip Pallmann ([email protected])
Roger L. Berger & Jason C. Hsu (1996) Bioequivalence trials, intersection-union tests and equivalence confidence sets. Statistical Science, 11(4), 283–319.
Eve Bofinger (1992) Expanded confidence intervals, one-sided tests, and equivalence testing. Journal of Biopharmaceutical Statistics, 2(2), 181–188.
Lawrence D. Brown, George Casella, J. T. Gene Hwang (1995) Optimal confidence sets, bioequivalence, and the limacon of Pascal. Journal of the American Statistical Association, 90(431), 880–889.
George Casella & Jiunn T. Hwang (1983) Empirical Bayes confidence sets for the mean of a multivariate normal distribution. Journal of the American Statistical Association, 78(383), 688–698.
Victor Chew (1966) Confidence, prediction, and tolerance regions for the multivariate normal distribution. Journal of the American Statistical Association, 61(315), 605–617.
Jason C. Hsu & Roger L. Berger (1999) Stepwise confidence intervals without multiplicity adjustment for dose-response and toxicity studies. Journal of the American Statistical Association, 94(446), 468–482.
Jason C. Hsu, J. T. Gene Hwang, Hung-Kung Liu, Stephen J. Ruberg (1994) Confidence intervals associated with tests for bioequivalence. Biometrika, 81(1), 103–114.
Willi Maurer, Ludwig A. Hothorn, Walter Lehmacher (1995) Multiple comparisons in drug clinical trials and preclinical assays: a priori ordered hypotheses. In: Joachim Vollmar (editor), Biometrie in der Chemisch-Pharmazeutischen Industrie, vol. 6, pp. 3–18. Fischer-Verlag, Stuttgart, Germany.
Philip Pallmann & Thomas Jaki (2017) Simultaneous confidence regions and intervals for multivariate bioequivalence. Submitted to Statistics in Medicine.
Donald J. Schuirmann (1987) A comparison of the two one-sided tests procedure and the power approach for assessing the equivalence of average bioavailability. Journal of Pharmacokinetics and Biopharmaceutics, 15(6), 657–680.
Yu-Ling Tseng (2002) Optimal confidence sets for testing average bioequivalence. Test, 11(1), 127–141.
Yu-Ling Tseng & Lawrence D. Brown (1997) Good exact confidence sets for a multivariate normal mean. The Annals of Statistics, 25(5), 2228–2258.
Weizhen Wang, J. T. Gene Hwang, Anirban DasGupta (1999) Statistical tests for multivariate bioequivalence. Biometrika, 86(2), 395–402.
# bootkern not included so far
csetMV for (simultaneous) confidence regions for normal mean and variance.
## Not run: # Example 1: simultaneous 90% confidence intervals for trivariate data trivar <- mvtnorm::rmvnorm(n=20, mean=rep(0.05, 3), sigma=toeplitz(c(0.05, 0.04, 0.03))) colnames(trivar) <- c("AUCinf", "AUCt", "Cmax") tost <- cset(dat=trivar, method="tost", alpha=0.1) summary(tost) # Example 2: simultaneous 90% confidence regions for bivariate data bivar <- mvtnorm::rmvnorm(n=20, mean=rep(0.05, 2), sigma=toeplitz(c(0.05, 0.04))) colnames(bivar) <- c("AUC", "Cmax") hotelling <- cset(dat=bivar, method="hotelling", alpha=0.1) summary(hotelling) plot(hotelling, main="90% Hotelling Region") limacon <- cset(dat=bivar, method="limacon.asy", alpha=0.1) summary(limacon) plot(limacon, main="90% Limacon Region") tseng <- cset(dat=bivar, method="tseng", alpha=0.1) summary(tseng) plot(tseng, main="90% Tseng Region") ## End(Not run)## Not run: # Example 1: simultaneous 90% confidence intervals for trivariate data trivar <- mvtnorm::rmvnorm(n=20, mean=rep(0.05, 3), sigma=toeplitz(c(0.05, 0.04, 0.03))) colnames(trivar) <- c("AUCinf", "AUCt", "Cmax") tost <- cset(dat=trivar, method="tost", alpha=0.1) summary(tost) # Example 2: simultaneous 90% confidence regions for bivariate data bivar <- mvtnorm::rmvnorm(n=20, mean=rep(0.05, 2), sigma=toeplitz(c(0.05, 0.04))) colnames(bivar) <- c("AUC", "Cmax") hotelling <- cset(dat=bivar, method="hotelling", alpha=0.1) summary(hotelling) plot(hotelling, main="90% Hotelling Region") limacon <- cset(dat=bivar, method="limacon.asy", alpha=0.1) summary(limacon) plot(limacon, main="90% Limacon Region") tseng <- cset(dat=bivar, method="tseng", alpha=0.1) summary(tseng) plot(tseng, main="90% Tseng Region") ## End(Not run)
Computes boundaries of (simultaneous) confidence regions for the mean and variance of a normal distribution using different methods.
csetMV(dat, n, method, alpha=0.1, scale="var", steps=500)csetMV(dat, n, method, alpha=0.1, scale="var", steps=500)
dat |
A vector of numeric values assumed to follow a normal distribution. Not required for |
n |
A numeric value giving the sample size. Only required for |
method |
A character string specifying the method to be used. See details for available |
alpha |
A numeric value giving the type I error level to be controlled. Default is |
scale |
A character string specifying whether the variance ( |
steps |
An integer setting the initial number of steps for the search algorithm. Default is |
Available methods are: mood for the classical region described in Mood (1950); large for the large-sample approximation region described in section 4.1 of Arnold & Shavelle (1998); plugin for a plug-in variant of the large-sample approximation region described in section 4.2 of Arnold & Shavelle (1998); pluginF for the plug-in variant of the large-sample approximation region described in section 4.3 of Arnold & Shavelle (1998) using an asymptotic F distribution as in Douglas (1993); lrt for the likelihood ratio test region described in section 4.4 of Arnold & Shavelle (1998); cheng.iles for the region described in Cheng & Iles (1983); min.area for the minimum-area region described in Frey et al. (2009).
An object of class JOCMV.
Warning: please use with care! Some of the functionality has not yet been thoroughly tested.
Philip Pallmann ([email protected])
Barry C. Arnold & Robert M. Shaville (1998) Joint confidence sets for the mean and variance of a normal distribution. The American Statistician, 52(2), 133–140.
R. C. H. Cheng & T. C. Iles (1983) Confidence bands for cumulative distribution functions of continuous random variables. Technometrics, 25(1), 77–86.
J. B. Douglas (1993) Confidence regions for parameter pairs. The American Statistician, 47(1), 43–45.
Jesse Frey, Osvaldo Marrero, Douglas Norton (2009) Minimum-area confidence sets for a normal distribution. Journal of Statistical Planning and Inference, 139(3), 1023–1032.
Alexander M. Mood (1950) Introduction to the Theory of Statistics. McGraw-Hill, New York, NY.
cset for (simultaneous) confidence regions and intervals around multivariate normal means.
## Not run: # Simultaneous 90% confidence regions for the mean and variance or sd of univariate normal data univar <- rnorm(n=50) moodvar <- csetMV(dat=univar, method="mood", alpha=0.1, scale="var") summary(moodvar) plot(moodvar) moodsd <- csetMV(dat=univar, method="mood", alpha=0.1, scale="sd") summary(moodsd) plot(moodsd) ## End(Not run)## Not run: # Simultaneous 90% confidence regions for the mean and variance or sd of univariate normal data univar <- rnorm(n=50) moodvar <- csetMV(dat=univar, method="mood", alpha=0.1, scale="var") summary(moodvar) plot(moodvar) moodsd <- csetMV(dat=univar, method="mood", alpha=0.1, scale="sd") summary(moodsd) plot(moodsd) ## End(Not run)
JOC and JOCMV
Generic functions for summarising and plotting objects of class JOC or JOCMV.
## S3 method for class 'JOC' plot(x, equi=log(c(0.8, 1.25)), axnames=NULL, main=NULL, xlim=log(c(0.77, 1.3)), ylim=log(c(0.77, 1.3)), col="black", convexify=FALSE, ...) ## S3 method for class 'JOC' print(x, digits=max(3, getOption("digits") - 4), ...) ## S3 method for class 'JOC' summary(object, digits=max(3, getOption("digits") - 4), ...) ## S3 method for class 'JOCMV' plot(x, axnames=NULL, main=NULL, xlim=NULL, ylim=NULL, col="black", ...) ## S3 method for class 'JOCMV' print(x, digits=max(3, getOption("digits") - 4), ...) ## S3 method for class 'JOCMV' summary(object, digits=max(3, getOption("digits") - 4), ...)## S3 method for class 'JOC' plot(x, equi=log(c(0.8, 1.25)), axnames=NULL, main=NULL, xlim=log(c(0.77, 1.3)), ylim=log(c(0.77, 1.3)), col="black", convexify=FALSE, ...) ## S3 method for class 'JOC' print(x, digits=max(3, getOption("digits") - 4), ...) ## S3 method for class 'JOC' summary(object, digits=max(3, getOption("digits") - 4), ...) ## S3 method for class 'JOCMV' plot(x, axnames=NULL, main=NULL, xlim=NULL, ylim=NULL, col="black", ...) ## S3 method for class 'JOCMV' print(x, digits=max(3, getOption("digits") - 4), ...) ## S3 method for class 'JOCMV' summary(object, digits=max(3, getOption("digits") - 4), ...)
x |
An output object of class |
object |
An output object of class |
digits |
A numeric value giving the number of significant digits to be printed. |
equi |
A numeric vector of length 2 specifying the equivalence region (lower and upper equivalence threshold) to be shaded in grey. When set to |
axnames |
A vector of two character strings giving the x and y axis labels. For |
main |
A character string giving the plot title. Default is |
xlim |
A numeric vector of length two specifying the plotting range on the x-axis. Default is |
ylim |
A numeric vector of length two specifying the plotting range on the y-axis. Default is |
col |
A character string specifying the colour of the plotted region or intervals. |
convexify |
A logical specifying whether the convex hull around a non-convex region should be plotted instead of the region itself. Ignored unless |
... |
Further plotting arguments to be passed to methods. Type |
print and summary summarise the estimates and confidence set boundaries of an object of class JOC or JOCMV that was created with cset or csetMV, respectively. plot displays a (simultaneous) confidence region or intervals when applied to an object of class JOC or JOCMV.
An on-screen summary or graphical display.
Warning: please use with care! Some of the functionality has not yet been thoroughly tested.
Philip Pallmann ([email protected])
Philip Pallmann & Thomas Jaki (2017) Simultaneous confidence regions and intervals for multivariate bioequivalence. Submitted to Statistics in Medicine.
cset and csetMV for computing (simultaneous) confidence regions and intervals.
## Not run: # Example 1: simultaneous 90% confidence region for bivariate data bivar <- mvtnorm::rmvnorm(n=100, mean=rep(0.05, 2), sigma=diag(2) * 0.05) hotelling <- cset(dat=bivar, method="hotelling", alpha=0.1) summary(hotelling) plot(hotelling, main="90% Hotelling Region") # Example 2: simultaneous 90% confidence region for the mean and variance of univariate normal data univar <- rnorm(n=50) moodvar <- csetMV(dat=univar, method="mood", alpha=0.1, scale="var") summary(moodvar) plot(moodvar, main="90% Mood Region") ## End(Not run)## Not run: # Example 1: simultaneous 90% confidence region for bivariate data bivar <- mvtnorm::rmvnorm(n=100, mean=rep(0.05, 2), sigma=diag(2) * 0.05) hotelling <- cset(dat=bivar, method="hotelling", alpha=0.1) summary(hotelling) plot(hotelling, main="90% Hotelling Region") # Example 2: simultaneous 90% confidence region for the mean and variance of univariate normal data univar <- rnorm(n=50) moodvar <- csetMV(dat=univar, method="mood", alpha=0.1, scale="var") summary(moodvar) plot(moodvar, main="90% Mood Region") ## End(Not run)
Computes the actual size of a intersection union test procedure that corresponds to a (1 – alpha) confidence set.
iutsize(p, n, alpha=0.1, sim=1e6)iutsize(p, n, alpha=0.1, sim=1e6)
p |
An integer giving the number of dimensions. |
n |
An integer giving the sample size. |
alpha |
A numeric value specifying the type I error level to be controlled. Default is 0.1. |
sim |
An integer giving the number of simulations to be carried out. Default is 1 million. |
A (1 – alpha) confidence set can be used to derive a two one-sided tests (TOST) procedure (Schuirmann 1987) whereby type I error rate control is ensured at level alpha due to the intersection union principle (Berger 1982). The actual test size, however, is often substantially lower than alpha i.e., the approach is conservative. It is well known for the one-dimensional case that the TOST corresponding to a (1 – alpha) confidence interval has size (1 – alpha/2). This function computes the achieved test size with dimension p and n according to the formula on p. 399 of Wang et al (1999).
A numeric value giving the actual size of the test.
Philip Pallmann ([email protected])
Roger L. Berger (1982) Multiparameter hypothesis testing and acceptance sampling. Technometrics, 24(4), 295–300.
Donald J. Schuirmann (1987) A comparison of the two one-sided tests procedure and the power approach for assessing the equivalence of average bioavailability. Journal of Pharmacokinetics and Biopharmaceutics, 15(6), 657–680.
Weizhen Wang, J. T. Gene Hwang, Anirban DasGupta (1999) Statistical tests for multivariate bioequivalence. Biometrika, 86(2), 395–402.
# For p=1 we get the well-known result that the 90% CI corresponds to the TOST at 5%: #iutsize(p=1, n=20) # With increasing dimension the test gets conservative: #iutsize(p=2, n=20) #iutsize(p=3, n=20) # For p>1 the conservativeness also depends on sample size: #iutsize(p=2, n=10) #iutsize(p=2, n=1000)# For p=1 we get the well-known result that the 90% CI corresponds to the TOST at 5%: #iutsize(p=1, n=20) # With increasing dimension the test gets conservative: #iutsize(p=2, n=20) #iutsize(p=3, n=20) # For p>1 the conservativeness also depends on sample size: #iutsize(p=2, n=10) #iutsize(p=2, n=1000)
Pharmacokinetic data from a study on the bioequivalence of a test and a reference formulation of ticlopidine hydrochloride in 24 healthy male volunteers, using a randomised crossover design (Marzo et al. 2002).
data("marzo")data("marzo")
A data frame with 24 observations on the following 8 variables.
VolunteerA numeric vector giving the volunteer ID.
SequenceA factor with levels RT and TR specifying the sequence a volunteer was randomised to (R=reference, T=test).
Cmax_TA numeric vector of the maximum concentration (Cmax) with the test product.
Cmax_RA numeric vector of the maximum concentration (Cmax) with the reference product.
AUC_TA numeric vector of the area under the concentration-time curve (AUC) from zero to the last observed time point with the test product.
AUC_RA numeric vector of the area under the concentration-time curve (AUC) from zero to the last observed time point with the reference product.
AUCinf_TA numeric vector of the area under the concentration-time curve AUC) from zero to infinity with the test product.
AUCinf_RA numeric vector of the area under the concentration-time curve AUC) from zero to infinity with the reference product.
The pharmacokinetic parameters (Cmax and AUC) were calculated using a non-compartmental approach. The data were taken from Tables I and II of Marzo et al. (2002).
Antonio Marzo, Lorenzo Dal Bo, Antonio Rusca, Pierangelo Zini (2002) Bioequivalence of ticlopidine hydrochloride administered in single dose to healthy volunteers. Pharmacological Research, 46(5), 401–407.
Philip Pallmann & Thomas Jaki (2017) Simultaneous confidence regions and intervals for multivariate bioequivalence. Submitted to Statistics in Medicine.
data(marzo) ## An example analysis of Cmax assuming log-normality # Difference of log(Cmax) marzo$deltalogCmax <- log(marzo$Cmax_T) - log(marzo$Cmax_R) # Estimated mean treatment effect with SE mean(marzo$deltalogCmax) sd(marzo$deltalogCmax) / sqrt(nrow(marzo)) # Two one-sided test (TOST) p-values t.test(x=marzo$deltalogCmax, alternative="less", mu=log(1.25))$p.value t.test(x=marzo$deltalogCmax, alternative="greater", mu=log(0.80))$p.value # 90% TOST confidence interval t.test(x=marzo$deltalogCmax, conf.level=0.9)$conf.int[1:2]data(marzo) ## An example analysis of Cmax assuming log-normality # Difference of log(Cmax) marzo$deltalogCmax <- log(marzo$Cmax_T) - log(marzo$Cmax_R) # Estimated mean treatment effect with SE mean(marzo$deltalogCmax) sd(marzo$deltalogCmax) / sqrt(nrow(marzo)) # Two one-sided test (TOST) p-values t.test(x=marzo$deltalogCmax, alternative="less", mu=log(1.25))$p.value t.test(x=marzo$deltalogCmax, alternative="greater", mu=log(0.80))$p.value # 90% TOST confidence interval t.test(x=marzo$deltalogCmax, conf.level=0.9)$conf.int[1:2]
Translates the confidence level of a joint 100(1 – alpha)% confidence ellipse into that of the corresponding marginal confidence interval when projecting the ellipse's boundary onto the axes. Also does the "inverse operation" i.e., calculates the confidence level of a joint confidence ellipse so that its perpendicular shadows onto the axes are 100(1 – alpha)% confidence intervals.
translate(level=0.95, ddf, direction)translate(level=0.95, ddf, direction)
level |
A numeric value giving the confidence level. |
ddf |
An integer specifying the denominator degrees of freedom. Setting this to |
direction |
A character string indicating what is to be computed. Choose either |
Setting direction="ci2cr" calculates the confidence level of a confidence interval generating ellipse (CIGE) whose perpendicular shadows onto the axes are 100(1 – alpha)% confidence intervals with a marginal confidence level (1 – alpha) as specified in level; see p. 205 of Fox (2008).
On the other hand, setting direction="cr2ci" computes the marginal confidence level of the intervals obtained by projecting a joint 100(1 – alpha)% confidence ellipse with (1 – alpha) as specified in level; see p. 254 of Monette (1990). These marginal intervals can be viewed as including a Scheffe penalty (Scheffe 1953).
For ddf=0 the F-distribution used for calculating the confidence levels is replaced with an asymptotic chi-square distribution.
A numeric value giving the calculated confidence level.
Philip Pallmann ([email protected])
John Fox (2008) Applied Linear Regression and Generalized Linear Models. Second Edition. SAGE, Thousand Oaks, CA.
Georges Monette (1990). Geometry of multiple regression and interactive 3-D graphics. In: John Fox & J. Scott Long (eds.) Modern Methods of Data Analysis. SAGE, Newbury Park, CA.
Henry Scheffe (1953) A method for judging all contrasts in the analysis of variance. Biometrika, 40(1–2), 87–104.
# Get CIGE level translate(0.95, ddf=1, "ci2cr") translate(0.95, ddf=9999, "ci2cr") translate(0.95, ddf=0, "ci2cr") # Get Scheffe CI level translate(0.95, ddf=1, "cr2ci") translate(0.95, ddf=9999, "cr2ci") translate(0.95, ddf=0, "cr2ci")# Get CIGE level translate(0.95, ddf=1, "ci2cr") translate(0.95, ddf=9999, "ci2cr") translate(0.95, ddf=0, "ci2cr") # Get Scheffe CI level translate(0.95, ddf=1, "cr2ci") translate(0.95, ddf=9999, "cr2ci") translate(0.95, ddf=0, "cr2ci")
Data from a study in quality control assessing the breaking strengths of 20 wire connections between a semiconductor wafer and a terminal post (King 1971).
data("wires")data("wires")
A data frame with 20 observations on the following 2 variables.
StrengthA numeric vector giving the strength at which the connection failed.
FailureA factor with levels b and w specifying whether the bond or wire failed.
The data were taken from Table 4.1 of Nelson (1982).
James R. King (1971) Probability Charts for Decision Making. Industrial Press, New York, NY.
Wayne B. Nelson (1982) Applied Life Data Analysis. Wiley, Hoboken, NJ.
## Not run: data(wires) # Simultaneous 90% confidence regions for the mean and variance plot(csetMV(wires$Strength, method="mood"), main="Mood") plot(csetMV(wires$Strength, method="lrt"), main="LRT") ## End(Not run)## Not run: data(wires) # Simultaneous 90% confidence regions for the mean and variance plot(csetMV(wires$Strength, method="mood"), main="Mood") plot(csetMV(wires$Strength, method="lrt"), main="LRT") ## End(Not run)