Slice-based regressors for physiological noise correction of 3D-EPI functional data

AFNI version info (afni -ver): AFNI_21.1.01 (Apr 9 2021)

Dear AFNI community,

My fMRI data were collected with a 3D-EPI sequence (Philips 7T scanner), with simultaneous recording of pulse and respiration signals. I am now trying to perform physiological noise correction by

  1. computing regressors via RetroTS.py, and
  2. performing image correction following the steps of the ricor block from afni_proc.py (3dDetrend, 3dDeconvolve, 3dREMLfit, 3dSynthesize, and 3dcalc).

Because the acquisition is 3D-EPI, I don't have slices in the conventional sense. For this reason, I first tried to compute slice-based regressors with RetroTS.py using -n 1 (number of slices = 1):

python RetroTS.py \
-phys_file sub-001_run-01_physio.tsv \
-phys_json sub-001_run-01_physio.json \
-n 1 \
-v 1.32 \
-prefix sub-001_run-01

However, 3dREMLfit failed because the number of regressors did not match the expected number, which corresponded to dim3 of the nifti image.

3dREMLfit \
-input sub-001_run-01_bold.nii.gz \
-matrix sub-001_run-01.ricor.xmat.1D \
-Obeta sub-001_run-01.ricor.betas \
-Oerrts sub-001_run-01ricor.errts \
-slibase_sm sub-001_run-01.ricor_det.1D
++ 3dREMLfit: AFNI version=AFNI_21.1.01 (Apr  9 2021) [64-bit]
++ Authored by: RWCox
*+ WARNING: Unknown regressor ordering. If the regressors
   were made via 'RetroTS' you are probably okay.
** interleave_cols: invalid nint=98 (ny=13)

When I ran RetroTS.py with a number of slices equal to dim3, 3dREMLfit succeeded.

At this point, my questions are:

  1. To bypass the 3dREMLfit error while considering that I do not have slices, would it be valid to compute regressors only for a single slice (-n 1 in RetroTS.py) and then replicate them across all slices along dim3?
  2. If I had to force the creation of regressors for each slice, which dimension should define the slice count in my case (i.e., is dim3 always the correct choice)?
  3. Related to question 2, I noticed that 3dREMLfit automatically matches slice counts between slibase regressors and dim3, and that the header of the detrended regressors (output of 3dDetrend) reports nifti orientation as:
#  ni_axes     = "R-L,A-P,I-S"

Does this mean AFNI assumes the I-S axis is stored as dim3? Is the anatomical orientation important in the RETROICOR procedure?

Any clarification on these questions would be very much appreciated!
Thank you,
Valeria

Hi, Valeria-

Great to have physio data. In the past few years, we have updated our physio-based regressor generation, by making a new program called physio_calc.py; please see it's help here or with physio_calc.py -hview. This updates and replaces the RetroTS.py functionality.

The new program has a lot more flexible functioning, better QC of intermediate steps, and checks for consistency with EPI time series properties, as well. I would strongly encourage moving to the newer program. Here are poster descriptions of it from OHBM-2024 and OHBM-2025 poster about it.

You might need to update your AFNI to get it, but since that computer's version is over 4 years old, that is probably a good thing to get many other updates from the meantime, too.

The physio_calc.py output can be put directly into afni_proc.py.

It would be better to troubleshoot this with the newer program (if the problem actually persists, in that case, that is).

--pt

Hi Paul,

Thank you very much for your quick reply and for pointing me to the updated resource!

I have tried creating my regressors with physio_calc.py using the -dset_epi option, and I still get a slibase with 98 slices, so I assume this is the only way I can perform my noise correction.

Thanks also for sharing the posters - very helpful indeed!

Best,
Valeria

Hi, Valeria-

What is the output of this for your EPI dset:

3dinfo -orient -n4 DSET_EPI

? The -n4 will output the 4 dimensional sizes of your data (the first three are spatial, and the last it number of volumes).

The 3rd dimension is assumed to be your number of slices. The dataset orientation is also output in the above command---the 3rd letter there should be the acquisition slice "letter", which would be "I" or "S" for axial slice acquisitions.

