How do we get resting state fMRI time series for a collection of regions?

[isukrit]

Hi, I want to get simulated rs-fMRI time series for a number of regions in the same brain but there is no function defined for it. The 'simTSrestingstate' function gives us time series for a single region. Do you want us to run the same function k times for k regions and is that equivalent to rs-fMRI for the all regions in the same brain?

Thanks!

[HBossier]

Hi @isukrit,

Unfortunately, simulating resting state time series for multiple regions is not yet implemented (as far as I can tell, I am not the main developer). It is possible to use the same function k times, however you will end up with regions that are not correlated.

One workaround I see at the moment is to generate k time series without any noise. Then define your own variance-covariance matrix between the k regions, draw noise from a multivariate distribution and add this noise on top of the raw time series. For instance, say you want to generate rs-fMRI time series for two regions R1 and R2. We use a simple correlation between both regions of 0.75 and add Gaussian white noise with a signal to noise ratio of 1 (if we take a base signal of 100, we then get sigma = 100). Then I would do the following in R:

# Scanning parameters
nscan <- 100
base <- 100
TR <- 2
sigma <- 100

# Variance-covariance matrix between the regions (correlation of 0.75)
covM <- matrix(c(sigma, 75, 75, sigma), nrow = 2, ncol = 2, byrow = TRUE)
covM
     [,1] [,2]
[1,]  100   75
[2,]   75  100

# Raw time series
R1 <- simTSrestingstate(nscan = nscan, base = base, TR = TR, noise = 'none')
R2 <- simTSrestingstate(nscan = nscan, base = base, TR = TR, noise = 'none')

# White (correlated) noise
whiteNoise <- MASS::mvrnorm(n = nscan, mu = c(0, 0), Sigma = covM)

# Final time series for both regions
restRegion <- cbind(R1, R2) + whiteNoise

# Check correlation
cor(restRegion[,1], restRegion[,2])
[1] 0.7687157

# Plot time series
par(mfrow = c(1,2))
plot(restRegion[,1], type = 'l', main = 'Region 1 - rest',
     ylab = '', xlab = 'Scan')
plot(restRegion[,2], type = 'l', main = 'Region 2 - rest',
     ylab = '', xlab = 'Scan')

And generalise this for multiple regions.

[isukrit]

Very interesting! I will try this out. Additionally, you could use Cholesky Decomposition to achieve the same as well.