ashr (version 0.9)

ashci: Credible Interval Computation for the ash object

Description

Given the ash object return by the main function ash, this function computes the corresponding credible interval of the mixture model.

Usage

ashci(a, level = 0.95, betaindex, lfsrcriteria = 0.05, tol = 1e-05, maxcounts = 100, shrinkingcoefficient = 0.9, trace = FALSE, ncores = FALSE)

Arguments

a

the fitted ash object

level

the level for the credible interval, (default=0.95)

betaindex

a vector consisting of locations of betahat where you would like to compute the credible interval

lfsrcriteria

a scalar, in which the function would autoselect betahats based on lfsr value smaller than lfsrcriteria when index is not supplied. Setting it to 1 would compute credible interval for all observations.

tol

the desired accuracy, default value is 1e-5.

maxcounts

a positive integer used as the maximum iterations to carry additional optimization on a refined interval,when you see Actual cdf deviates from (1-level)/2,(1+level)/2 by too much, try to increase this number. (default is 100)

shrinkingcoefficient

A positive real number smaller than 1 (default is 0.9), used to shrink to search interval for lower and upper bound of the credible interval

trace

a logical variable denoting whether some of the intermediate results of iterations should be displayed to the user. Default is FALSE.

ncores

Whether to use parallel computing, defaults to FALSE, user could specify number of cores they would like to use. Further, if user does not specify and we have over 1000 entries to compute, then the function would perform parallel computation using number of CPU cores on the current host.

Value

A matrix, with first column being the location,second column being lfsr, 3rd column being posterior mean,4th and 5th column being the lower bound and upper bound for the credible interval.

Details

Uses default optimization function and perform component-wise credible interval computation. The computation cost is linear of the length of betahat.

Examples

beta = c(rep(0,20),rnorm(20))
sebetahat = abs(rnorm(40,0,1))
betahat = rnorm(40,beta,sebetahat)
beta.ash = ash(betahat, sebetahat)

CImatrix=ashci(beta.ash,level=0.95)
print(CImatrix)
print(CImatrix[order(CImatrix[,2]),]) # Sorted according to the lfsr

CImatrix1=ashci(beta.ash,level=0.95,betaindex=c(1,2,5))
CImatrix2=ashci(beta.ash,level=0.95,lfsrcriteria=0.1)
#CImatrix3=ashci(beta.ash,level=0.95, betaindex=c(1:length(beta)),ncores=4)
print(CImatrix1)
print(CImatrix2)
#print(CImatrix3)

##A larger example
#beta = c(rep(0,1000),rnorm(1000))
#sebetahat = abs(rnorm(2000,0,1))
#betahat = rnorm(2000,beta,sebetahat)
#beta.ash = ash(betahat, sebetahat)
#CImatrix4 = ashci(beta.ash,level=0.95, betaindex=c(1:length(beta)),ncores=4)