Setting SPM defaults
SPM has many default options, and mostly they're sensible. There are a few places where you can go very wrong by using the default values without thinking of them, however! We'll go through a few that are important to consider.
General Recommendations
- Standardize the default options for your analysis group (i.e., lab)
- Make sure that all computers used for analysis use the same defaults, to avoid unexpected results
- Don't rely solely on the defaults - some studies will require different choices, or selection of choices based on characteristics of the experimental design and data.
- Check the output of your SPM models, especially the mask of which voxels are analyzed and the filtering/autocorrelation options saved in the SPM.mat file in your first-level model output.
Which files to edit
There are a couple of options for setting defaults. First, you can edit the file spm_defaults.m in the copy of SPM you're running. Try:
Second, you can create a file called spm_my_defaults.m somewhere on your Matlab path.
What to change
Here is an overview of what we recommend changing. You'd edit these lines in spm_defaults.m, or include them in spm_my_defaults.m. We'll cover why to do this below.
global defaults % only need to add for the spm_my_defaults.m option
%==========================================================================
defaults.stats.fmri.hpf = 180; % tor changed; 128 can be too low; should choose yourself % 128;
defaults.stats.fmri.cvi = 'none'; % tor changed. 'AR(1)', 'None', 'FAST'
%==========================================================================
defaults.mask.thresh = -Inf; % tor changed; Was 0.8. 0.8 performed implicit masking and made explicit masking impossible (may be fixed?)
%==========================================================================
defaults.stats.maxmem = 2^33; % tor modified. 2^33 = 8.59 GB. Was 2^28; depends on your memory.
defaults.stats.topoFDR = 0; % tor modified. Was 1 for "topological FDR" on peaks.
Why?
Let's go through line by line and discuss why.
High-pass filtering
defaults.stats.fmri.hpf = 180; % tor changed; 128 can be too low; should choose yourself % 128;
The default is 128-sec high pass filter, which means that any design-related variance (aka "real signal") that is slower than 1/128 Hz = 0.0078 Hz will be filtered out. If you have a block design with 25 sec on/off alternating periodicity, this is fine. But if you have:
- on/off blocks that are 40 sec long (80 sec periodicity)
- 3 blocks that are 20 sec long and alternate pseudorandomly, so that 2 blocks of the same type may appear next to each other
- Long trials, e.g., a working memory study with a 20+ sec delay interval and a response probe like this one, or pain epochs separated by 30-60 sec like those included in this study
- Other designs with energy (variance) at slower periodicities (> about 50-60 sec)
...then a substantial proportion of your true signal will be filtered out!!
If you have a randomized event-related design in which multiple events may be clumped together, the default may be problematic! You can check this by looking at the cumulartive variance in your design at different frequencies and comparing it to the energy of the noise at those frequencies (e.g., scn_spm_design_check.m in CANlab Core Tools).
...but in general, 180 sec is a safer default.
Autocorrelation estimation
defaults.stats.fmri.cvi = 'None'; % tor changed. 'AR(1)', 'None', 'FAST'
Estimating the autocorrelation using a noise AR model or similar can improve efficiency (i.e., power) if the noise model is accurate. But if the noise model is misspecified (inaccurate) or its parameters are poorly estimated (imprecise), then it can be worse than doing nothing.
The default AR(1) model in SPM is a linear approximation of an AR(1) model that assumes autocorrelation is uniform everywhere in the brain (at least in regions that count as important). It's not. And short, noisy time series are not well suited to precisely estimating AR-model parameters.
In our experience, we've seen AR models perform very poorly in fMRI data.
The good news is that group analyses are valid and unbiased even if no autocorrelation model is used at the first level. This is because the betas (regression slopes) are unbiased estimators -- i.e., they are not systematically "wrong", and your group model will thus not be wrong. They simply lose efficiency (are more variable) in the absence of a good noise model.
Our default choice for SPM AR modeling is "none" (don't do it). You might consider using it if you:
- Want to make inferences about single subjects' activation and need to do this to have a hope of obtaining valid single-subject P-values
- Have lots of high-quality time series data per person (e.g., lots of within-peson data)
...but in these cases, consider using diferent software for this, like FSL.
Masking
defaults.mask.thresh = -Inf; % tor changed; Was 0.8. 0.8 performed implicit masking and made explicit masking impossible (may be fixed?)
By default, SPM masks out voxels that are below 20% of the global mean image value for a given subject -- that is, values that are not very intense in the images.
However, inhomogeneity in the raw values of fMRI images exists for lots of reasons. For example, some sequences tend to produce brighter values in the center of the brain (e.g., spiral sequences).
And even if overall signal intesity is low (e.g., in areas subject to susceptibility artifacts like the ventromedial prefrontal cortex, medial temporal lobes, and amygdala), if it correlates meaningfully with your task design in a way that is consistent across participants, there might be important true effects there. Low overall intensity makes it less likely you'll obtain strong statistical results, but does not invalidate positive statistical findings.
In SPM, voxels that are missing (not analyzed) for one subject are also missing in any SPM group analysis containing that subject. If different subjects have different missing voxels, you can end up failing to even analyze important parts of the brain in your group analysis. We've seen cases where the entire cortex was missing from analysis!
Some toolboxes don't have this problem and analyze available subjects at each voxel, like the CANlab robust regression (see also here) and mediation toolboxes. But even so, it's not a good idea to exclude participants based on low overall signal intensity. The covariation with your task is what counts. You can examine which voxels you analyzed by looking at the mask.nii image in the SPM first-level model output.
Memory usage
defaults.stats.maxmem = 2^33; % tor modified. 2^33 = 8.59 GB. Was 2^28; depends on your memory.
SPM was created when very little memory was available, so only part of a slice was loaded into memory at once. Now, we can easily load a whole subject's time series into memory on a standard laptop. Increase this to load more data in at once and speed up anaysis. Set it so that it uses much, but not all, of your computer's RAM.
False discovery rate correction type
defaults.stats.topoFDR = 0; % tor modified. Was 1 for "topological FDR" on peaks.
False discovery rate control is a great and simple way to eliminate most false positives while preserving sensitivity. It is always as liberal or more liberal than familywise error rate correaction, and provides a better picture of the overall pattern of what is active and what is not by reducing Type II errors.
The original formulation was voxel-wise FDR, which provides a voxel-level P-value threshold. This is straightforward. Later, Chumbley and Friston proposed topological FDR control (see also here), which controls the false discovery rate of "peaks". The basic idea is that signal is everywhere (a continuous field over brain space, with no truly null voxels), so the question is not which voxels are active, but which local "peaks" (topological high-points) are significant. ("Peaks" are different from "clusters" or blobs). Consider a mountain-peak analogy, where statistical strength is elevation and geographical location is brain space. Suppose that condition A is "pre-earthquake" and condition B is "post-earthquake". The topological question becomes "is there a mountain peak in condition A that is higher than in condition B". This is compared with the standard voxel-wise question of "is the ground higher in condition A than condition B".
Chumbley and Friston take an interesting philosophical stance; but if the premise that "signal is everywhere and there are no false positives" is accepted, it's unclear what the value is of assessing local activation peaks at all. It is also less straightforward in terms of interpreting which of a series of measures (voxels) passes a given statistical threshold. So, being pragmatists who prefer to stick to inferences about local signal strength rather than inferred topological features (at least for standard brain mapping purposes), we prefer to use the voxel-wise method.