knitr::opts_chunk$set(echo = TRUE)

R Markdown

Synthetic simulations to assess issues with pbj F-test.

library(RandomFields)
library(RNifti)
library(pbj)
library(parallel)

x = 1:100
y = 1:100
n=50
model = RMexp(var=5, scale=10)
mask = updateNifti(matrix(1, nrow=length(x), ncol=length(y)) )
nsim=2

cfts = c(0.1, 0.25)
pvals = as.data.frame(matrix(NA, nrow=nsim, ncol=length(cfts)))
names(pvals) = paste0('cft', cfts)

pvals =mclapply(1:nsim, function(sim){
  images = RFsimulate(model, x = x, y=y, n=n)
  images = updateNifti(simplify2array(lapply(images@data, matrix, nrow=length(x))))

  data = data.frame(x = rnorm(n), Winv = runif(n))
  statmap = lmPBJ(images, form = ~ x + I(x^2) + I(x^3),
                  formred = ~ x, mask = mask,
                  data = data, Winv = data$Winv)
  pbjres = pbjSEI(statmap, cfts.s = cfts, nboot = 500)
  sapply(pbjres[grep('cft', names(pbjres))], function(x) min(x$pvalues))
} )


simonvandekar/NIsim documentation built on Oct. 12, 2020, 5:06 p.m.