| Type: | Package |
| Title: | Simulate fMRI Data |
| Version: | 0.2-14 |
| Date: | 2023-10-18 |
| Depends: | R (≥ 3.1.1), deSolve |
| Description: | Generates functional Magnetic Resonance Imaging (fMRI) time series or 4D data. Some high-level functions are created for fast data generation with only a few arguments and a diversity of functions to define activation and noise. For more advanced users it is possible to use the low-level functions and manipulate the arguments. See Welvaert et al. (2011) <doi:10.18637/jss.v044.i10>. |
| License: | GPL-2 | GPL-3 [expanded from: GPL (≥ 2)] |
| LazyLoad: | yes |
| Repository: | CRAN |
| Date/Publication: | 2023-10-18 12:50:02 UTC |
| Packaged: | 2023-10-18 09:49:35 UTC; tabelow |
| NeedsCompilation: | yes |
| Author: | Marijke Welvaert [aut], Joke Durnez [ctb], Beatrijs Moerkerke [ctb], Yves Rosseel [ctb], Karsten Tabelow [ctb, cre], Geert Verdoolaege [ctb] |
| Maintainer: | Karsten Tabelow <karsten.tabelow@wias-berlin.de> |
Functions to Generate fMRI Data Including Activated Data, Noise Data and Resting State Data
Description
The package allows users to generate fMRI time series or 4D data. Some high-level functions are created for fast data generation with only a few arguments and a diversity of functions to define activation and noise. For more advanced users it is possible to use the low-level functions and manipulate the arguments.
Author(s)
Marijke Welvaert with contributions from Joke Durnez, Beatrijs Moerkerke, Yves Rosseel, Karsten Tabelow, and Geert Verdoolaege
Maintainer: Karsten Tabelow <karsten.tabelow@wias-berlin.de>
References
Welvaert, M., Durnez, J., Moerkerke, B., Verdoolaege, G. and Rosseel, Y. (2011). neuRosim: An R Package for Generating fMRI Data. Journal of Statistical Software, 44(10), 1–18
Examples
## Generate fMRI time series for block design
design <- simprepTemporal(totaltime=200, onsets=seq(1,200,40),
durations=20, TR=2, effectsize=1, hrf="double-gamma")
ts <- simTSfmri(design=design, SNR=1, noise="white")
plot(ts, type="l")
## Generate fMRI slice for block design with activation in 2 regions
design <- simprepTemporal(totaltime=200, onsets=seq(1,200,40),
durations=20, TR=2, effectsize=1, hrf="double-gamma")
region <- simprepSpatial(regions=2, coord=list(c(32,15),c(57,45)),
radius=c(10,7), form="sphere")
out <- simVOLfmri(design=design, image=region, dim=c(64,64),
SNR=1, noise="none")
plot(out[32,15,], type="l")
Calculates a discrete Gaussian smoothing kernel (adopted from AnalyzeFMRI)
Description
Calculates a simple, discrete Gaussian smoothing kernel of a specifice size given the covariance matrix of the Gaussian.
Usage
GaussSmoothKernel(voxdim=c(1,1,1), ksize=5, sigma=diag(3,3))
Arguments
voxdim |
The dimensions of each voxel. |
ksize |
The size (in voxels) of the kernel with which to filter the independent field. |
sigma |
The covariance matrix of the Gaussian kernel. |
Value
An array of dimension (ksize,ksize,ksize) containing the smoothing kernel.
Author(s)
J. L. Marchini
See Also
Examples
a <- GaussSmoothKernel(voxdim=c(1,1,1), ksize=5, sigma=diag(1,3))
Simulate a GRF (adopted from AnalyzeFMRI)
Description
Simulates a Gaussian Random Field with specified dimensions and covariance structure.
Usage
Sim.3D.GRF(d, voxdim, sigma, ksize, mask=NULL, type=c("field","max"))
Arguments
d |
A vector specifying the dimensions of a 3D or 4D array. |
voxdim |
The dimensions of each voxel. |
sigma |
The 3D covariance matrix of the field. |
ksize |
The size (in voxels) of the kernel with which to filter the independent field. |
mask |
A 3D mask for the field. |
type |
If |
Details
The function works by simulating a Gaussian r.v at each voxel location and the smoothing the field with a discrete filter to obtain a field with the desired covariance structure.
Value
mat |
Contains the simulated field if |
max |
The maximum value of the simulated field |
Author(s)
J. L. Marchini
See Also
Examples
d <- c(64, 64, 21)
FWHM <- 9
sigma <- diag(FWHM^2, 3) / (8 * log(2))
voxdim <- c(2, 2, 4)
msk <- array(1, dim = d)
field <- Sim.3D.GRF(d = d, voxdim = voxdim, sigma = sigma,
ksize = 9, mask = msk, type = "max")
Balloon model
Description
Generates the BOLD signal based on the Balloon model of Buxton et al. (2004).
Usage
balloon(stim, totaltime, acc, par=list(), verbose=TRUE)
Arguments
stim |
Vector representing the presence/absence (1-0 coding) of a stimulus/activation in seconds. |
totaltime |
Total duration of stimulus vector in seconds. |
acc |
Microtime resolution of stimulus vector in seconds. |
par |
List representing the parameters of the Balloon model. The list should contain the following:
|
verbose |
If |
Details
Based on the provided stimulus boxcar function, a neural activation function is generated that enters the Balloon model to generate a BOLD response. The microtime resolution ensures a high-precision generation of the response. More details can be found in Buxton et al. (2004).
Value
Vector representing the values of the BOLD signal for the given stimulus vector and microtime resolution.
Author(s)
G. Verdoolaege, M. Welvaert
References
Buxton, RB, Uludag, K, Dubowitz, DJ and Liu, TT (2004). Modeling the hemodynamic response to brain activation. NeuroImage, 23, S220-S233.
See Also
Examples
s <- rep(rep(0,10), rep(1,10), 5)
T <- 100
it <- 0.1
out <- balloon(s, T, it)
#takes a couple of seconds due to solving of the differential equations
Double-gamma Haemodynamic reponse function
Description
Specifies a double-gamma variate haemodynamic response function for the given time vector and parameters.
Usage
canonicalHRF(x, param = NULL, verbose = TRUE)
Arguments
x |
Time vector in seconds. |
param |
List of parameters of the haemodynamic response function. The list should contain the following:
|
verbose |
If |
Value
Vector representing the values of the function for the given time vector and parameters.
Author(s)
M. Welvaert
References
[1] Friston, KJ, Fletcher, P, Josephs, O, Holmes, AP, Rugg, MD and Turner, R (1998). Event-related fMRI: Characterising differential responses. NeuroImage, 7, 30-40.
[2] Glover, GH (1999). Deconvolution of impulse response in event-related BOLD fMRI. NeuroImage, 9, 416-429.
See Also
Examples
t <- 1:100
out <- canonicalHRF(t, verbose=FALSE)
Single Gamma Haemodynamic response function.
Description
Specifies a Gamma variate haemodynamic response function for the given time vector and FWHM.
Usage
gammaHRF(x, FWHM = 4, verbose = TRUE)
Arguments
x |
Time vector in seconds. |
FWHM |
Full Width Half Maximum of the Gamma variate function. |
verbose |
If |
Value
Vector representing the values of the function for the given time vector and FWHM.
Author(s)
M. Welvaert
References
Buxton, RB, Uludag, K, Dubowitz, DJ and Liu, TT (2004). Modeling the hemodynamic response to brain activation. NeuroImage, 23, S220-S233.
See Also
Examples
t <- 1:100
out <- gammaHRF(t, verbose=FALSE)
Generate low frequency drift
Description
Generates a low-frequency drift dataset with specified dimensions and frequency.
Usage
lowfreqdrift(dim, freq = 128, nscan, TR, template, verbose = TRUE)
Arguments
dim |
A vector specifying the dimensions of the image. |
freq |
The frequency of the drift in seconds. |
nscan |
The number of scans in the dataset. |
TR |
The repetition time in seconds. |
template |
An array representing the anatomical structure or mask with dimensions equal to dim. |
verbose |
Logical indicating if warnings should be printed. |
Details
The function generates low-frequency drift based on a basis set of cosine functions. The result is an array with specified dimensions and frequency.
Value
An array containing the drift with dimensions specified in dim.
Author(s)
Y. Rosseel, M. Welvaert
References
Friston et al. (2007). Statistical Parametric Mapping: The analysis of functional brain images. Academic Press.
See Also
temporalnoise, systemnoise, physnoise, tasknoise, spatialnoise
Examples
d <- c(10,10,10)
freq <- 80
nscan <- 100
TR <- 2
out <- lowfreqdrift(d, freq, nscan, TR, verbose=FALSE)
Generate physiological noise
Description
Generates a physiological noise dataset with specified dimensions and standard deviation. The physiological noise is defined as noise caused by heart beat and respiratory rate.
Usage
physnoise(dim, nscan, TR, sigma, freq.heart = 1.17,
freq.resp = 0.2, template, verbose = TRUE)
Arguments
dim |
A vector specifying the dimensions of the image. |
nscan |
The number of scans in the dataset. |
TR |
The repetition time in seconds. |
sigma |
The standard deviation of the noise. |
freq.heart |
The frequency in Hz of the heart beat. |
freq.resp |
The frequency in Hz of the respiratory rate. |
template |
An array representing the anatomical structure or mask with dimensions equal to dim. |
verbose |
Logical indicating if warnings should be printed. |
Details
The function generates physiological noise. Heart beat and respiratory rate are defined as sine and cosine functions with specified frequencies. Additional Gaussian noise creates variability over voxels. The result is a noise dataset with specified dimensions and desired standard deviation.
Value
An array containing the noise with dimensions specified in dim and nscan.
Author(s)
M. Welvaert
See Also
temporalnoise, lowfreqdrift, systemnoise, tasknoise, spatialnoise
Examples
d <- c(10,10,10)
sigma <- 5
nscan <- 100
TR <- 2
out <- physnoise(d, nscan, TR, sigma, verbose=FALSE)
The Rice Distribution
Description
Density and random generation for the Rician distribution
Usage
rrice(n, vee, sigma) Arguments
n |
number of observations. Must be a positive integer of length 1. |
vee |
non-centrality parameter of the distribution. Must be a positive integer of length 1. |
sigma |
scale parameter of the distribution. Must be a positive integer of length 1. |
Details
See VGAM for more details on the parameters and the formula of the probability density function.
Value
Random deviates for the given number of observations.
Author(s)
T.W. Yee
Examples
x <- rrice(n=10,vee=2,sigma=1)
Simulate fMRI time series
Description
Simulates an fMRI time series for the specified design and noise type.
Usage
simTSfmri(design = list(), base=0, nscan = NULL, TR = NULL, SNR=NULL,
noise = c("none", "white", "temporal", "low-frequency",
"physiological", "task-related", "mixture"), type = c("gaussian", "rician"),
weights, verbose = TRUE, rho = 0.2, freq.low = 128, freq.heart = 1.17,
freq.resp = 0.2, vee=1)
Arguments
design |
List generated by |
base |
Baseline value of the time series. |
nscan |
Number of scans. |
TR |
Repetition time in seconds. |
SNR |
Signal-to-noise ratio of the time series. |
noise |
Type of noise (white is default). |
type |
If |
weights |
If |
verbose |
Logical indicating if warnings should be returned. |
rho |
If |
freq.low |
If |
freq.heart |
If |
freq.resp |
If |
vee |
If |
Value
A vector representing the fMRI time series.
Author(s)
M. Welvaert
See Also
Examples
design <- simprepTemporal(totaltime=200, onsets=seq(1,200,40),
durations=20, effectsize=1, TR=2, hrf="double-gamma")
ts <- simTSfmri(design=design, SNR=1, noise="white")
plot(ts, type="l")
Simulate fMRI resting state time series
Description
Synthesizes a single time series x representing resting state activity. The fluctuation frequencies f are limited to a square passband 0.01 Hz <= f <= 0.1 Hz. TR is the repetition time (needed to compute the passband limits), expressed in seconds. N is the required number of samples (needs not be a power of 2).
Usage
simTSrestingstate(nscan, base=0, TR, SNR=NULL, noise = c("none", "white",
"temporal", "low-frequency", "physiological", "mixture"),
type = c("gaussian", "rician"), weights, verbose = TRUE, rho = 0.2,
freq.low = 128, freq.heart = 1.17, freq.resp = 0.2, vee=1)
Arguments
nscan |
Number of scans. |
base |
Baseline value of the time series. |
TR |
Repetition time in seconds. |
SNR |
Signal-to-noise ratio of the time series. |
noise |
Type of noise (white is default). |
type |
If |
weights |
If |
verbose |
Logical indicating if warnings should be returned. |
rho |
If |
freq.low |
If |
freq.heart |
If |
freq.resp |
If |
vee |
If |
Value
A vector representing the resting state time series
Author(s)
J. Durnez, G. Verdoolaege, M. Welvaert
References
[1] C.G. Fox, Computers & Geoscience, Vol. 13, pp. 369-374, 1987.
[2] M. Fukunaga, Magnetic Resonance Imaging, Vol. 24, pp. 979-992, 2006.
See Also
Examples
out <- simTSrestingstate(nscan=50, TR=2, SNR=1, noise="none")
plot(out, type="l")
Simulate 3D or 4D fMRI data
Description
Simulates a 3D or 4D fMRI dataset for the specified design and with activation in the specified regions.
Usage
simVOLfmri(design = list(), image = list(), base=0, dim, nscan = NULL,
TR = NULL, SNR=NULL, noise = c("none", "white", "temporal",
"spatial", "low-frequency", "physiological", "task-related",
"mixture"), type = c("gaussian", "rician"),
spat = c("corr", "gaussRF", "gammaRF"), weights, verbose = TRUE,
rho.temp = 0.2, rho.spat = 0.75, freq.low = 128,
freq.heart = 1.17, freq.resp = 0.2, FWHM = 4, gamma.shape = 6,
gamma.rate = 1, vee=1, template)
Arguments
design |
List generated by |
image |
List generated by |
base |
Baseline of the data. Should be a single number or an array with the same dimensions as in |
dim |
Dimensions of the image space. |
nscan |
Number of scans for noise images. |
TR |
Repetition time for noise images. |
SNR |
Signal-to-noise ratio. |
noise |
Type of noise, default is white. |
type |
If |
spat |
If |
weights |
If |
verbose |
Logical indicating if warning should be printed. |
rho.temp |
If |
rho.spat |
If |
freq.low |
If |
freq.heart |
If |
freq.resp |
If |
FWHM |
If |
gamma.shape |
If |
gamma.rate |
If |
vee |
If |
template |
An array representing the anatomical structure or mask with dimensions equal to dim. |
Value
A 3D or 4D array specifying the values for each voxel in the data.
Author(s)
M. Welvaert
See Also
simTSfmri, simprepTemporal, simprepSpatial
Examples
design <- simprepTemporal(totaltime=200, onsets=seq(1,200,40),
durations=20, TR=2, effectsize=1, hrf="double-gamma")
region <- simprepSpatial(regions=2, coord=list(c(32,15),c(57,45)),
radius=c(10,7), form="sphere", fading=TRUE)
out <- simVOLfmri(design=design, image=region, dim=c(64,64),
SNR=1, noise="none")
plot(out[32,15,], type="l")
image(1:64, 1:64, out[,,10], col = grey(0:255/255))
Prepare spatial structure of the data
Description
Prepare a list defining the necessary parameters to specify the spatial structure of the activation data.
Usage
simprepSpatial(regions, coord, radius = NULL,
form = c("cube", "sphere", "manual"), fading = 0)Arguments
regions |
Number of activated regions. |
coord |
List of coordinates specifying the xyz-coordinates. |
radius |
If form=cube or sphere, the distance between the center and the edge, if form=manual, the number of voxels in each region. |
form |
The form of the activated regions. |
fading |
Decay rate between 0 and 1. 0 means no fading, while 1 results in the fastest decay. |
Value
A list with the necessary arguments to be used in simVOLfmri.
Author(s)
M. Welvaert
See Also
simVOLfmri, simprepTemporal, specifyregion
Examples
coord <- list(c(3,3,3),c(6,6,6))
radius <- c(1,2)
out <- simprepSpatial(2, coord, radius, form="cube", fading=0.2)
Prepare temporal structure of the data
Description
Prepare a list defining the necessary parameters to specify the temporal structure of the activation data.
Usage
simprepTemporal(totaltime, regions = NULL, onsets, durations,
TR, effectsize, accuracy=0.1,
hrf = c("gamma", "double-gamma", "Balloon"),
param = NULL)
Arguments
totaltime |
Duration of the experiment. |
regions |
Number of regions. If not specified, it is assumed that all regions have the same design matrix. |
onsets |
List or vector representing the onsets of the stimulus in seconds. |
durations |
List or vector representing the durations of the stimulus in seconds. |
TR |
Repetition time in seconds. |
effectsize |
List or number representing the effectsize in each condition. |
accuracy |
Microtime resolution in seconds. |
hrf |
Haemodynamic response function (double-gamma is default) |
param |
Vector, matrix or array representing the parameters of the haemodynamic response function. |
Value
A list with the necessary arguments to be used in simVOLfmri or simTSfmri.
Author(s)
M. Welvaert
See Also
simVOLfmri, simTSfmri, simprepSpatial, specifyregion
Examples
ncond <- 2
os <- list(c(20,60),c(15,35))
d <- list(20, 10)
effect <- list(7,10)
total <- 80
TR <- 2
out <- simprepTemporal(total, onsets=os, durations=d, TR=TR,
effectsize=effect, hrf="double-gamma")
Generate spatially correlated noise
Description
Generates a spatially correlated noise dataset with specified dimensions and standard deviation.
Usage
spatialnoise(dim, sigma, nscan, method = c("corr", "gammaRF", "gaussRF"),
type=c("gaussian","rician"), rho = 0.75, FWHM = 4, gamma.shape = 6,
gamma.rate = 1, vee=1, template, verbose = TRUE)
Arguments
dim |
A vector specifying the dimensions of the image. |
sigma |
The standard deviation of the noise. |
nscan |
The number of scans in the dataset. |
method |
Method specifying the type of spatial correlation. Default is |
type |
Type of distribution if |
rho |
If |
FWHM |
If |
gamma.shape |
If |
gamma.rate |
If |
vee |
If |
template |
An array representing the anatomical structure or mask with dimensions equal to dim. |
verbose |
Logical indicating if warnings should be printed. |
Details
The function generates spatially correlated noise. When method=="corr", AR(1) voxelwise correlations are introduced.
If method=="gaussRF" of method=="gammaRF", respectively a Gaussian Random Field or a Gamma Random Field is created. The result is a noise array with specified dimensions and desired standard deviation.
The generation of the random fields is based on the function Sim.3D.GRF from J.L. Marchini in the package AnalyzeFMRI.
Value
An array containing the noise with dimensions specified in dim and nscan.
Author(s)
J. Durnez, B. Moerkerke, M. Welvaert
See Also
temporalnoise, lowfreqdrift, physnoise, tasknoise, systemnoise, Sim.3D.GRF
Examples
d <- c(10,10,10)
sigma <- 5
nscan <- 100
rhospat <- 0.7
out <- spatialnoise(d, sigma, nscan, method="corr", rho=rhospat, verbose=FALSE)
Generate design matrix.
Description
Generates a design matrix to be used as a model for the simulated activation.
Usage
specifydesign(onsets, durations, totaltime, TR, effectsize, accuracy=0.1,
conv = c("none", "gamma", "double-gamma", "Balloon"),
cond.names = NULL, param = NULL)
Arguments
onsets |
List or vector representing the onsets in seconds. |
durations |
List or vector representing the durations in seconds. |
totaltime |
Duration of the experiment in seconds. |
TR |
Repetition time in seconds. |
effectsize |
List or number representing the effectsize in each condition. |
accuracy |
Microtime resolution in seconds. |
conv |
Should the design matrix be convoluted, default is none. |
cond.names |
Optional names for the conditions. |
param |
Parameters of the haemodynamic response function. See |
Value
A matrix specifying the design.
Author(s)
M. Welvaert
See Also
specifyregion,gammaHRF,canonicalHRF,balloon
Examples
os <- list(c(20,60),c(15,35))
d <- list(20, 10)
total <- 80
TR <- 2
out <- specifydesign(os, d, total, TR, effectsize=list(2,5), conv="double-gamma")
Generate activation image
Description
Generates an image with activated regions for specified dimensions. The regions are defined by their center and radius or can be entered manually.
Usage
specifyregion(dim, coord, radius = NULL,
form = c("cube", "sphere", "manual"),
fading = 0)
Arguments
dim |
Dimensions of the image space. |
coord |
Coordinates of the activated region, if |
radius |
If |
form |
The form of the activated region. Default is |
fading |
Decay rate between 0 and 1. 0 means no fading, while 1 results in the fastest decay. |
Value
An array representing the activation image with specified regions.
Author(s)
M. Welvaert
See Also
specifyregion,gammaHRF,canonicalHRF,balloon
Examples
d <- c(10,10,10)
coord <- c(3,3,3)
radius <- 1
out <- specifyregion(d, coord, radius, form="sphere")
Generate a stimulus boxcar function.
Description
Generates a stimulus boxcar vector for the specified time duration and microtime resolution based on the user-defined onsets and durations.
Usage
stimfunction(totaltime, onsets, durations, accuracy)
Arguments
totaltime |
Total time of the design in seconds. |
onsets |
Vector representing the onsets of the stimulus in seconds. |
durations |
Vector representing the durations of the stimulus in seconds. |
accuracy |
Microtime resolution in seconds. |
Details
If duration is a single number, it is assumed that all stimulus onsets have the same duration.
Value
A vector in microtime resolution specifying the stimulus boxcar function in 1-0 coding.
Author(s)
M. Welvaert
See Also
Examples
total <- 100
os <- c(1, 21, 41, 61, 81)
d <- 10
out <- stimfunction(total, os, d, 0.1)
Generate system noise
Description
Generates a system noise dataset with specified dimensions and standard deviation. The noise can be either Gaussian or Rician distributed.
Usage
systemnoise(dim, nscan, type=c("gaussian","rician"), sigma, vee, template,
verbose = TRUE)
Arguments
dim |
A vector specifying the dimensions of the image. |
nscan |
The number of scans in the dataset. |
type |
Distribution of system noise. Default is gaussian. |
sigma |
The standard deviation of the noise. |
vee |
If |
template |
An array representing the anatomical structure or mask with dimensions equal to dim. |
verbose |
Logical indicating if warnings should be printed. |
Value
An array containing the noise with dimensions specified in dim and nscan.
Author(s)
M. Welvaert
See Also
temporalnoise, lowfreqdrift, physnoise, tasknoise, spatialnoise
Examples
d <- c(10,10,10)
sigma <- 5
nscan <- 100
out <- systemnoise(d, nscan, type="rician", sigma, verbose=FALSE)
Generate task-related noise
Description
Generates a Gaussian noise dataset with specified dimensions and standard deviation only when a task is performed or activation is present.
Usage
tasknoise(act.image, sigma, type=c("gaussian","rician"), vee=1)
Arguments
act.image |
Array defining where and when activation is present. |
sigma |
Standard deviation of the noise. |
type |
Distribution of task-related noise. Default is gaussian. |
vee |
If |
Details
The function generates random Gaussian noise for those voxels in the dataset that show activation. The result is a noise array with specified dimensions and desired standard deviation.
Value
An array containing the noise.
Author(s)
M. Welvaert
See Also
temporalnoise, lowfreqdrift, physnoise, systemnoise, spatialnoise
Examples
d <- c(10,10,10)
sigma <- 5
nscan <- 100
act <- array(rep(0, prod(d)*nscan), dim=c(d,nscan))
act[2:4,2:4,2:4,c(20:30,40:50,60:70)] <- 1
out <- tasknoise(act, sigma)
Generate temporally correlated noise
Description
Generates an autoregressive noise dataset with specified dimensions and standard deviation.
Usage
temporalnoise(dim, nscan, sigma, rho = 0.2, template, verbose = TRUE)
Arguments
dim |
A vector specifying the dimensions of a 2D or 3D array. |
nscan |
The number of scans in the dataset. |
sigma |
The standard deviation of the noise. |
rho |
The autocorrelation coefficients. The length of the vector determines the order of the autoregressive model. |
template |
An array representing the anatomical structure or mask with dimensions equal to dim. |
verbose |
Logical indicating if warnings should be printed. |
Value
An array containing the noise with dimensions specified in dim.
Author(s)
J. Durnez, B. Moerkerke, M. Welvaert
See Also
systemnoise, lowfreqdrift, physnoise, tasknoise, spatialnoise
Examples
d <- c(10,10,10)
sigma <- 5
nscan <- 100
rho <- c(0.3,-0.7)
out <- temporalnoise(d, nscan, sigma, rho, verbose=FALSE)