| Title: | FAB Confidence Intervals |
|---|---|
| Description: | Frequentist assisted by Bayes (FAB) confidence interval construction. See 'Adaptive multigroup confidence intervals with constant coverage' by Yu and Hoff <DOI:10.1093/biomet/asy009> and 'Exact adaptive confidence intervals for linear regression coefficients' by Hoff and Yu <DOI:10.1214/18-EJS1517>. |
| Authors: | Peter Hoff and Chaoyu Yu |
| Maintainer: | Peter Hoff <[email protected]> |
| License: | GPL-3 |
| Version: | 0.2 |
| Built: | 2025-02-02 02:39:03 UTC |
| Source: | https://github.com/pdhoff/fabci |
Compute emprirical Bayes estimates of the error variance and distribution of the regression coefficients.
ebayes_est(y, X, emu = FALSE, dof = min(50, round(0.5 * (dim(X)[1] - dim(X)[2]))))ebayes_est(y, X, emu = FALSE, dof = min(50, round(0.5 * (dim(X)[1] - dim(X)[2]))))
y |
a numeric vector of data |
X |
a design matrix |
emu |
(logical) estimate mean of coefficient (TRUE) or assume it is zero (FALSE)? |
dof |
degrees of freedom to use for the t-quantiles (the remainder go to adaptive estimation of the prior) |
This function computes the adaptive FAB confidence interval for each coefficient in a linear regression model.
A list (s,sigma2,tau2,mu) where
s an estimate of the error standard deviation
sigma2 an estimate of the error variance, independent of s
tau2 an estimate of the coefficient variance, independent of s
mu an estimate of the coefficient mean, independent of s
Peter Hoff
Compute the adaptive FAB t-intervals for the coefficients of a regression model.
fabregCI(y, X, alpha = 0.05, dof = min(50, round(0.5 * (dim(X)[1] - dim(X)[2]))), verbose = TRUE)fabregCI(y, X, alpha = 0.05, dof = min(50, round(0.5 * (dim(X)[1] - dim(X)[2]))), verbose = TRUE)
y |
a numeric vector of data |
X |
a design matrix |
alpha |
the type I error rate, so 1-alpha is the coverage rate |
dof |
degrees of freedom to use for the t-quantiles (the remainder go to adaptive estimation of the prior) |
verbose |
logical, print progress or not |
This function computes the adaptive FAB confidence interval for each coefficient in a linear regression model.
A matrix where each row corresponds to the interval and OLS estimate of a coefficient.
Peter Hoff
Computation of a 1-alpha FAB t-interval
fabtCI(y, psi = c(0, 100, 1, 2), alpha = 0.05)fabtCI(y, psi = c(0, 100, 1, 2), alpha = 0.05)
y |
a numeric vector with at least two non-missing values |
psi |
a length-four vector of hyperparameters for the prior |
alpha |
the type I error rate, so 1-alpha is the coverage rate |
A FAB interval is the "frequentist" interval procedure
that is Bayes optimal: It minimizes the prior expected
interval width among all interval procedures with
exact 1-alpha frequentist coverage. This function computes
the FAB t-interval for the mean of a normal population with an
unknown variance, given a user-specified prior distribution
determined by psi. The prior is that the population mean
and variance are independently distributed as normal and
inverse-gamma random variables.
Referring to the elements of psi
as mu, t2, s20, nu0, the prior is determined as follows:
mu is the prior expectation of the mean
t2 is the prior variance of the mean
the population variance is inverse-gamma(nu0/2,nu0 s20/2)
Peter Hoff
y<-rnorm(10) fabtCI(y,c(0,10,1,5)) fabtCI(y,c(0,1/10,1,5)) fabtCI(y,c(2,10,1,5)) fabtCI(y,c(0,1/10,1,5))y<-rnorm(10) fabtCI(y,c(0,10,1,5)) fabtCI(y,c(0,1/10,1,5)) fabtCI(y,c(2,10,1,5)) fabtCI(y,c(0,1/10,1,5))
Computation of a 1-alpha FAB t-interval using z-optimal spending function
fabtzCI(y, s, dof, alpha = 0.05, psi = list(mu = 0, tau2 = 1e+05, sigma2 = 1))fabtzCI(y, s, dof, alpha = 0.05, psi = list(mu = 0, tau2 = 1e+05, sigma2 = 1))
y |
a numeric scalar, a normally distributed statistic |
s |
a numeric scalar, the standard error of y |
dof |
positive integer, degrees of freedom for s |
alpha |
the type I error rate, so 1-alpha is the coverage rate |
psi |
a list of parameters for the spending function, including
|
n<-10 y<-rnorm(n) fabtzCI(mean(y),sqrt(var(y)/n),n-1) t.test(y)$conf.intn<-10 y<-rnorm(n) fabtzCI(mean(y),sqrt(var(y)/n),n-1) t.test(y)$conf.int
Computation of a 1-alpha FAB z-interval
fabzCI(y, mu, t2, s2, alpha = 0.05)fabzCI(y, mu, t2, s2, alpha = 0.05)
y |
a numeric scalar |
mu |
a numeric scalar |
t2 |
a positive numeric scalar |
s2 |
a positive numeric scalar |
alpha |
the type I error rate, so 1-alpha is the coverage rate |
A FAB interval is the "frequentist" interval procedure
that is Bayes optimal: It minimizes the prior expected
interval width among all interval procedures with
exact 1-alpha frequentist coverage. This function computes
the FAB z-interval for the mean of a normal population with an
known variance, given a user-specified prior distribution
determined by psi. The prior is that the population mean
is normally distributed.
Referring to the elements of psi
as mu, t2, s2, the prior and population variance are
determined as follows:
mu is the prior expectation of the mean
t2 is the prior variance of the mean
s2 is the population variance
Peter Hoff
y<-0 fabzCI(y,0,10,1) fabzCI(y,0,1/10,1) fabzCI(y,2,10,1) fabzCI(y,0,1/10,1)y<-0 fabzCI(y,0,10,1) fabzCI(y,0,1/10,1) fabzCI(y,2,10,1) fabzCI(y,0,1/10,1)
Estimate across-group heterogeneity of means and variances
hhetmodel(y, g)hhetmodel(y, g)
y |
a numeric vector of data |
g |
a group membership vector, of the same length as y |
This function estimates parameters in a hierarchical model for normally distributed groups, where the across-group model for means is normal and the across group model for variances is inverse-gamma.
A vector (mu,t2,s20,nu0), where
mu is the mean of the group means
t2 is the variance of the group means
the the distribution of group variances is inverse-gamma(nu0/2,nu0 s20/2)
Peter Hoff
Estimate across-group heterogeneity of means
hhommodel(y, g, group, p1)hhommodel(y, g, group, p1)
y |
a numeric vector of data |
g |
a group membership vector, of the same length as y |
group |
the index of the group |
p1 |
number of groups used to pool sample variance |
This function estimates parameters in a hierarchical model for normally distributed groups, where the across-group model for means is normal and the variance is the same across groups.
A vector (s2,df,muw,t2w,s2w), where
s2 is the pooled variance
df is the degree of freedom of the t-quantiles
muw is the estimate mean of the group means
t2w is the estimate variance of the group means
s2w is the estimate within-group variance
Chaoyu Yu
Computation of 1-alpha FAB t-intervals for heteroscedastic multigroup data.
multifabCI(y, g, alpha = 0.05)multifabCI(y, g, alpha = 0.05)
y |
a numeric vector of data |
g |
a group membership vector, of the same length as y |
alpha |
the type I error rate, so 1-alpha is the coverage rate |
For each group j, this function computes an estimate of the parameters in a hierarchical model for means and variances from data other than group j, and uses this information to construct a FAB t-interval for group j. These intervals have 1-alpha frequentist coverage, assuming within-group normality.
Peter Hoff
## -- simulated data p<-10 ; n<-10 y<-rnorm(n*p) ; g<-rep(1:p,n) ## -- more interesting data takes longer # data(radon) ; y<-radon[,2] ; g<-radon[,1] ## -- FAB t-intervals FCI<-multifabCI(y,g) ## -- UMAU t-intervals ybar<-tapply(y,g,mean) ; ssd<-tapply(y,g,sd) ; n<-table(g) qtn<-cbind( qt(.025,n-1), qt(.975,n-1) ) UCI<-sweep(sweep(qtn,1,ssd/sqrt(n),"*"),1,ybar,"+") mean( (UCI[,2]-UCI[,1])/(FCI[,2]-FCI[,1]) , na.rm=TRUE)## -- simulated data p<-10 ; n<-10 y<-rnorm(n*p) ; g<-rep(1:p,n) ## -- more interesting data takes longer # data(radon) ; y<-radon[,2] ; g<-radon[,1] ## -- FAB t-intervals FCI<-multifabCI(y,g) ## -- UMAU t-intervals ybar<-tapply(y,g,mean) ; ssd<-tapply(y,g,sd) ; n<-table(g) qtn<-cbind( qt(.025,n-1), qt(.975,n-1) ) UCI<-sweep(sweep(qtn,1,ssd/sqrt(n),"*"),1,ybar,"+") mean( (UCI[,2]-UCI[,1])/(FCI[,2]-FCI[,1]) , na.rm=TRUE)
Computation of 1-alpha FAB t-intervals for homoscedastic multigroup data.
multifabCIhom(y, g, alpha = 0.05, prop = 0.5)multifabCIhom(y, g, alpha = 0.05, prop = 0.5)
y |
a numeric vector of data |
g |
a group membership vector, of the same length as y |
alpha |
the type I error rate, so 1-alpha is the coverage rate |
prop |
the proportion of groups to obtain the sample variance estimate |
For each group j, this function computes an estimate of the parameters in a hierarchical model for means using data from other groups, and uses this information to construct a FAB t-interval for group j. These intervals have 1-alpha frequentist coverage, assuming within-group normality and that the within group variance is the same across groups.
Chaoyu Yu
## -- simulate the data mu = 0; sigma2 = 10; tau2 = 1; p =100; theta = rnorm(p,mu,sqrt(tau2)) ns = round(runif(p,2,18)) Y=c() for(i in 1:p){ d2 = rnorm(ns[i],theta[i],sqrt(sigma2)) d1 = rep(i,ns[i]) d = cbind(d1,d2) Y = rbind(Y,d)} y = Y[,2] g = Y[,1] ## -- FAB t-intervals FCI = multifabCIhom(y,g) ## -- UMAU t-intervals ybar<-tapply(y,g,mean) ; ssd<-tapply(y,g,sd) ; n<-table(g) qtn<-cbind( qt(.025,n-1), qt(.975,n-1) ) UCI<-sweep(sweep(qtn,1,ssd/sqrt(n),"*"),1,ybar,"+") mean( (UCI[,2]-UCI[,1])/(FCI[,2]-FCI[,1]) , na.rm=TRUE)## -- simulate the data mu = 0; sigma2 = 10; tau2 = 1; p =100; theta = rnorm(p,mu,sqrt(tau2)) ns = round(runif(p,2,18)) Y=c() for(i in 1:p){ d2 = rnorm(ns[i],theta[i],sqrt(sigma2)) d1 = rep(i,ns[i]) d = cbind(d1,d2) Y = rbind(Y,d)} y = Y[,2] g = Y[,1] ## -- FAB t-intervals FCI = multifabCIhom(y,g) ## -- UMAU t-intervals ybar<-tapply(y,g,mean) ; ssd<-tapply(y,g,sd) ; n<-table(g) qtn<-cbind( qt(.025,n-1), qt(.975,n-1) ) UCI<-sweep(sweep(qtn,1,ssd/sqrt(n),"*"),1,ybar,"+") mean( (UCI[,2]-UCI[,1])/(FCI[,2]-FCI[,1]) , na.rm=TRUE)
Radon levels in 919 homes from 85 Minnesota counties
data(radon)data(radon)
A numeric matrix
http://www.stat.columbia.edu/~gelman/arm/software/
Compute Bayes optimal spending function
sfabz(theta, psi, alpha = 0.05)sfabz(theta, psi, alpha = 0.05)
theta |
value of theta being tested |
psi |
a list of parameters for the spending function, including
|
alpha |
level of test |
This function computes the value of s that minimizes the acceptance probability of a biased level-alpha test for a normal population with known variance, under a specified prior predictive distribution.
Peter Hoff
Compute the usual t-intervals for the coefficients of a regression model
umauregCI(y, X, alpha = 0.05)umauregCI(y, X, alpha = 0.05)
y |
a numeric vector of data |
X |
a design matrix |
alpha |
the type I error rate, so 1-alpha is the coverage rate |
This function computes the 'usual' uniformly most accurate unbiased confidence interval for each coefficient in a linear regression model.
A matrix where each row corresponds to the interval and OLS estimate of a coefficient.
Peter Hoff