pbjInference: Performs (semi)Parametric Bootstrap Joint ((s)PBJ) Inference

In simonvandekar/pbj: Parametric Bootstrap Joint Testing Procedures for Neuroimaging

Performs (semi)Parametric Bootstrap Joint ((s)PBJ) Inference

Description

Performs (semi)Parametric Bootstrap Joint ((s)PBJ) Inference

Usage

pbjInference(
  statMap,
  statistic = mmeStat,
  null = TRUE,
  nboot = 5000,
  rboot = function(n) {
     (2 * stats::rbinom(n, size = 1, prob = 0.5) - 1)
 },
  method = c("wild", "permutation", "nonparametric"),
  runMode = c("cdf", "bootstrap"),
  progress = FALSE,
  rdata_rds = NULL,
  cft_s = NULL,
  cft_p = NULL,
  cft_chisq = NULL,
  mc.cores = 1,
  mc.preschedule = FALSE,
  ...
)

Arguments

statMap	statMap object as obtained from lmPBJ function.
statistic	A user specified function that takes an RNifti image object and computes a particular statistic of interest. There are several provided in the pbj package. See referenced functions below.
null	Perform bootstrap under the null or alternative hypotheses?
nboot	Number of bootstrap samples to use.
rboot	Function for generating random variables. Should return an n vector. Defaults to Rademacher random variable.
method	Character, method to use for resampling procedure. Wild bootstrap, permutation, or nonparametric
runMode	character, that controls output. cdf returns the empirical CDFs, bootstrap returns the bootstrapped statistics as a list.
progress	Logical indicating whether to track progress with a progress bar.
rdata_rds	A character string specifying a .rdata or .rds file to write the output. If specified, resampling is run as background process.
cft_s	A vector of robust effect size index (RESI) cluster forming thresholds
cft_p	A vector of p-value cluster forming thresholds
cft_chisq	A vector of p-value cluster forming thresholds
mc.cores	Number of cores to parallelize bootstrap.
mc.preschedule	Logical, where to preschedule jobs see [parallel()]. FALSE is substantially faster.
...	arguments passed to statistic function.

Details

This function runs resampling-based methods to perform inference on topological features of neuroimaging statistics. The topological feature is determined by the 'statistic' function. Several exist in the pbj package (see below). The statistic function should take an image as input and compute some topological feature from the image. returned as a list. Multiple topological features can be returned, as in [mmeStat()] and [cluster()]. #' To use default methods the statistic must have a logical 'rois' argument that outputs an integer valued image identifying where each topological features is located.

Value

Returns the statMap object, with a pbj object added. If runMode=='cdf', the first element is the observed statistic value, and the subsequent elements are the CDFs and ROIs, used for computing adjusted p-values and plotting. If runMode=='bootstrap', the first element is the observed statistic value and the second is a list of the boostrap values.

See Also

[mmeStat()], [maxima()], and [cluster()] for statistic functions. See [lmPBJ()] to create statMap objects. See [image.statMap()], and [table.statMap()] for producing summaries of the results.

Examples

# loading example data
library(pain21)
pain = pain21()
pdata = pain$data

# fitting regression of images onto study sample size,
# weights proportional to study sample size
pbjModel2 = lmPBJ(images=pdata$images, form=~n, formred=~1,
                  W = pdata$n, mask=pain$mask, data=pdata)
table.statMap(pbjModel2, 'maxima')


# p-value thresholding for cluster extent/mass
# not run, parallel processing
# pbjModelAll <- pbjInference(pbjModel2, nboot=50, cft_p=0.01,
#                             CEI=TRUE, CMI=TRUE, max=TRUE,
#                             mc.cores=ceiling(parallel::detectCores()/2))
pbjModelAll <- pbjInference(pbjModel2, nboot=5, cft_p=0.01,
                            CEI=TRUE, CMI=TRUE, max=TRUE)
head(table.statMap(pbjModelAll, method = 'maxima'))
head(table.statMap(pbjModelAll, method = 'CEI', cft_p=0.01))
# returns the first threshold used
head(table.statMap(pbjModelAll, method = 'CMI', cft_p = 0.01))
# RESI effect size thresholding for cluster extent/mass
pbjModel2 = pbjInference(pbjModel2, nboot=2, cft_s=c(0.1, 0.25), CMI=TRUE)
head(table.statMap(pbjModel2, method = 'CEI', cft_s=0.1))
head(table.statMap(pbjModel2, method = 'CEI', cft_s=0.25))
# Inference on local maxima
pbjModel2 = pbjInference(pbjModel2, nboot=2, statistic = maxima)
# Cluster extent inference
pbjModelCEI = pbjInference(pbjModel2, nboot=2, statistic = cluster, cft_s=0.1)

# not run
# outfile = paste0(tempfile(), '.rds')
# run in parallel in the background
# bgInfo <- pbjInference(pbjModel2, rdata_rds = outfile, nboot=50, cft_p=0.01,
# CEI=TRUE, CMI=TRUE, max=TRUE, mc.cores=ceiling(parallel::detectCores()/2))
# statMap = readRDS(outfile)

simonvandekar/pbj documentation built on Nov. 3, 2023, 9:33 a.m.