Were your slices not acquired axially, or does the orientation and matrix dimension information not match the order of slice acquisition?

--pt

Hi Paul,

Here is the output:

3dinfo -orient -n4 sub-001_run-01_bold.nii.gz
AIL 112 112 98 232

physio_calc.py correctly estimates regressors along the expected dimension for slices (the 3rd), but in my case this dimension does not correspond to actual slices, since the data were acquired with a 3D-EPI sequence. Would this affect the validity of regressors?

Thanks for your help!

Hi, Valeria-

The dataset orientation order should reflect the acquisition order, with the slice plane being the third one; that should also be consistent with the number of slice timing values, if you have slice timing available as well.

--pt

Hi Paul,

Sorry if I wasn’t clear enough previously. I believe the logic you described applies to 2D EPI, where one acquires a slice at a time until the full volume is covered. In my case, since I used a 3D acquisition, the entire volume is sampled at once. The slices only appear after the reconstruction, but they are not directly acquired. For this reason, I also don’t have any slice timing information.

That is why I was concerned about how to apply slice-based regressors here, given that the acquisition itself does not provide a true slice-by-slice structure. Does this make sense to you?

EDIT:
I looked more closely at the slice-based regressors generated with:

python physio_calc.py \
-phys_file sub-001_run-01_physio.tsv \
-phys_json sub-001_run-01_physio.json \
-dset_epi sub-001_run-01_bold.nii.gz \
-prefix sub-001_run-01

and it turns out that all regressors are identical across slices. This matches the dset_slice_times in the slibase JSON, which is a vector of 98 zeros. I think this resolves my issue, because the slice-based regressors are compatible with the dim3 dimension (= 98) of the nifti image, but the information is the same across slices.

Yeah, I don't believe you can use physio_calc (or any RETROICOR approach) with 3D-EPI. When I have 3D-EPI I remove the tshift block from my afni_procs, too. There's probably still a way you can use those physio measurements but it may involve some kind of down-sampling to the TR to use as GLM regressors of non-interest but I'm not sure how to do that well. You might have to go with ROI-based cleaning approaches like aCOMPCOR or ANATICOR. I'm responding to this because: (a) gov shutdown and (b) i'm really interested in this question and want to know if a physics/pulse seq expert chimes in. A lot of the 3D-EPI users I know are in the laminar world (super high resolution) and therefore their data is in the "thermal-noise-dominated regime" where physio correction isn't as critical so this issue doesn't always come up. -Sam

Hi Sam,

Thank you very much for your reply!

Some time ago, I ran some comparisons between performing physiological noise correction before any other preprocessing step using slice-based regressors versus after preprocessing using physio data as regressors of no interest in a GLM. It was quite some time ago, so I don’t recall all the details, but we observed better results when the physio correction was performed first, in terms of both ICA components and beta estimates of regressors of interest. That’s why we decided to adopt that approach.

I understand that slice-based regressors are not designed for 3D-EPI, which indeed suggests that other approaches might be more appropriate here. However, running the noise correction after preprocessing, including the signals as regressors of no interest, seems somewhat suboptimal to me, since they would correlate with motion regressors, and removing physiological noise before preprocessing should also help with better overall preprocessing.

I also believe that using slice-based regressors with the same estimates across slices (i.e., varying only across volumes) would make sense, as one is effectively regressing out the same information volume-by-volume. The results I obtained in the previous days with physio_calc regressors seem reasonable (at least to me).

I’d also be very interested in the opinion of a physics/pulse sequence expert, let’s see if someone joins this conversation!

Hi Valeria,

I have some experience with physiological noise correction, although I have not yet used physio_calc.py. (I have my own very old scripts to perform RETROICOR and just haven't switched to physio_calc.py yet).

If the 3D-EPI acquisition is acquired in a single shot (one RF excitation), then it should be possible to apply RETROICOR, just with the same timing across all slices. If the 3D-EPI is made up of multiple excitations (e.g. one RF pulse for each phase encoding step), then it would be difficult to apply RETROICOR since each voxel will have essentially been acquired at multiple time points relative to the cardiac cycle. This latter 3D-EPI is often referred to as a multishot 3D-EPI.

-Rasmus Birn

2 Likes