ftest.pow: 'ftest' power function

In olli0601/abc.star: abc.star

Description

Compute the power of the one-sample multivariate equivalence test for population means of multivariate normal summary values with unknown population variance.

Usage

ftest.pow(rho, tau, n.of.y, p, alpha = 0.01, support = c(0, Inf),
  log = FALSE, norm = 1)

Arguments

rho	Vector of quantiles
tau	Upper boundary point of the equivalence region
p	Number of variables
alpha	Level of the equivalence test
support	Support of the truncated power function (vector of dimension 2).
log	If TRUE, the power function is returned on the log scale.

Note

The power function can be truncated to support.

References

http://arxiv.org/abs/1305.4283

Examples

# power function of the F-test, to test equality of means for multivariate
# normal samples with unknown covariance matrix

#set number of variables (i.e. summary statistics)
p <- 3
#set number of simulations
n <- 100

#calculate power for fixed equivalence value
tau <- 1.2
rho <- seq(0, .5, length = 1024)
ftest.pow(rho, tau, n = n, p = p)

# power increases as size of equivalence region increases but power function
# flattens out as equivalence region gets large
tmp <- lapply(c(0.05, 0.1, 0.2, 0.3), function(tau)
		{
			data.table(tau = as.factor(tau), rho = rho, power = ftest.pow(rho, tau, n, p, alpha = 0.01))
		})
tmp	<- do.call('rbind', tmp)
pp <- ggplot(tmp, aes(x = rho, y = power, colour = tau, group = tau)) + geom_line() + labs(y = 'Power\n(ABC acceptance probability)')
print(pp)

# power increases as number of simulations increase
tau 	<- 0.2
rho		<- seq(0, .3, length = 1024)
tmp  	<- lapply(c(25, 50, 100, 200, 400), function(n)
		{
			data.table(n = as.factor(n), rho = rho, y = ftest.pow(rho, tau, n, p, alpha = 0.01), d='power')
		})
tmp		<- do.call('rbind', tmp)
pp 		<- ggplot(tmp, aes(x = rho, y = y, colour = n, group = n)) + geom_line() + labs(y = 'Power\n(ABC acceptance probability)')
print(pp)

# add likelihood density to last power plot
t2.x	<- 0.25
tmp		<- rbind(tmp, data.table(n=n, rho=rho, y=ftest.sulkl(rho, t2.x, n, p, norm = 1, support= c(0,Inf), log=FALSE), d='prtl.lkl'))
pp 		<- ggplot(tmp, aes(x = rho, y = y, colour = n, linetype=d, group = interaction(n,d))) + geom_line() + labs(y = 'Power\n(ABC acceptance probability)')
print(pp)

olli0601/abc.star documentation built on May 24, 2019, 12:53 p.m.