gof_power: Power estimation of goodness-of-fit tests.

View source: R/gof_power.R

gof_powerR Documentation

Power estimation of goodness-of-fit tests.

Description

Find the power of various goodness-of-fit tests.

Usage

gof_power(
  pnull,
  rnull,
  ralt,
  param_alt,
  phat = function(x) -99,
  dnull = function(x) -99,
  TS,
  TSextra,
  With.p.value = FALSE,
  alpha = 0.05,
  Ranges = matrix(c(-Inf, Inf, -Inf, Inf), 2, 2),
  nbins = c(5, 5),
  minexpcount = 5,
  rate = 0,
  SuppressMessages = FALSE,
  maxProcessor,
  B = 1000
)

Arguments

pnull

function to find cdf under null hypothesis

rnull

function to generate data under null hypothesis

ralt

function to generate data under alternative hypothesis

param_alt

vector of parameter values for distribution under alternative hypothesis

phat

=function(x) -99, function to estimate parameters from the data, or -99

dnull

=function(x) -99, density function under the null hypothesis, if available, or -99 if missing

TS

user supplied function to find test statistics

TSextra

list provided to TS (optional)

With.p.value

=FALSE does user supplied routine return p values?

alpha

=0.05, the level of the hypothesis test

Ranges

=matrix(c(-Inf, Inf, -Inf, Inf),2,2), a 2x2 matrix with lower and upper bounds, if any, for chi-square tests

nbins

=c(5, 5), number of bins for chi square tests.

minexpcount

=5 minimal expected bin count required for chi square tests.

rate

=0 rate of Poisson if sample size is random, 0 if sample size is fixed

SuppressMessages

=FALSE, should informative messages be shown?

maxProcessor

maximum of number of processors to use, 1 if no parallel processing is needed or number of cores-1 if missing

B

=1000 number of simulation runs

Details

For details on the usage of this routine consult the vignette with vignette("MDgof","MDgof")

Value

A numeric matrix of power values.

Examples

# All examples are run with B=10 and maxProcessor=1 to pass CRAN checks.
# This is obviously MUCH TO SMALL for any real usage.
# Power of tests if null hypothesis specifies a bivariate standard normal 
# distribution but data comes from a bivariate normal with different means, 
# without parameter estimation.
rnull=function() mvtnorm::rmvnorm(100, c(0, 0))
ralt=function(p) mvtnorm::rmvnorm(100, c(p, p))
pnull=function(x) {
  if(!is.matrix(x)) return(mvtnorm::pmvnorm(rep(-Inf, 2), x))
  apply(x, 1, function(x) mvtnorm::pmvnorm(rep(-Inf, 2), x))
}
gof_power(pnull, rnull, ralt, c(0, 1), B=10, maxProcessor = 1)
# Same as above, but now with density included
dnull=function(x) {
  if(!is.matrix(x)) return(mvtnorm::dmvnorm(x))
  apply(x, 1, function(x) mvtnorm::dmvnorm(x))
}
gof_power(pnull, rnull, ralt, c(0, 1), dnull=dnull, B=10, maxProcessor = 1)
# Power of tests when null hypothesis specifies a bivariate normal distribution, 
# with mean parameter estimated, wheras data comes from a t distribution
rnull=function(p) mvtnorm::rmvnorm(100, p)
ralt=function(df) mvtnorm::rmvt(100, sigma=diag(2), df=df)
pnull=function(x,p) {
  if(!is.matrix(x)) return(mvtnorm::pmvnorm(rep(-Inf, 2), x, mean=p))
  apply(x, 1, function(x) mvtnorm::pmvnorm(rep(-Inf, 2), x, mean=p))
}
dnull=function(x, p) {
  if(!is.matrix(x)) return(mvtnorm::dmvnorm(x, mean=p))
  apply(x, 1, function(x) mvtnorm::dmvnorm(x, mean=p))
}
phat=function(x) apply(x, 2, mean)
gof_power(pnull, rnull, ralt, c(50, 5), dnull=dnull, phat=phat, B=10, maxProcessor = 1)
# Example of a discrete model, with parameter estimation
# Under null hypothesis: X~Bin(10, p), Y|X=x~Bin(5, 0.5+x/100)
# Under alternative hypothesis: X~Bin(10, p), Y|X=x~Bin(5, K+x/100)
rnull=function(p=0.5) {
  x=stats::rbinom(1000, 10, p)
  y=stats::rbinom(1000, 5, 0.5+x/100)
  MDgof::sq2rec(table(x, y))
}
ralt=function(K=0.5) {
  x=stats::rbinom(1000, 10, 0.5)
  y=stats::rbinom(1000, 5, K+x/100)
  MDgof::sq2rec(table(x, y))
}
pnull=function(x, p) {
  f=function(x) sum(dbinom(0:x[1], 10, p[1])*pbinom(x[2], 5, 0.5+0:x[1]/100))
  if(!is.matrix(x)) x=rbind(x)
  apply(x, 1, f)
}
phat=function(x) {
  tx=tapply(x[,3], x[,1], sum)
  mean(rep(as.numeric(names(tx)), times=tx))/10
}
gof_power(pnull, rnull, ralt, c(0.5, 0.6), phat=phat, B=10, maxProcessor = 1)

MDgof documentation built on Feb. 13, 2026, 1:06 a.m.