> A suggested future addition, to make this more straightforward, is to have a 3dLM version of 3dLMEr (or even better, allow robust linear modeling as well).

In fact 3dMVM does allow voxel-wise covariate, but the limitation is that currently it can only take one voxel-wise covariate (you have two in this case). And it also does robust regression.

> For 2) Remove the population mean of 'brain' at the voxel level and call it BRAIN;
> Is that a 3dcalc command in which mean is the mean of all brain voxel b-values?

The mean/average should be performed across all the subjects for each voxel separately. You can use 3dMean or 3dTstat to compute the voxel-wise mean, and then use 3dcalc to subtract the mean from each subject's data.

Gang