equivtest: Multivariate equivalence testing

In ecopower: Power Estimates and Equivalence Testing for Multivariate Data

Multivariate equivalence testing

Description

equivtest takes in a copula model fitted to data and a matrix of effect sizes to execute a a multivariate equivalence test.

Usage

## S3 method for class 'cord'
equivtest(
  object,
  coeffs,
  term = NULL,
  object0 = NULL,
  stats = NULL,
  test = "LR",
  nsim = 999,
  ncores = detectCores() - 1,
  show.time = TRUE
)

equivtest(
  object,
  coeffs,
  term = NULL,
  object0 = NULL,
  stats = NULL,
  test = "LR",
  nsim = 999,
  ncores = detectCores() - 1,
  show.time = TRUE
)

Arguments

object	objects of class cord, typically the result of a call to cord.
coeffs	Coefficient matrix for a manyglm object that characterises the size of effects to be simulated. See effect_alt for help in producing this matrix.
term	Name of predictor of interest in quotes. Defaults to NULL, see details.
object0	object of class cord that specifies the null hypothesis. Defaults to NULL, see details.
stats	Statistics simulated under the null hypothesis. Optional, defaults to NULL. If not NULL, equivtest will not simulate test statistics and use the stats specified.
test	Test statistic for computing p-value. Defaults to "LR".
nsim	Number of simulations for p-value estimate to be based upon. Defaults to 999.
ncores	Number of cores for parallel computing. Defaults to the total number of cores available on the machine minus 1.
show.time	Logical. Displays time elapsed. Defaults to TRUE.

Details

equivtest takes a cord object and a coefficient matrix coeffs which specifies an effect size of interest to perform an equivalence test.

First, marginal parameters of the data are obtained from a manyglm object. Next, a copula model is fitted using cord to estimate the factor analytic covariance structure of the data. The cord function uses two factors by default. The p-value is then obtained by parsing the cord object into extend, nsim times with an effect size specified by coeffs.

The test statistics are simulated under the hypothesis that the effect size equals a certain threshold. The p-value is computed as the proportion of times the simulated test statistics are less than the observed statistic. Equivalence is declared if the estimated effect is less than the threshold.

equivtest can handle any user-defined null hypothesis, so only the fitted null model (object0) or the predictor of interest (term) needs to be specified. If both object0 and term are NULL, equivtest will automatically set the predictor of interest as the last term in the fitted object model or drop the only term in the model to obtain the intercept model.

Simulations are computed in parallel using the "socket" approach, which uses all available cores minus 1 for clustering to improve computation efficiency. Using 1 less than the number of available cores for your machine (detectCores()-1) is recommended to leave one core available for other computer processes.

Value

Equivalence test results, and;

p	p-value;
stat_obs	observed statistic;
stats	simulated statistics.

Functions

• equivtest(): Multivariate equivalence testing

See Also

effect_alt

Examples

library(ecoCopula)
library(mvabund)
data(spider)
spiddat = mvabund(spider$abund)
X = data.frame(spider$x)

# Specify increasers and decreasers
increasers = c("Alopacce", "Arctlute", "Arctperi", "Pardnigr", "Pardpull")
decreasers = c("Alopcune", "Alopfabr", "Zoraspin")

# Equivalence test for continuous predictor at effect_size=1.5
fit.glm = manyglm(spiddat~bare.sand, family="negative.binomial", data=X)
threshold = effect_alt(fit.glm, effect_size=1.5,
       increasers, decreasers, term="bare.sand")
fit.cord = cord(fit.glm)
equivtest(fit.cord, coeffs=threshold, term="bare.sand", nsim=99, ncores=2)

# Equivalence test for categorical predictor with 4 levels at effect_size=1.5
X$Treatment = rep(c("A","B","C","D"),each=7)
fit_factors.glm = manyglm(spiddat~Treatment, family="negative.binomial", data=X)
threshold = effect_alt(fit_factors.glm, effect_size=1.5,
       increasers, decreasers, term="Treatment")
fit_factors.cord = cord(fit_factors.glm)
equivtest(fit_factors.cord, coeffs=threshold, term="Treatment", nsim=99, ncores=2)

# Specify object0
object0.glm = manyglm(spiddat~1, family="negative.binomial")
object0.cord = cord(object0.glm)
equivtest(fit_factors.cord, coeffs=threshold, object0=object0.cord, nsim=99, ncores=2)

ecopower documentation built on Aug. 24, 2023, 5:09 p.m.