Singularity for an FMRI time series across subjects


I simulated fmri data via R’s neurosim. For 20 subjects, there was a 200-second experiment in which the onset (only a single stimulus) occurred at every 30 seconds, starting at 1s. TR=1, and I used a 64x64 grid, so there are 4096 voxels… I then convolved these onsets with the canonical HRF to get the evoked stimulus, and ran a mixed effects model: Bold=EV+subject, where subject is random. Unfortunately, some voxels, not all, return a singularity error. Does anyone know why I would get singularity issues in a design such as this? Thank you.

Ahoi hoi @JayBro94,

full disclosure: I’m by far no expert re these things,
thus not sure how helpful my comments might be, sorry.

Would you mind providing more information re your EV and voxels in question (e.g. does their simulated signal vary from those of voxels that don’t return the error, are they spatially adjacent/neighboring voxels)? IIRC, one would usually get this sort of error in cases where the predictor variables have a prominent linear relationship (e.g. multicollinearity) or the data doesn’t support a given model (which should then be reduced in complexity/be more parsimonious), eh?

Cheers, Peer

Hi @PeerHerholz ,

Thanks for the response. My onset times are every 30 seconds, starting at the 1st second. I have 20 subjects, and for each subject, the effect size varies from 0 to 4. I convolved this stimulus with the canonical hrf to get my EV.

I am finding that this error occurs for voxels that are near the region with the stimulus. But voxels that are far from the stimulus don’t give me an error. I did add some noise to the data to mitigate this, but it didn’t seem to help. Below is my R code of the simulation for reference. Thanks!

for (i in 1:num_sub){
design ← simprepTemporal(totaltime=time_len, onsets=seq(1,time_len,30),
durations=duration, TR=tr, effectsize=runif(1,0,4), hrf=“double-gamma”)
region ← simprepSpatial(regions=2, coord=list(c(runif(1,29,35),runif(1,13,15)),c(runif(1,55,57),runif(1,43,47))),
radius=c(7,4), form=“sphere”, fading=0.5)
out ← simVOLfmri(design=design, image=region, dim=c(64,64),
SNR=runif(1,0.5,3), noise=“mixture”,weights=c(1/6,1/6,1/6,1/6,1/6,1/6),rho.temp=0.2,rho.spat=0.2)

Just to clarify, the code above only simulates the times series data. After simulating this data, I create a vector for the EVs, which are the same for all subjects. So, I get a vector in which the first 20 elements are the same, then the next 20 are the same, and so on. This vector has (20 subjects * 200 time points) 4000 elements. I then run a mixed model as follows:

lme(unlist.get.nam…~1+EV_rep, random=~1 +EV_rep | idnum, data=get(nam),

This gives singularity for all voxels that are near the activated region.

Nevermind, I figured it out. The mixed model I was running was a random slopes model, and there was hardly any variation in the EV between subjects. If I run a random intercept model, it works. Thanks!