psadil
September 23, 2022, 12:31am
1
I’m trying to ensure that I’m describing fMRIPrep’s implementation of aCompCor correctly, but I am unsure about one step in the CSF mask generation. The boilerplate description of my run (fMRIPrep version 20.2.3, using FreeSurfer) says
For aCompCor, three probabilistic masks (CSF, WM and combined CSF+WM) are generated in anatomical space. The implementation differs from that of Behzadi et al. in that instead of eroding the masks by 2 pixels on BOLD space, the aCompCor masks are subtracted a mask of pixels that likely contain a volume fraction of GM. This mask is obtained by dilating a GM mask extracted from the FreeSurfer’s aseg segmentation, and it ensures components are not extracted from voxels containing a minimal fraction of GM. Finally, these masks are resampled into BOLD space and binarized by thresholding at 0.99 (as in the original implementation)
I think that the relevant dilation and extraction are done here (fmriprep/confounds.py at 469384ea47f1a7105705df0382ad3e6ddf116eca · nipreps/fmriprep · GitHub ), where I see a subtraction of gray from white matter and gray from the combined mask. But I’m not seeing gray subtracted from CSF.
When fMRIPrep generates the CSF mask for aCompCor, is the (dilated) gray matter mask subtracted from the CSF mask (e.g., somewhere else)?
Thanks!
jsein
September 23, 2022, 11:53am
2
Very good question! Indeed, looking at the code it seems that no GM subtraction is done to the CSF mask! Subtraction is only done to WM and combined WM+CSF masks in this code.
psadil
September 23, 2022, 1:15pm
3
Sounds good! Thanks for double-checking.
With the subsequent thresholding and binarization, I would believe that the subtraction would only barely change the CSF mask (if at all?). Also, the lack of a reduction of that mask’s size matches the Behzadi et al. implementation (they didn’t erode the CSF mask).
jsein
September 23, 2022, 7:08pm
4
Some history about the discussions that preceded this implementation in FMRIPREP:
nipreps:master
← oesteban:fix/2129-compcor-masks-reliability
opened 11:46PM - 20 May 20 UTC
This PR revises the implementation of CompCor masks in an attempt to make it clo… ser to the original proposal and at the same time address some recurring problems of *fMRIPrep*'s tCompCor implementation.
Finally, with the more careful resampling of prior knowledge from the anatomical scan, this refactor should also make the aCompCor components more run-to-run repeatable.
## aCompCor
The massaging of CompCor masks is now done in anatomical space where it is more precise, and a careful resampling to BOLD space follows.
The implementation deviates from Behzadi et al. Their original implementation thresholded the CSF and the WM partial-volume masks at 0.99 (i.e., 99% of the voxel volume is filled with a particular tissue), and then binary eroded that 2 voxels:
> Anatomical data were segmented into gray matter, white matter, and CSF partial volume maps using the FAST algorithm available in the FSL software package (Smith et al., 2004). Tissue partial volume maps were linearly interpolated to the resolution of the functional data series using AFNI (Cox, 1996). In order to form white matter ROIs, the white matter partial volume maps were thresholded at a partial volume fraction of 0.99 and then eroded by two voxels in each direction to further minimize partial voluming with gray matter. CSF voxels were determined by first thresholding the CSF partial volume maps at 0.99 and then applying a threedimensional nearest neighbor criteria to minimize multiple tissue partial voluming. Since CSF regions are typically small compared to white matter regions mask, erosion was not applied.
This particular procedure is not generalizable to BOLD data with different voxel zooms as the mathematical morphology operations will be scaled by those. Also, from reading the excerpt above and the tCompCor description, I (@oesteban)
believe that they always operated slice-wise given the large slice-thickness of their functional data.
Instead, *fMRIPrep*'s implementation deviates from Behzadi's implementation on two aspects:
* the masks are prepared in high-resolution, anatomical space and then projected into BOLD space; and,
* instead of using binary erosion, a dilated GM map is generated -- thresholding the corresponding PV map at 0.05 (i.e., pixels containing at least 5% of GM tissue) and then subtracting that map from the CSF, WM and CSF+WM (combined) masks. This should be equivalent to eroding the masks, except that the erosion only happens at direct interfaces with GM.
When the probseg maps provene from FreeSurfer's ``recon-all`` (i.e., they are discrete), binary maps are *transformed* into some sort of partial volume maps by means of a Gaussian smoothing filter with sigma adjusted by the size of the
BOLD data.
## tCompCor
In the case of *tCompCor*, this commit removes the heavy erosion of the brain mask because 1) that wasn't part of the original proposal by Behzadi et al., and 2) the erosion was the potential source of errors from numpy complaining that it can't take from an empty axis of an array.
> Based on these results, we chose a 2% threshold (∼20–30 voxels per slice) as a reasonable empirical threshold that effectively identified voxels with the highest fractional variance of physiological noise.
Although they do the calculation slice-wise, this commit rolls tCompCor back to calculate the 2% threshold on the whole-brain mask.
I would appreciate a deep review from @rciric & @jdkent to ensure my interpretation of Behzadi's paper is correct. I'm also tagging @effigies for a more top-level check of the changes. Further feedback from others is of course very much appreciated.
Resolves: #2129.
References: #2052.
This behavior of not reducing the CSF mask matches indeed the Behzadi et al. implementation:
Since CSF regions are typically small compared to white matter regions mask, erosion was not applied.